m‘ h‘ -. ~ muwmx‘.‘
I v A: .

 

, ‘ “‘ I
‘W M ,......,...........
. J. h
. '8]

g.-. hw'h-‘
“h .— .

A.
,2 h-...."‘"""“’“

.
is L N:

.. m ~-2~3"u
‘.\:L§E§“I "r" ‘N: ~44:

7 “VII; d~1fi p,.
@274 ~*-::‘ '3
2%;

,5”..—
'3‘. '
m,»

iwéfu aw

- 4‘. Lv~
’ . ->
“db—“IMEJ it...“
:u‘nms-k

‘ | ‘M agnmf-l‘
4.. \ ~ '
.... J

lyb‘:-?~’{

La - 1

~ a. V“ ..

%“~:“- .( '

”3&1“:
m. “a”?

r

‘ .44 v Jl".-:.-v

u ..

«pr «4
v 1.: ~41- \Iw ‘
: 4 .(N .‘rhf'ww
kiwi. " ‘ -;g' ‘3
0...

r

av

:31:
4’; at; t
-“‘u If:

. vim]. '-

“-:lK/‘I~q l. t
nun .~-
, -: ~N

’s

‘ a ll ‘1
if $65."

i
.-‘~‘L

'5‘. .4.

«u

t‘ 25‘

r. .:.‘.,
. u. 4 . - 4 3
,,. . Jim '4

J
. "’r" 'v
”‘1- ‘. “.s. ..
"1:22;... a.
..'.'.':':’:”';.
. . .
,.

H,..;.’...p,,.

JUNE“
' 4

 

;
w in u ‘3‘...“41 ‘

,m

r . -~ - 7:3,", '32-.
x . J” .,...- .. é:
.I;..§:..'..J . 4- , , 3?;

m ”I. 1

.'.r. ’l.."..';"
macs-47 a.
l :4” 4M4

"1' '31-...
xfiflmh
1f.- '3,

 

:5" -

‘ .31; W

;.~:,.::.:E;"W'J’:flm" '. 4,}.-.wr51flrjgg
L inf.“ ’1?" 3w:

.3‘

IGAN TATE ITY LIBRARIES

lllllllllllllll llllllm l I I

31293 00906 3516

 

 

 

 

 

 

 

 

 

 

 

 

This is to certify that the

dissertation entitled

Acoustical Imaging by Broadband Signals

presented by

Nathan Nai-Hsien Wang

has been accepted towards fulfillment
of the requirements for

Ph . D degree in Electrical Engi :‘ac:sz'23:ir:~3

{gr/W

[Majogjrofessor

 

 

[hm October 21, 1991

 

MSU is an Affirmative Action/Equal Opportunity Institution 0-12771

1 l
LIBRARY
Michigan State
University

 

—r——

 

I

PLACE IN RETURN BOX to remove this checkout from your record.
TO AVOID FINES return on or before due due.

 

DATE DUE DATE DUE DATE DUE

l l

ILLJLJ
”—ljm
‘___ll_ c i,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ACOUSTICAL IMAGING BY BROADBAND
SIGNALS

By

Nathan N ai-H sien Wang

A DISSERTATION

Submitted to
Michigan State University
in partial fulfillment of the requirements

for the degree of

DOCTOR OF PHILOSOPHY

Department of Electrical Engineering

1991

620} 6 / "/U

ABSTRACT

ACOUSTICAL IMAGING BY BROADBAND
SIGNALS

By

Nathan N ai-H sien Wang

The conventional method for measuring attenuation and velocity of material gen-
erally involves tedious and /or repetitious procedures. Using a narrow acoustic pulse
(broad frequency spectrum) for transducer excitation, it is shown that the reflection
coefficients of different materials are frequency dependent. Measuring the reflection
coefficient as a function of frequency, the velocity-density product and attenuation-
density ratio can be obtained with high precision. This technique can be extended
to an n-layered structure. A number of different materials are identified based on
the dispersion of a narrow acoustic pulse. In addition to the theoretical derivation,
the acoustical imaging by broadband signals are constructed. Using hierarchical clus-
tering techniques, the imaging is further characterized. The techniques discussed in

this dissertation apply not only to nondestructive evaluation, but also to biological

diagnoses.

Copyright by
Nathan Nai-Hsien Wang
1991

To my parents and my Wife

iv

 

ACKNOWLEDGEMENTS

I would like to express my sincere appreciation to Dr. Bong Ho and Dr. H.
Roland Zapp, my advisors, for their guidance and support. I want to thank Dr. John
R. Deller and Dr. Joseph C. Gardiner for their help. Most importantly, I would like
to thank my parents, Mrs. and Mr. Chiou-Yueh and Chin-Chi Wang, and my wife

Ye-Min Nei. Because of their encouragement and support, I was able to write this
dissertation.

TABLE OF CONTENTS

LIST OF TABLES viii
LIST OF FIGURES ix
1 INTRODUCTION 1
1.1 Overview ..................................
1.2 Chapter Descriptions ........................... 3

2 MEASUREMENTS AND APPLICATIONS OF ULTRASONIC

SIGNALS 4
2.1 Measurements of Ultrasonic Velocity and Attenuation ......... 4
2.1.1 Measurement of Ultrasonic Velocity ............... 5
2.1.2 Measurement of Ultrasonic Attenuation Constant ....... 7
2.2 Nondestructive Evaluation of Materials by Ultrasound ......... 11
2.2.1 Nondestructive Evaluation .................... 11
2.2.2 Nondestructive Evaluation by C-scan .............. 12
2.2.3 Resolution in C-scans ....................... 14
2.3 Medical Ultrasonics and Tissue Characterization ............ 20
2.3.1 Acoustical Imaging ........................ 20
2.3.2 Tissue Characterization .............. ' ....... 21
3 BASIC PROPERTIES OF ULTRASOUND 25
3.1 Review of Basic Ultrasound Principles ................. 25
3.2 The Acoustic Plane Wave Equation ................... 26
3.3 Transmission and Reflection Coefficients ................ 31
3.4 Broadband Signal Probing Techniques ................. 34

4 DETERMINATION OF VELOCITY AND ATTENUATION BY
BROADBAND SIGNALS 37
4.1 Damping and Attenuation ........................ 38

vi

4.2 Derivations under the Assumption that the Attenuation Coefficient Is

Proportional to Frequency Squared ................... 42

4.3 Derivations under the Assumption that the Attenuation Coefficient Is
Proportional to the I-th Power of Frequency .............. 54
5 EXPERIMENTAL PROCEDURE AND RESULTS 59
5.1 Experimental Procedure ......................... 60
5.2 Experimental Results ........................... 63
5.3 Two Dimensional Imaging by Using Velocity-Density Product . . . . 72

6 ULTRASONIC IMAGES AND TISSUE CHARACTERIZATION

BY HIERARCHICAL CLUSTERING TECHNIQUES 79
6.1 Hierarchical Clustering .......................... 80
6.2 Material Identification by Ultrasonic Images .............. 87
6.3 Tissue Characterization by Ultrasonic Broadband Signals ....... 92
6.3.1 Ultrasonic Imaging of Tissues .................. 92
6.3.2 Feature Selection ......................... 93
6.3.3 Image Enhancement ....................... 96
6.3.4 Tissue Characterization by Hierarchical Clustering ...... 100
6.3.5 Data Reduction .......................... 102
6.3.6 Image Reconstruction ....................... 102
6.3.7 Cluster Validation ........................ 104

7 CONCLUSIONS AND DIRECTIONS FOR FUTURE WORK 110

7.1 Conclusions ................................ 110
7.2 Directions for Future Work ........................ 111
A PROGRAM LIST - SCANNING AND COLOR DISPLAY 112
BIBLIOGRAPHY 144

vii

3.1

5.1

6.1

6.2

LIST OF TABLES

Approximate values of ultrasonic velocities of various media ..... 27

Experimental results of measuring velocity and attenuation by broad-

band signals ................................ 71
The values of CPCC by different clustering methods for the 2-layered
model ................................... 88
The values of CPCC by different clustering methods for a section of
human brain with a hemorrhaged tumor ................ 108

viii

3.1
3.2

5.1

5.2
5.3
5.4
5.5
5.6
5.7
5.8
5.9
5.10
5.11

5.12
5.13

6.1

6.2

6.3

6.4

6.5
6.6

LIST OF FIGURES

Infinite homogeneous medium with uniform applied pressure ..... 28

Transmission and reflection at an interface ............... 32

Computer simulation of the reflected Gaussian pulse from four different
materials .................................. 61

The block diagram of the experimental arrangement for reflected mode 62

The procedure of the sampling program ................. 64
The procedure for data processing .................... 65
The incident signal and reflected signals from test materials ..... 67

The spectra of the incident signal and reflected signals from test materials 68
The reflection coefficients of test materials within their 3-dB bandwidth 69
The results of regression with 90% prediction confidence interval . . . 70

The block diagram of the scanning system ............... 73
A 2-layered 2-dimensional sample model ................ 75
Imaging of the 2-layered 2-dimensional sample model using velocity-

density product data ........................... 76
A 3-layered 2-dimensional sample model ................ 77

Imaging of the 3-layered 2-dimensional sample model using velocity-
density product data ........................... 78

Clustering result of a simple 2-layered model by single link clustering
method with 3 clusters .......................... 89

Clustering result of a simple 2-layered model by complete link cluster-

ing method with 2 clusters ........................ 89
Clustering result of a simple 2-layered model by Ward’s clustering

method with 2 clusters .......................... 90
Clustering result of a simple 2-layered model by UPGMA clustering

method with 2 clusters .......................... 90
Comparison of CPCC by the different methods for the 2-layered model 91
A section of human brain with hemorrhaged tumor .......... 95

ix

6.7

6.8

6.9

6.10

6.11

6.12

6.13

6.14

6.15

6.16

6.17

6.18

6.19

6.20

Ultrasonic imaging constructed using the total energy of the reflected
signals ...................................
Central frequency, peak frequency, and bandwidth of a typical ultra-
sonic spectrum ..............................
Ultrasonic image constructed using the peak frequency of the reflected
signals ...................................
Ultrasonic image constructed using the central frequency of the re-
flected signals ...............................
Ultrasonic image constructed using the 3-dB bandwidth of the reflected
signals ...................................
Ultrasonic image constructed using the cross correlation of the incident
signal and the reflected signals ......................
The original image obtained from a composite material with a hole in
the center .................................
The enhanced image by interpolation ..................
Projection of five features in a two-dimensional space .........
Clustering result of a section of human brain with a hemorrhaged tumor
by single link clustering method with 8 clusters ............
Clustering result of a section of human brain with a hemorrhaged tumor
by complete link clustering method with 4 clusters ...........
Clustering result of a section of human brain with a hemorrhaged tumor
by Ward’s clustering method with 9 clusters ..............
Clustering result of a section of human brain with a hemorrhaged tumor
by UPGMA clustering method with 3 clusters .............
Comparison of CPCC by the different methods for a section of human
brain with a hemorrhaged tumor ....................

95

97

98

98

99

99

105

105

106

106

CHAPTER 1

INTRODUCTION

1 .1 Overview

Ultrasonic techniques are becoming increasingly important in medicine [1, 2], non-
destructive evaluation [3], and many other applications [4]. Most ultrasonic imaging
systems use echo returns from boundaries of different acoustic impedances to show
material properties related to boundary variations. It is the purpose of the research
reported here to show that the material (or tissue) characteristics can be extracted
from the shape of the echo pulse when an extremely narrow pulse excitation is used.

Techniques for the measurement of ultrasonic velocity and attenuation of a
medium have existed for a long time, but none permits data over a wide frequency
range to be acquired rapidly. Usually a sequence of measurements is made at discrete
frequencies over the required range [1, 5]. This is time consuming, and a disadvantage
if the system being studied is likely to be changing with time [6]. The broadband
video pulse technique provides the best solution for this kind of problem.

The video pulse signal contains a broadband of frequencies. Each frequency com-
ponent propagates through the material with a different velocity and attenuation
depending upon the target and material properties. When all the reflected frequency

components are received by the receiving antenna (transducer), the output will be

 

a video pulse which resembles a spread version of the incident pulse. The reflection
of each spectral component is determined by the characteristics of the medium from
which it is reflected. If the bandwidth is sufficiently broad, many frequency compo-
nents can be obtained. Therefore, the frequency dependent velocity and attenuation
can be measured by a single pulse.

Although the video pulse technique has been used for a variety of applications
in target identification by electromagnetic waves [7, 8], it has not been applied to
acoustic probing in spite of the advantage of lower operation frequencies and lower
wave attenuation, since the attenuation coefficient is directly related to frequency.
In addition to the lower operating frequency, another attractive feature is that the
velocity of an acoustic wave ( about 1500 m / sec in water) is five orders of magnitude
lower than that of the electromagnetic wave in free space. Based on the relationship
between frequency and wavelength for any propagating wave, u = f /\ , the acoustic
wave should have much shorter wavelength which in turn provides superior range
resolution.

Since the acoustic wave is a pressure wave, the properties of the acoustic wave
are slightly different from the electromagnetic wave. For electromagnetic wave prop-
agation, the dominant characteristics of the medium are the conductivity and the
dielectric constant, while for ultrasonic waves, the propagation and reflection depend
on the medium density, elasticity, and viscosity.

Theoretically, the characteristics of a medium can be determined if the reflected
(or transmitted) and incident pulses are known. Conversely, if the incident pulse and
the characteristics of the medium are known, it is possible to calculate the reflected
(or transmitted) pulse shape.

Instead of using discrete frequencies for echo measurements, a broad bandwidth
acoustic pulse can be utilized to detect the dispersive property of a material for iden-

tification purposes. Each frequency component propagates through the material with

a different velocity and attenuation depending upon the properties of the medium.
Consequently, the velocity and attenuation constant can be measured by comparing
the pulse shape of the echo return to that of the known incident pulse. The acoustical

imaging by broadband signals can also be constructed.

1.2 Chapter Descriptions

Chapter 2 reviews the methods to measure velocity and attenuation by ultra-
sound and applications of ultrasonics. Chapter 3 provides the derivation of acoustic
waves and states the basic properties of ultrasound. Chapter 4 provides theoreti-
cal derivations of reflections with broadband signals and develops the velocity and
attenuation models for multilayer structures. Chapter 5 gives the experimental pro-
cedure to verify the derivations developed in Chapter 4. In Chapter 6, a section
of human brain with a hemorrhaged tumor is used for experimental measurements.
The acoustical images using various features are constructed and using hierarchical
clustering methods the image is characterized. Chapter 7 provides a conclusion and

some suggestions for future work.

CHAPTER 2

MEASUREMENTS AND
APPLICATIONS OF
ULTRASONIC SIGNALS

2.1 Measurements of Ultrasonic Velocity and At-
tenuation

The measurements of ultrasonic velocity and attenuation are very important in
ultrasonic applications and many methods for measuring these properties have been
developed. Although each method has its own advantages and disadvantages, the
choice of optimum method depends on the nature of the particular material and the
desired measurement accuracy. In the following sections, we will describe some ap-

proaches for measuring ultrasonic velocity and attenuation.

2.1.1 Measurement of Ultrasonic Velocity

Pulse-Superposition Method

The pulse-superposition method allows a very accurate measurement of sound
velocity [4, 9, 10]. The radio-frequency (rf) pulses excited by a transducer, are sent to
a specimen, from which multiple echoes return. Each succeeding echo is constructively
added to the earlier echo of a given pulse by controlling the pulse repetition rate
based on the reciprocal of the travel time in the target. The technique is called
pulse superposition because the ultrasonic pulses are literally superposed on each
other. The high degree of accuracy arises from the fact that the method is capable
of precisely measuring the time between cycles. The velocity is twice the thickness of

the target times the inverse time separation between pulses.

Sing-Around Method

The sing-around method [4, 10, 11] uses two transducers to measure the velocity
of sound by placing one transducer at each end of the test material. One transducer
is used as the transmitter the other as the receiver. The receiving ultrasonic pulse
transducer generates a trigger signal for the transmitting transducer in order to es-
tablish a pulse repetition rate, (PRR). The velocity can therefore be measured by

the following equation:

2) _ length x PRR
— l—e >< PRR

 

R5 length x PRR (2.1)

where PRR is the pulse repetitions rate per second, length is the length of the test
material, and e is a correction factor for delays in the transducers, in the coupling
between the transducers and the material, and in the electrical components. The

measurement of sound velocity by this method is of moderately high accuracy.

Interferometer Method

The interferometer [4] is a continuous-wave (CW) device which has been used
for accurately measuring velocity and attenuation of sound in liquids and gases in
which standing waves can be sustained. The transducer is fixed at one end of a
fluid column with a movable rigid reflector at the other end. The reflector is moved
by a micrometer adjustment mechanism. For a fixed frequency as the reflector is
moved, the reflected wave generates in phase and out of phase components which add
constructively or destructively with the transmitted wave. The half wavelength of
the particular sound frequency is determined by the distance the micrometer moves
between two successive maxima. The accuracy of the measurement thus depends
on the accuracy of the micrometer readings, the parallelism between reflector and

transducer surfaces, and the accuracy of the frequency determination.

Resonance Method

The resonance method is similar to the interferometer method. It can be applied
to measure the velocity of sound in gases, liquids, or solids. The resonance method
uses either one fixed transducer and one fixed reflector, or two transducers located
by a known distance apart. By changing the frequency of the ultrasonic wave, the
successive resonances can be determined. The difference between two successive res-
onant frequencies in a nondispersive medium is equal to the fundamental resonant
frequency of the medium. From the resonance frequency difference it is possible to

determine the sound velocity from the relationship for a single transducer:

1) = 2lAf (2.2)

where l is the distance between the transducer and the reflecting surface and A f is

the difference between successive resonant frequencies. The accuracy of the resonance
method depends on the accuracy of the frequency determination and the measurement

of the distance between reflecting surfaces.

Other Methods of Interest

In addition to the above methods, there are several methods in use for measuring
the velocity of acoustic waves in special cases [12]-[16]. For examples, McSkimin [17]
and Papadakis [10] introduced a method which is particularly good for rf measure-
ments in thin specimens such as rare single crystals and thin sheet metal. Dameron
used an inhomogeneous media model, instead of the traditional layered model, to de-
termine the acoustic velocity in tissue [12], and Ophir et al. presented a pulse-echo

beam-tracking method and transmission method to measure the speed of sound in

tissue [18]- [21].

2.1.2 Measurement of Ultrasonic Attenuation Constant

The ultrasonic wave propagating in the z—direction can be defined as:
A(t, z) = Age-o" cos(Kz — wt) (2.3)

where t is time, A is the amplitude, A0 is the amplitude at t = 0 and z = 0, K is the
propagation constant, and a is the attenuation constant.

The attenuation constant a can be determined from:

a _ 111(Al/A2)

— (22 _ 21) (2.4)

where A1 and A2 are the amplitudes measured at position 21 and 22 at times of

maximum amplitude. The unit of the attenuation constant is nepers per unit length.

There are many methods developed to measure the attenuation constant and some

will be outlined below.

Spectral-Difference Approach
Kuc introduced the spectral-difference method to measure the attenuation in bi-
ological tissue [1], [22]-[24]. If it is assumed that tissue acts like a linear system, the

two-way power transfer function, [H ( f )|2, is given by:

PF(f)

|H(f)|”= Pm)

 

(2.5)

where Pp( f) and PN( f) are the spectra of the reflected signals from the far and near
surfaces respectively.

Since the acoustic attenuation coefficient of most soft biological tissue has a linear-
with-frequency attenuation characteristic, the attenuation coefficient at frequency f,

a( f), can be expressed as:

0(f) = flf (dB/cm) (2-6)

where ,3 is the slope of the attenuation versus frequency.

The power transfer function of biological tissue can thus be expressed as:
]H(f)]2 = 10-a(f)(2D)/10 = 10—0.2fifD (2.7)

where 20 is the additional path length through the tissue traveled by the far surface
reflection.

If we use the log-spectra, we have

1010g10|H(f)|2 = 1010g10PF(f) — 1010810 P~(f)

= 4st (2.8)

Therefore, the value of B can be found by

‘10108101H(f)l2
2Df
if.

= 20 dB/(cm - MHz) (29)

 

where Sf is the slope of the attenuation with respect to frequency of the log—spectral

difference.

Spectral-Shift Method

In addition to the spectral-difference method, Kuc introduced the spectral-shift
method to measure the attenuation in biological tissue [22]- [24]. The spectral-shift
approach requires that the propagating pulse have a Gaussian-shaped spectrum, while
the spectral-difference method did not require a specific spectral form.

Using the same notation as in section 2.1.2, the spectrum of the reflected signal

from the near surface, PN, will have a Gaussian shape given by:
PM) = CNe-U-f~>’/B’ (2.10)

where CN is a constant, B is a measure of the pulse bandwidth, and fN is the centroid
of the spectrum, or central frequency.

The two-way power transfer function IH ( f )|2 is given by:
|H(f)|2 = 6‘4”" (2.11)

where 6 is the slope of the attenuation coefficient having units of nepers per

centimeter-megahertz. The spectrum from the far surface reflection, Pp( f), is:

PF(f) = Cpe‘("’F’2/B’ (2.12)

10

where CF is a constant, and

fp = fN — ZBDB2 (2.13)

The spectrum of the far surface reflection still has the same Gaussian form as
the spectrum of the near surface reflection, but it has been shifted to a lower center-
frequency. Therefore, the value of B can be determined from the frequency shift in

the central frequency:

_fN"'fF

Since 1 Np/(cm-MHz) = 8.686 dB/(cm-MHz), 3 can be expressed as follows:

_ 4-343(f1v - fF)
— D32

 

fl dB/(cm-MHz) (2.15)

Other Methods and Applications

Many papers discuss the measurement of the attenuation constant in different me-
dia [25]-[29]. Since the attenuation constant is one of the most important properties
of materials (or tissue), many applications use the attenuation constant to classify
the medium, such as whether livers are normal or abnormal [30]-[34]. For example,
Matsuzawa et al. developed a technique for measuring both velocity and attenuation
of sound in highly absorptive materials [25], Parker and Waag proposed a method to
measure the ultrasonic attenuation from B-scan images in biological applications [31],
and Wang et al. used the video pulse technique to investigate velocity and attenuation

in homogeneous models [35]- [40] as discussed in Chap. 4.

11
2.2 Nondestructive Evaluation of Materials by

Ultrasound

2.2.1 Nondestructive Evaluation

Nondestructive evaluation (NDE) of materials by ultrasound has experienced
tremendous growth in recent years, especially for composite materials [4, 41, 42].
The development of reliable test methods for the nondestructive evaluation of com-
posites is a challenging task. In the past ten to fifteen years, a great deal of research
has focused on the development of quantitative active ultrasonic testing as well as
passive acoustic emission techniques.

Measurements by ultrasound can provide a great deal of information on the prop-
erties of materials. Measurement techniques utilizing ultrasonic waves are especially
attractive because of the direct connection between the characteristics of the wave
propagation and the mechanical properties of a material. The ultrasonic waves can
propagate through a material, carrying the information of the samples out without
damaging the samples.

In general, a successful ultrasonic materials characterization procedure requires
the solution of two problems [42]. The first is related to the reliable detection of
the ultrasonic properties, such as amplitude, arrival time, velocity, attenuation, or
spectral feature, etc. The second is related to the evaluation of the waveform pa-
rameter used to recover the material properties. One can determine the material
properties from some particular waveform parameters by establishing a correlation

between them.

12

2.2.2 Nondestructive Evaluation by C-scan

Time-of-Flight C-scans

Conventional ultrasonic C-scans display the amount of energy reflected from some
feature in the sample as a function of two dimensions at some fixed depth [3]. Com-
mon usage of the C-scan is to locate delamination and other defects, such as porosity,
or cracks in materials. Therefore, conventional C-scanning is widely used for quality
control during laminate manufacture. However, conventional C-scanning does not
provide information on the depthwise distribution of damage. The consequences of
having a defect localized between two plies may be quite different from having a
uniform distribution of defects throughout the laminate, so that the depth informa-
tion may be quite important. The depth information can be obtained by recording
the position of defects in time by an ultrasonic A-scan and using this time-of-flight

information to construct a three dimensional C—scan.

Computer-Automated Measurement and Control

Computer-automated measurement and control (CAMAC) is an instrumentation
and interface standard that provides both hardware and software flexibility [43].
By use of the standard CAMAC system a relatively easy interfacing of the motor
controllers to the transient recorder for complete computer control is possible. The

advantages of the CAMAC system are:

1. Many modules are available, such as analog-to—digital converters and motor

controllers.
2. The interface software can be written in a high-level language.

3. New custom applications are easily established by changing modules and writing

simple interface routines in high level languages.

13

C-scan by High Speed Sampling

At present, C-scan ultrasonic imaging is used routinely in examining composite
materials [41]. It is well known that there are many attractive features for such a
nondestructive evaluation technique. By using the C-scan it is possible to precisely
identify voids and defects at various depths inside the material. In practice the
following drawbacks limit the application for detecting flaws and defects at various

depths inside the material [41]:

l. The range resolution is restricted by the detected pulse width. Improvement can

be made by using large amplitude narrow pulses, which however is technically

difficult to achieve.

2. It is impossible to obtain a complete impact damage display of a specific fiber
layer by the constant-depth scan (C-scan) due to the fact that the fibers at the

damage site have been displaced.

3. To detect flaws, all echo returns from the material interior must be processed.
This requires a great deal of processing time and thus impedes the effort of

rapidly detecting flaws over a large object such as the wing of an aircraft.

By recording the detailed signal waveform and performing some signal processing,
the shortcomings mentioned above can be removed. Specifically, the shortcomings
can be addressed as discussed below.

Theoretically, the range resolution of an ultrasonic imaging system can approach
the operating wavelength. Most commercially available imaging systems offer sub-
stantially less resolution than is achievable based on the theoretical limit. By using
pulse amplitude detection of echo returns to identify the acoustic boundaries, an un—

certainty in range is introduced because of the wide pulse width and the requirement

14

to exceed a threshold setting in the receiver system. The precision in locating bound-
aries is further deteriorated by the timing jitter inherent in the transmitted pulse
signal.

The difficulty of retrieving the phase information from the rf carrier is due to the
high sampling rate required to satisfy the Nyquist criteria. The sampling rate should
be 4 or 5 times greater than the signal frequency, or more than 10 MHz for the system
utilized in this research. A very high speed real time recording system is necessary
to achieve the high sampling rate. A system was designed which involved a very high
speed sampling circuit, a memory unit to store the sampled data and the appropriate
signal processing software to identify the boundaries. In conjunction with the precise
time reference from a fixed spatial reference, the echo returns can be orderly sorted
by the software displayed. The unique feature of this signal processing scheme is that
the data points can be lined up by layer structure rather than at constant depth from
the transducer. Since the retrieved data are stored in groups by layer as well as by
the spatial location from the transducer, a three-dimensional contour plot of each

layer becomes possible.

2.2.3 Resolution in C-scans

One important topic for C-scan images is the resolution. Resolution refers to the
ability to separate echoes from two discontinuities located very close together [4].
To increase the dynamic range usually requires digital rather than analog techniques
[2]. The resolution of a C-scan image can be divided into range resolution and lateral

resolution. Each of these are considered below.

Range Resolution

The range resolution is the ability to distinguish two different echoes in depth.

Since the time axis represents the depth from which the ultrasonic signals were ef-

15

fected, the problem of range resolution translates to the ability to distinguish signals
along the time axis. When a second echo arrives at the receiver before the end of the
first these signals will overlap and make it very difficult to distinguish between the
two signals. Typically the range resolution is limited by the transmitting pulse width
and operating frequency [4]. Pulses with long pulse width will cause poor resolution,
while short pulses produce high resolution. In order to develop a high resolution
system, a broad bandwidth, low-Q transducer is necessary. Many researchers have
tried to increase range resolution by different methods, for example, Lees [44] used
the peak ratios of ultrasonic signals to determine the thickness of thin layers. By
measuring the peak ratios of the thin layers of known thicknesses, the thickness of
unknown layers could be determined by interpolation. This method provides a simple
approximation of the layer thicknesses.

In addition to the methods discussed above, another interesting method to increase

the range resolution is by deconvolution, as discussed below.

Increasing Range Resolution by Deconvolution

Although the range resolution is related to the ultrasonic wavelength, it is still
possible to increase the range resolution beyond the wavelength limitation. The key
point is to carefully record the phases of the signals. If it is possible to distinguish
the signals of different phases, the resolution is limited only by the identification of
phase difference. A powerful tool to assist in phase discrimination and thus to help
to increase the range resolution is signal deconvolution.

For a linear time invariant system, deconvolution can be used in many applica—
tions. There are many papers which discuss the applications of deconvolution to
ultrasonic signals [45]-[49]. Beretsky and Farrell used frequency domain deconvolu-
tion to improve ultrasonic imaging [45]. Deconvolution in the frequency domain is

straightforward and easy to implement. Assume an incident ultrasonic signal, a:(t),

16

and a reflected signal, y(t) with an impulse response of the test material given by

h(t). If the reflection process is linear, the relationship between x(t) and y(t) is:

W) = $(t)*h(t)

= /: h(r)x(t — 7')d7' (2-16)

where (at) is the operation of convolution.
If the Fourier transformations of a:(t), y(t), and h(t) are X(w), Y(w), and H(w),
respectively, the above relationship between the incident signal and the reflected signal

can be written as:

Y(w)

H‘“) = m

(2.17)

Therefore, the impulse response of the system can be easily obtained by finding

the inverse Fourier transformation of H (w):

h(t) = f“{HWH
-1 2122
s {W} (2...)

Although this method is straightforward and well known, it has a major difficulty,
when the frequency spectrum, X (w), is close to zero at some frequencies, which results
in the denominator of Eq. 2.17 being close to zero. This will cause substantial error
in the calculation of H (w) Thus the spectral technique is not usually practical.

Steiner et al. introduced a generalized cross-correlator to improve resolution in de-

lay estimation measurements [47]. Yamada presented an on-line deconvolution for the

17

high resolution ultrasonic pulse-echo measurements using a narrow-band transducer
[48]. The on-line deconvolution filter is optimized so as to maximize the resolution
capability subject to the constraints of small processing noise.

Silvia summarized many methods of deconvolution [49]. For a discrete-time
linear time-invariant system, the deconvolution problem can be solved by a matrix

formulation. Assume the sequence of the incident signal is:

93(2) #0 i=0...N—1 (2.19)

a:(i) =0 i<00riZN

and the sequence of the impulse response is:

h(j) #0 j=0...M-1 (2.20)

h(j) =0 j<Oorj2M

and the sequence of the reflected signal is:

y(k) ,eo k=0...N+M—2 (2.21)

h(k) =0 k<00rk2N+M——2

The relationship between :L'(n) and y(n) can be written as:

M = .(n) * h(n)
M—l
= Zh(k)a:(n—k) n=0...N+llI—2 (2.22)

k=0

18

In matrix format, the convolution can be written as:
y = X’H (2.23)

where
31(0)

yll) (2.24)

 

 

y(N+M—2)

q

.2:(0) a'(—1) :r(—M +1)

X = :z:(1) a:(0) . . . a:(—M + 2) (2.25)

 

 

_x(N+M—2) :z:(N+M—1) .r(N—1)_

and

2(0)

’H = hi1) (2.26)

 

 

_h(M—1)_

yisa(N+M-1)bylmatrix,Xisa(N+M—1)bmeatrix,and'HisaM

by 1 matrix. Therefore, the impulse response sequence ’H can be obtained from:
”H = X ‘13? (2.27)

To solve the above equation is usually very difficult because the sizes of y, X,
and ’H are usually quit large. Fortunately, many elements in ’H can be set to zero by

checking the reflected signals y(n). This will greatly reduce the size of matrices and

19

make the equation solvable.
If there are only K nonzero elements in ’H, where K is less than or equal to M,

then Eq. 2.23 can be rewritten as:

y = XrHr (2.28)

where ’H, is a reduced matrix from ’H obtained by deleting some elements with zero
amplitude, and X, is a reduced matrix from X obtained by deleting some columns
corresponding to the zeros in ’H. Thus, y is still a (N + M -— 1) by 1 matrix, X, is a
(N + K — 1) by K matrix, and ’H, is a K by 1 matrix.

In Eq. 2.28, there are K unknowns, but we have (N + M — 1) equations. Usually
no exact solution can be found since the number of equations is greater than the
number of unknowns. However, The solution with the minimum least squares error

may be obtained by [50]:

H, = (aft) " x332 (2.29)

where X; is the transpose of X,.

Lateral Resolution

Lateral resolution is the ability to distinguish two different objects in the trans-
verse spatial space. Because of the poor focus property of ultrasonic transducers,
the lateral resolution is usually low and directly related to the size and frequency
of the transducer. To improve resolution, Hundt and Trautenberg used a special
digital filter to increase the lateral resolution of B-scan images [51]. The filter is
constructed from measured signal amplitudes across the transmitter beam providing
a lateral resolution improvement of 50 percent.

Sato, et al., used a homomorphic transform and deconvolution to increase the

20

sector-scan image [46]. Yokota and Sato developed a super-resolution ultrasonic imag-
ing system by using adaptive focusing [52]. By their method, a set of reflected wave
field data from an object was acquired by repeating the transmission and reception
for all pairs of transducers on the array. From these data the object—dependent adap-
tive focusing is performed to illuminate the desired points on the object in order to

suppress the signal.

2.3 Medical Ultrasonics and Tissue Characteri-
zation

Ultrasound has been widely used in medical applications since the 1970’s and the
potential for ultrasound utilization is now recognized [2]. Medical ultrasonics is in
a period of rapid growth and on the verge of making significant impacts on clinical
medicine. Due to fast application specific digital circuit techniques and powerful com-
puters, the applications of ultrasound becomes practical. Many ultrasonic systems

and automatic tissue classification systems have recently been developed.

2.3.1 Acoustical Imaging

Ultrasound Volumetric Imaging System

Shattuck et al. introduced a parallel processing technique, called explososcan,
for high speed ultrasound imaging with linear phased arrays [53]. This technique
allows the data acquisition rate to increase by a factor of four by the simultaneous
acquisition of four B-mode image lines from each individual broadband transmitted
pulse. By using the higher data rate, averaging of the image data to reduce noise is

practical.

21

Smith et al. has built a real-time volumetric ultrasound imaging system for med-
ical diagnosis at Duke University [54, 55]. The system uses pulse-echo phased array
principles to steer a 2 dimensional array transducer of 289 elements in a pyramidal
scan format. By using the explososcan technique [53], the system can produces 4992

scan lines at a rate of approximately 8 frames / sec by their parallel processing scheme.

Ultrasonic Tomography

In order to obtain quantitative images representing certain tissue properties, ultra-
sonic computerized tomography has been utilized by Greenleaf [56, 57], Ermert [58],
and others.

Ermert and Rohrlein built an ultrasonic reflection mode computerized tomography
imaging system by utilizing a conventional B-scanner in a multi-view operation mode
[58]. The superposition of several B-scans from different aspect angles of one cross-
sectional area and subsequent 2—dimensional processing procedures used in X-ray

CT-scanning produced high quality images of the tissue samples.

Another Method of Interest

In addition to the above ultrasonic systems, Richard used time-gain correction
methods to enhance ultrasonic B-scan images [59]. This technique achieves signal
multiplication by first logarithmically compressing the complete dynamic range of the

received echo.

2.3.2 Tissue Characterization

The literature contains many papers on tissue characterization by pattern recog-
nition and segmentation techniques. The objective in this research are to develop
automatic detection systems. Usually the choice of classification method is depen—

dent on the object of interest. Feature selection can improve the results impressively.

22

Tissue Characterization by Bayesian Classifier

Momenan et al. applied pattern recognition techniques to ultrasound tissue char-
acterization [60]. They used the Bayesian classifier to determine the feasibility of
classifying tissue type as well as determining some properties of the feature spaces.
An unsupervised clustering technique, hypothesis, was discussed also. Four features
are used in their paper, the average distance between regularly positioned specular
scatterers, the ratio of specular to diffuse backscatter intensities, a measure of the
variability in the specular component normalized by the diffuse contribution, and the
slope of the attenuation coefficient as a function of frequency. The first three features
can be obtained directly from the statistics of the squared envelope (intensity data),

and the fourth feature can be obtained from the frequency dependent attenuation.

Detection of Tumors by Mammograms

Besides the ultrasonic images, Brzakovic et al. presented an automated system for
detection and classification of particular types of tumors in digitized mammograms
[61]. They used two stages to analyze mammograms, first the system identifies pixel
groupings which may correspond to tumors. Second, the detected pixel groupings are

subjected to classification. The features used in their work are:

1. area - total number of pixels enclosed in an extracted region.

2. shape descriptor - least mean square difference between the signatures of the

tumor model and the extracted region.

3. edge distance variation descriptor - mean square distance between the edge of

the tumor model and the locally detected edge from the extracted region.

4. edge intensity variation - local intensity variation along the edges of the ex-

tracted region.

23

They used a classification hierarchy to detect tumors. The first measured feature,
which is the easiest feature to measure, can determine whether the samples are non-
tumor or tumor. It is only necessary to measure a second feature when the output
of the first classification is tumor. The Bayer classifier will be used to perform the
second classification. The classification by a third feature is only used if the output
of the second feature is abnormal.

Lin et. a1 presents a method for detecting a type of breast tumor, or circumscribed
mass from mammograms [62]. It relies on a combination of criteria used by experts,
including the shape, brightness contrast, and uniform density of the tumor areas.
They used modified median filtering to enhance mammogram images and template
matching to detect breast tumors. In the template matching step, suspicious areas
are picked by thresholding the cross-correlation values and a percentile method is
used to determine a threshold for each film. There are two possible approaches to
enhancing mammographic features. One is to increase the contrast of suspicious areas

and second is to remove background noise.

Detection of Lung Tumors

Hashimoto et al. used a two-feature classification technique for the detection of
the edges of x-ray images of lung tumors [63]. The features in this classification
technique are the outputs of various forms of gradient and Laplacian operators. They
found three two-feature classifiers that are statistically superior to a single feature

detector based on the gradient alone in the presence of noise.

Segmentation of Magnetic Resonance Imaging
Magnetic resonance imaging is a powerful tool in medical diagnosis. Bomans et
al. applied three dimensional segmentation to magnetic resonance imaging [64].

Raya presented a rule-based low-level segmentation system that can automatically

24

identify the space occupied by different structures of the brain via Magnetic Resonance

Imaging [65].

Other Methods of Interest

In addition to the imaging and detection procedures outlined above, other ap-
proaches are possible. For example, Trivedi et al. used a two-dimensional (gradient,
density) feature space for the segmentation of medical images [66], Klingler et al. used
mathematical morphology to classify the echocardiographic images [67], Michael and
Nelson introduced HANDX, a model-based computer vision system for automatic
segmentation of bones from digital hand radiographs [68], and Rosenfeld reviews
some basic computer vision techniques and speculates about their possible relevance

to the modeling of human visual processes [69].

CHAPTER 3

BASIC PROPERTIES OF
ULTRASOUND

This chapter provides the theoretical analysis for the ultrasonic wave equations
and acoustic broadband signal probing techniques [35, 70, 71]. In section 3.1,
some related ultrasound principles are reviewed. The acoustic plane wave equation is
derived in section 3.2. In section 3.3, the transmission coefficient and the reflection
coefficient are determined from the boundary conditions at interfaces. The properties

of the acoustic video pulse will be explained in section 3.4.

3.1 Review of Basic Ultrasound Principles

Ultrasound is the name given to those waves with frequencies above 20 KHz (above
the audible range). Ultrasonic waves consist of propagating periodic disturbances in
an elastic medium where the particles of the medium vibrate about their equilibrium
positions either perpendicular to or parallel to the direction of propagation. Since
the longitudinal vibrations are dominant and of lower propagation losses only these
variations are considered. The propagation of the resulting motion-strain effects away

from the source results in a longitudinal compression wave that transmits mechanical

25

26

or acoustic energy away from the source. The vibratory motion of the medium is
strongly dependent on the ultrasound frequency and on the state of the medium.
Table 3.1 provides the values of velocity of longitudinal sound waves in various
materials [4, 71]. Velocity in solids are the highest, those in liquids and biological
soft tissues are lower, and those in gases are lowest.

The relationship between frequency and wavelength of a wave is given by

A = (3.1)

u
7
where u is the propagation velocity, f the the frequency of the wave, and A is the
wavelength. Since the propagation velocity of ultrasound (1.5 X 103 m / sec in water)
is much lower than that of an electromagnetic wave (3 x 108m / sec in free space), the
wavelength of ultrasound is five orders of magnitude smaller for a given frequency. For
a frequency of 5 MHz, the wavelength of ultrasound is approximately 0.3 mm while
for an electromagnetic wave it is 60 meters. Because of the much shorter wavelength,

the use of ultrasound for detection will give much higher resolution.

3.2 The Acoustic Plane Wave Equation

In this section, we will investigate sound propagation in an infinite homogeneous
medium in equilibrium. The coordinate system is shown in Figure 3.1. If a force
is applied in the :r-direction, it will result in a uniform pressure p(a:o, t) to the y - 2
plane at a distance 270 from the origin. Due to this applied force, the particles will
experience a displacement in the x-direction. Particles at the location so will also be
displaced. If we consider a small differential volume element with incremental length

dz and area A in the plane of the applied force, the equilibrium volume, V, of the

Table 3.1. Approximate values of ultrasonic velocities of various media

 

 

 

L

Longitudinal Acoustic
Material Mass Density Velocity Impedance
leg/m3 m/sec x 106 kg/m2 - sec
Aluminum 2695 6350 17.1
Copper 8900 4700 42.0
Gold 19300 3240 62.6
Ice 900 3980 3, 6
Iron(steel) 7830 5950 46.6
Lead 11400 1960 22.3
Nickel 8800 5600 49.0
Plexiglass 1182 2670 3.17
Silver 10500 3700 39.0
Zinc 7100 4170 29.6
Acetone 790 1174 .93
Alcohol (ethyl) 790 1170 .92
Oil, SAE 20 870 1740 1.50
Water(20°C ) 1000 1480 1.5
Air (0°C) 1.293 332 428 x 10'6
Air (20°C) 1.21 340 411 x 10-6
Carbon dioxide (16.6°C) 1.935 277 536 x 10‘6
Hydrogen (199°C) 0.09 1265 115 x 10‘6
Nitrogen (19.9°C) 1.207 349 421 x 10'6
Oxygen (196°C) 1.38 327 451 x 10‘6

 

element will be

 

z
A \ ' \ '
—' é r— —"§+d§ r—
I I
13% .t) ]u(xo .t) E
—> -:-> :
l l
v i l
: :
0 x0 : xo-I—dx :
I I
: :3"3X
I I
I I

 

 

 

Figure 3.1. Infinite homogeneous medium with uniform applied pressure

As a result of this compressional force the differential volume changes to

V’ = A(d:c + d{) = A(d:c + git) (3.3)
where 015 represents a compressional displacement.
Therefore, the change in volume is
06
dV = V' — V = A —d .
(3.“) a4)

ConVentionally, the strain produced in the volume element is defined as the ratio of

29

the volume change to the original volume:

1L2:
V—Ba:

strain E (3.5)
According to Hooke’s law, the ratio of stress to strain in an elastic medium is constant.

This relationship can be expressed as

p = — g (3.6)
where k is the coefficient of elasticity. The negative sign in Eq. 3.6 results from a
positive pressure in the x-direction giving a strain in the negative x—direction.
We can define the particle velocity u as the time rate of change of particle dis-
placement. That is,
36

u(a:,t) = 52 (3.7)

The partial derivative of Eq. 3.6, with respect to time, gives:

8p __ 8 66 Bu

01'
an _ —10p
5; — TE (3.9)

If the applied pressure varies with time, the pressure magnitude will be a function
of distance, x, as well as of time, t. From Newton’s second law of motion, the

acceleration of a material element resulting from an applied pressure is:

3p Bu

0:: —pE‘ (3.10)

The negative sign in Eq. 3.10 results from a net acceleration to the right for a

30

negative spatial pressure gradient. In Eq. 3.10 p is the medium density.
Equations 3.9 and Eq. 3.10 are the coupled pressure particle velocity. By decou-

pling these equations, it is possible to obtain the acoustic plane wave equations

32p k 32p
and,
6221 k 3211
__ = __ 3.12
at2 p 62:2 ( )
The general solution for pressure and particle velocity are
p = p(0)ej(“’t”K”) (3.13)
n = u(0)ej(“’t‘K”) (3.14)

where K is the wave propagation constant given by

K = tum (3.15)

The wave number K is in general a complex quantity. It consists of the phase constant

5 and the attenuation constant a,
K = fl — ja (3.16)

From Eq. 3.10, 3.13, and 3.14 , the relationship between pressure and particle

velocity is

_ 92
p— K“ (3.17)

The characteristic acoustic impedance is defined as the ratio of the pressure to particle

31

speed. From Eq. 3.15, the acoustic impedance Z can be expressed as

= E (3.18)

Z K

QI'B

For a lossy medium, the acoustic impedance is a complex quantity. For a lossless
medium the attenuation constant a is zero, and the wave number reduces to K = B,
which gives a real acoustic impedance.

Under the lossless assumption, the phase velocity 12,, is

(3.19)

Therefore, the acoustic impedance for a lossless medium can be written as

Z = pvp (3.20)

3.3 Transmission and Reflection Coefficients

When an acoustic plane wave is incident to a boundary between two different
media, it will be partially reflected. The ratio of the characteristic impedances of the
two media determines the magnitude of the reflection coefficient and transmission
coefficient. Assume there is an acoustic plane wave propagating from medium 1 to
medium 2, as shown in Figure 3.2. Snell’s law states that the angles of incidence
and reflection are equal when the wavelength of the wave is small compared to the

dimensions of the reflector. That is,

32

From the relationship between 0; and 0,, the angle of transmission, can be obtained

 

from:
sin 0.- u;
= — .22
sin 9, U2 (3 )
1 r
9 i 9 r
Medium 1

 

 

Figure 3.2. Transmission and reflection at an interface

In the equilibrium state, the pressure on both sides of the boundary remains the

same in order to maintain a stationary boundary. This gives

p; + pr = p: (3.23)

Furthermore, the normal component of the particle velocity must be the same on

33

both sides of the boundary or else the two media will not remain in contact. Thus,

21,-cos 0.- — 11, cos 0,. = u, cos 6t

 

 

(3.24)
From Eqs. 3.17 and 3.24, we have
piKl cos 0.- er1c050, _ ptKg cos 91 (3 95)
P1 P1 P2 ..
Solving Eqs. 3.23 and 3.25, the pressure ratios are
Pr _ £1 cos 0; —— £21 cos 0; (3 26)
p.- — %c030;+I—;1c080¢ °
and
p, 275—? cos 0.- (3 27)
__ - K , 1? °
p. 71 cos 0, + 721 cos 0‘)
For normal incidence 0.- = 9: = 0, the pressure amplitude ratios reduce to
£1 _ £2
Reflection coefficient, R = K'- = p‘ ”2 (3.28)
Pi 5* + 51
121 pa
2.5.1.
Transmission coefficient, T = fl = p‘ (3.29)
Pi 5* + I—(Z
p1 p2

Using the acoustic impedance definition Z = wp/ K , the reflection coefficient and

transmission coefficient can be reduced to:

R_a-a

_ —— 3.30
Z2 + Z1 ( )

T 2Z2

=—— .1
Zz+Z1 (33)

Up to this point, the reflection and transmission coefficients were expressed in the

34

time domain. Using the Fourier Transform the reflection and transmission coefficients
can be expressed in the frequency domain as well.

The Fourier Transform of p(:z:, t) is defined as
P(:1:,w) = [00 p(x,t)e—j“’tdt (3.32)

Because the reflection coefficient R and the transmission coefficient T are independent

of time, the R and T can be written as

_P,(0,w)_ $11—12 _zz—z1

R______7(____ 3.33

H(0,w) %L+721 Z2+Zl ( )
ygl

Tzwzmzfl— (3.34)

3.4 Broadband Signal Probing Techniques

Reflection coefficients and transmission coefficients are crucial factors in deter-
mining wave propagation in media with different material properties. When the
incident wave is a narrow band signal, the reflection coefficient and the transmission
coefficient are usually constant. This is not true for a broadband signal. The ba-
sic principle of using video pulses is to exploit the broadband signal characteristics
to detect media variations. The broadband signal from the transducer is dispersed
by the frequency dependent reflection and transmission at the interfaces of different
acoustic impedances of different media. The energy of the reflected wave or of the
transmitted wave depends on the frequency dependent reflection coefficient or the
frequency dependent transmission coefficient.

A necessary requirement for using video pulses is that the waveform of the inci-

dent wave in both time and frequency must be well known. When the video pulse

35

impacts on the interface between the first and second media, the reflected wave con-
tains all frequency components, each modified by its own frequency dependent reflec-
tion coefficient. The reflected signal spectrum contains the information necessary to
characterize the unknown media.

The Fourier spectrum of the incident pressure wave is given by:
P;(a:,w) = / p,(a:, t)e’j“"dt (3.35)

If a: is set to zero at the interface of two media, the pressure in the frequency

domain at the interface is:

R-(OM) = }_{140,10} (336)

where f {} = ff; ~e‘j‘”‘dt is the Fourier transformation.

In section 3.3, it was shown that the reflection coefficient is independent of time,
so that the Fourier Transform of the reflected wave at the interface in the frequency
domain is given by:

Pr(0,w) = P,(0,w)R (3.37)

and the reflected pressure wave at the interface in the time domain is
pr(0,t) = .74 {P,(0,w)} = f‘1{P.-(0,w)R} (3.38)

where .7‘1{} = 31; [3°00 -ej“"‘dw is the inverse Fourier transform and R is the time
independent reflection coefficient.

A detected signal a distance d from the interface has the form:

Md, t) = Pr(0, t)6’jK’d (339)

36

The signal to be detected is thus given by:

p,(a:,t) = Real {7‘1 {P;(0,w)R} e’jK’d}

= Real {f'l {f{p;(0, 1)} R} e’jK‘d} (3.40)

The real part of the pressure is taken since only this signal can be detected ex-

perimentally.

CHAPTER 4

DETERMINATION OF
VELOCITY AND
ATTENUATION BY
BROADBAND SIGNALS

In this chapter the acoustic plane wave equations including attenuation is derived
in 4.1. Based on the wave equations, the measurement of velocity-density prod-
uct and attenuation-density ratio by broadband signals will be discussed from two
different conditions. The first assumes the attenuation coefficient is proportional to
the frequency squared. The theoretical derivation is shown in section 4.2. In order
to generalize the frequency dependent attenuation, the attenuation coefficient is as-
sumed proportional to the l-th power of frequency, where I is a real number between

1 and 2. The details of this condition are discussed in section 4.3

37

w
4.1 Damping and Attenuation

As discussed in section 3.2, the acoustic wave equation Eqs. 3.11 and 3.12 was
derived for the condition that the wave propagating through the medium suffers no
energy loss, that is, the acoustic wave propagates in the medium without attenuation.
Ideal materials of this kind do not exist, although weakly damped materials are often
approximated as ideal. Elastic damping usually depends on temperature, frequency,
and the type of vibration. At room temperature, acoustic losses in many materials
may be adequately described by adding a viscous damping term. In this section, the
damping force is used to obtain the attenuation constant a and phase constant 6.

In an ideal lossless medium, Hooke’s law describes the linear relationship between
pressure and strain:

pideal = —CS (41)

where c is a proportional constant, called elastic stiffness constant, while S is strain,

as defined in Eq. 3.5. The damping force in material can be put in the form:

65
pdamping = —n—5t_ (42)
where 17 is the material viscosity.
The total pressure becomes
05
P = pideal + pdamping = —(kS + 77—5?) (4.3)

For a one-dimensional problem, the strain is

5-85

_ 37: (4.4)

39

From Eqs. 3.7 and 4.4 we get

an 326 85

5; = 0331 = 791‘ (4‘5)

Taking a partial derivative with respect to time in Eq.4.3, and using Eqs. 4.4 and

4.5 we obtain:

 

0p _ Bu 02a
5 — _ (Caz: + "6161) (4'6)

From Eq. 3.10, the spatial derivative of pressure and the time derivative of particle

velocity are related by
3p Bu

Following the same procedure used in section 3.2 to decouple Eqs. 4.6 and 4.7, the

acoustic plane wave equations including attenuation become:

3219 c8219 17 83p

375 = gag; zaxzat (4'8)
and,
8211 c 0221 17 8323
372‘ - 53.72 + ;—axzat (4-9)
The general solutions for pressure and velocity are
p = p(0)ej(“’"K’) (4.10)
u = u(0)ej(“’“K‘-') (4.11)

where

K = ,6 -—ja (4.12)

40

Using these solutions, the following equations are obtained:

pwzp = 6K2}? + jnszp (4-13)

pa:2 = cK2 + 3'1]sz (4.14)

With K 2 = 62 — 206 j — a2, and both a and )6 real, the following equations are

obtained by equating the real and imaginary parts of Eq. 4.14:

C(fl2 — 012) + 20600) 2 pt.)2 (4.15)

nw(B2 — 02) — ‘2ch = 0 (4.16)

Solving Eqs. 4.15 and 4.16 for real 0 and 5, gives:

02_ _%(-____ -_2_.2__n:€w+ c2)[‘/1 +(17w/c)2 — 1] (4.17)

and
,62 = %( _2—":€w: Cz)[\/1+(nw/c)2 +1] (4.18)

Finally,

- _1__:P___(/——2 "2
a_w{2(n2w2+cz)[ 1+(nw/c) —1]} (4.19)

.. 1m 2 ”2
fl—w{2(n2w2 +c,,)1(/1+(w/c) +11} (4.20)

The above equations give the general solution for the attenuation constant a and
Phase constant )8. Eqs. 4.19 and 4.20 indicate that both the attenuation constant a

and the phase constant ,6 are functions of frequency. We can express 0 and 17 in terms

41

ofaandflas:

 

2 2 _ 2
= “(293+ (12;) (4.21)
and,
_ 2pwafl
77 _ -——(,B2 + 02V (4.22)

For a given frequency, the attenuation constant a is proportional to the square root
of the medium density. As the viscosity constant increases, the attenuation constant
also increases. If the elasticity stiffness constant increases, the attenuation constant
decreases.

The acoustic impedance of material is given by [70, 72]:

Z = wp/K (4.23)

where K = B — ja, is the complex wave propagation constant. As a result, for
lossy media, the acoustic impedance is a complex quantity as well. At the boundary

between media 2' and i + 1, for normal incidence, the reflection coefficient becomes:
ZH-l - Z;
ZH-l + Z;

(Ki/Pi) " (Kin/P141)
(Ki/Pi)+(Ki+1/Pi+1)

 

= (fli/Pi - Bi+1/Pi+l) — flat/Pi - 0i+1/Pi+1) (4.24)

(fli/Pi 'l' fli+1/Pi+1)—j(ai/Pi+ ai+l/Pi+l)

42

wherej = \/—1. Thus,

(Bi/Pi — IB£+1/P£+l)2(ai/Pi — 0i+1/Pi+1)2

IRiI = (3i/Pi + fli+1/Pi+1)2(ai/Pi + ag+1/p.‘+1)2

 

(4.25)

From above equation, the reflection coefficient is a function of the attenuation
constant, the phase constant, and the material density. Therefore, the reflection

coefficient will contain information of a, ,6, and p.

4.2 Derivations under the Assumption that the
Attenuation Coefficient Is Proportional to
Frequency Squared

In this section measurements of the velocity-density product and attenuation-
density ratio of materials for isotropic, homogeneous, linear, lossy media are devel-
oped. For typical materials and frequencies of operation ( where (963)” << 1 ), the

attenuation and phase constants can be simplified from Eqs. 4.19 and 4.20 to:

a w "7137/7/31 4142/3211"
"21637511 — $4511

2
mo
-2c_ P/C

aw2 (4.26)

22

22

43

a z wMu-§(1’§)W
4271741- $123)?)

w p/c

22

22

5w (4.27)

where a and b are constants for given materials,

(4.23)

9
III
M
n I:
‘b
\
O

and

b E p/c (4.29)

From Eqs. 4.26 and 4.27, when the frequency of the acoustic wave is low (313)” << 1,
the attenuation constant is proportion to w”, and the phase constant is proportion to
w. For very low frequencies, the attenuation approaches zero. For higher frequencies
the attenuation increases rapidly. The longitudinal acoustic wave velocity in isotropic

media becomes
w

fl

From the approximations of a and 3 in Eqs. 4.26 and 4.30, it is clear that

‘0: R3

(4.30)

o‘lr-d

the attenuation coefficient is proportional to frequency squared and the longitudinal

acoustic wave velocity is constant. Therefore, the following expression can be obtained

from 4.25:

|R{(w)|2 = [(bi/Pi) _ (bi+l/Pi+l )]2 + [(ag/pg) — (a,+1/p,-+1)]2w2

[Uh/m) + (b.,,/,.,,,)]2 + [(a./p.-) + (a.+./p.-.1)12w2 (4°31)

44

For various frequency components of the echo return, a set of simultaneous equations
in the form of Eq. 4.31 can be established. In order to demonstrate the procedure for

identifying unknown material, the following cases are considered:

Case 1: The ith medium is lossless (a; = 0)

For simplicity, consider water as the first medium, which is most often used in
acoustic imaging setups. Water is an excellent acoustic transmission medium with
practically no attenuation in a finite length for frequencies well into the GHz region.
Setting the water attenuation coefficient a.- = 0, Eq. 4.31 reduces to:

A + Btu2 — 1

. 2 =
IR.(w)| - A+Bw2+1 (4.32)

 

where A and B are frequency independent real coefficients:

 

 

= (bi/Pi)2 + (bin/P1402
A — 2l(bibe+1)/(Pipi+1)l (4°33)
___ (Cl.'-1_1/,0.'+1)2
B ‘ masts/(11.3.41)] (4'34)

Assume the measurement error of lit-(11.2).”2 is 61;, that is,

A+Bw§—1

R,- 2=
I (wkll A+Bw§+1 +61:

 

(4.35)

Equation 4.35 must now be solved for A and B by regression. Unfortunately, this
is a nonlinear model for which it is difficult to get a closed form solution. In order
to linearize the equation, assume the measurement error e is small, that is, 6% z 0.

Thus Eq. 4.35 can be rearranged as follows:

45

(A + Bwilll - |R1(wk)l2 + 0:12 = [1 - IRKMN2 - Chill + l1531'(wk)l2 + 6k] (4.35)

(A + 19¢v12.)[1-ll’i‘(cwc)lr"]2 + 264(4 + 1969130--|1ris(wk)|2]+(A + Bwiki

= [1 - IR.-(w1.)l2][1+ 1190001214” 64145314604)!2 — 61‘: (437)

For 6% z 0, this expression can be reduced to:

l A +BwZ][1-|Ri(wk)|212-[1-|Ri(wk)|2][1+lRi(wk)|2]

z ZCkIR,‘(wk)|2 — 2€k(A + sz)[l - IR..'(wk)I2] (4.38)

The second term in the right hand side of Eq. 4.38 can be further simplified to:

6k (A + mel - |R1(wk)|’]

€k(1 - |R£(Wk)|2)(1 + IR.-(w1.)l’ - 6k)
1- ll‘ldcwcll2 + 6::

 

= €k(1+lRi(wk)12)_ 62(1+1R£(wk)12 _ 6k)

1- |R1(wk)|2 + 6k

22

0:0 + IR.-(wk)|2) (439)

Therefore, a linear model can be constructed as follows:

46

l A +Bw§][1—IR.-(wk)|2]2 — [1 - IR.-(w1.)l2][1+ |R.~(wk)|2]
$3 zéklR,‘(wk)l2 — 2€k[1+|Ri(wk)l2]

: _2€k (440)

Therefore, equation 4.36 reduces to the linearized mode:

[A + 15’willl-IR.'(wk)|2l2 m [1- le'(wk)|2][1 + |R.(w,.)|2] - 26,, (4-41)

Assuming the random error 25 is normally distribution with standard deviation
0, the solutions of A and B, with the minimum sum-of-the-squares error of n mea-

surement data points, A and B, are [50, 73]:

A
as (XTX)-1XTY (4.42)

where X is a n by 2 matrix and Y is a n by 1 matrix as follows:

p

(141213312)? (1— 1344231321»: '

X = , _ (4.43)

 

 

_ (1 — IR.-(wn)|2)” (143442312242: .

47

' (1— |&(w1)|2)(1+|E(w1)|2) ]

Y = , (4.44)

 

 

( (1 — |R4(wn)|2)(1+ 1112442312) .
and X T is the transpose of X.
In general, for a linear regression model, if the random error e is normally dis-

tributed, the estimated values A and B will be unbiased, normally distributed esti-

mators of A and B respectively [50]. The variances of A and B will be:
Var(A) = 0,402 (4.45)

and

Var(B) = c1302 (4.46)

where cA and CB are the diagonal elements in (XTX)‘1, or

CA CAB
z (XTX)-1 (4.47)

03.4 CB
We must either know a in Eqs. 4.45 and 4.46 or possess a good estimate of a so
that the variances of A and B can be estimated. In general, the value a is unknown

so that the estimated standard deviation, 5', should be calculated as follows [50, 73]:

YTY—[A f3]XTY

S n—2

 

(4.48)

In order to find the 90 ‘70 confidence interval of A, a hypothesis H0 is constructed

48

as follows:

Ho : A = A0 (4.49)
where A0 is a specified value of A. The quantity TA, which is defined as:

A

A—Ao
S\/c_A

will have a Student’s tdistribution with n — 2 degrees of freedom. Therefore, the 90%

TA =

 

(4.50)

confidence interval for A is [50, 73]:

A j: 10,05,,._2s,/c—A (4.51)

where t0_05,n_2 is a t distribution with (n — 2) degrees of freedom and 0.05 refers to
the 90% confidence and @452 is the estimated variance of A.
By the same method for finding the confidence interval of A, the 90% confidence

interval for B is:

B :t to,05,n_2S‘/CB (4.52)

where 6352 is the estimated variance of B.

With the properties of the ith medium known, the material properties in the
(i+1)th medium, such as (12.41 pg“) and (01,-4.1 /p.-+1), can be determined. The velocity-
density product in the second medium can be obtained, using Eqs. 4.30 and 4.33, as

follows:

”Mi
[.4 2: (L42 —- 1]
The attenuation-density ratio, using Eqs. 4.26, 4.33, and 4.34, becomes

3‘“ = (517) (\f2 [4 :1: (M? - 1' 1‘3) (.122 (4.54)
i+l i i

vi+1Pi+1 = Pi+lbi+1 = (4-53)

 

 

49

The :1: sign in Eqs. 4.53 and 4.54 can be determined by comparing the phases of
the incident signal and the reflected signal. If they are in phase, which implies that
the impedance of the second medium is greater than that of the first medium, the
“minus” sign is used. Otherwise, the “plus” sign should be chosen.

After A and B are evaluated, the velocity-density product as well as the

attenuation-density ratio can be obtained from Eqs. 4.53 and 4.54.

Case 2: The ith medium is lossy (a.- aé 0)

In general, the attenuation constant of the (i+1)th medium is not equal to zero.
For such cases the procedure for obtaining velocity-density products and attenuation-
density ratios becomes more complicated. If (bg/pg), (ag/pg), and IR;(w)| are known,

the following equation can be constructed from Eq. 4.31:

=C+(D—H)u.:"’-1

 

 

 

 

IR‘MP " C + (D + H)w2 +1 (4'55)

where C, D, and H are the frequency independent real coefficients:
DE‘“;(z.>:..:;7z;:.<::sa>2
H :— 2:: (4.58)

Since there are only two unknown quantities,(b.~+1/p.-+1) and (a;+1/p.-+1), in Eq.
4.55, we can choose C and D as independent variables and H as a dependent variable

of C and D, expressed as:

50

 

_ “ix/273W i «film/Pi)" — a?
_ b.-(C 1: \F—‘Cz — 1)

Assume the measurement error of le(wk)12 in Eq. 4.55 is 5)., and 5 has a normal

H

 

(4.59)

distribution with zero mean. Therefore,

C+(D—H)w,2,-1
C+(D+H)w§+1

 

11160-0012 = + 51: (4-60)

The sum of the square error, E,(C, D), is:

Ei(CaD)

ll
0,
arm:

 

II
M:

C+(D—H)w,f—1) (4.61)

0124“”le _ C + (D — 11).): +1

For the minimum sum of the square error, we can set the derivatives of E;(C, D)

with respect to C and D to zero, respectively. Therefore,

813,-(0, D)
50

2:? (10+ (D -411)w13+112) '

C+(D—H)<.u2—1 2

 

 

and

51

313,-(0, D)
DD

= [2(1e+ef“5)n+n)-

C+(D—H)w2—1 . 2
(C + (D _ H)“; +1 " 11210001 ) (4.63)

 

 

Since the model in Eq. 4.60 is nonlinear, the solution of C, D, and H cannot be
obtained directly. However, we can solve for C, D, and H by Newton’s methods [74].
From Eqs. 4.62 and 4.63, two functions, f1(C, D) and f2(C, D) can be defined

03,-(0, D)

fl(C,D) = ac

 

_ 1.; ([0 + (D —4H)w£ + 1121'

C+(D—H)w2—1 2

 

and

0E,(C, D)

6D
_ " 4903 ,
- .§(n+e_nwg+nzl

C+(D—H)w2—1 2
(0+ (D _ H)»; +1 _ [127034)] ) (4.65)

f2(C, D) =

 

 

The solutions of C and D can be obtained when both f1(C, D) and f2(C, D) are equal
to zero. This result can be found by iteration. Initial values Co and D0 are given and

new values CH1 and Dk+1 obtained from Ck and Up, as follows:

52

C C C,D
”1 = " —J(C,.,D,.)-l M " k) (4.66)

Dk+1 DI: f2(Ck, Dk)

where J (Ck, D1.) is the Jacobian matrix of (Ck, Dk), and defined as:

aflgcp) 8f1§C,D)
(4.67)

J C ,D =
( k k) 3f2(C'.Dl 3f2(C,Q)
6C 3

C=Ck , D=Dk

The algorithm to find solutions of C and D can be described as follows:

Step 1: set k = 0, assuming Co and D0 are given initial values;

 

ag \/2Dk(Ck:f:‘/ CE-IXM/ng—a?

 

Step 2: H = 51(Cix/CE—1)
Step 3: Ck“ = C" _ J(Ck,Dk)‘1 f1(Ck,Dk)
Dk+1 Dk f2(Ck,Dk)

Step 42 1f]f1(Ck+1,Dk+1)[+ [f2(Ck+1,Dk+1)[ > E and k < kmax
then 16 = 10 +1; go to Step 2;
else STOP.

In the above algorithm, the iteration will continue until the value of
I f1(Ck+1,D,‘ +1)I + I f2(Ck+1,Dk+1)| is less than a small quantity E, at which point a
good approximate solution is obtained. To avoid the solutions from local minimum
squared error, some other initial values should be tried to make sure the solutions
have a global minimum. If the iterative number, k, is equal to the maximum allow-

able number of iterations, km”, and the value of |f1(Ck+1,Dk+,)[ + [f2(Ck+1, Dk+1)l

53

is greater than 5, the iteration may diverge so that a good approximation is not

obtained. If an unacceptable solution occures, other initial values should be tried.
When the solutions of C and D with minimum squared error, C and D, have

been found, the velocity-density product and the attenuation-density quotient in the

(i + 1)th medium can be expressed in terms of the quantities of the ith layer by:

 

 

 

 

vi+1Pi+1 = if“
3+1
= . "”0: (4.68)
[C :1: (E12 — 1]
CYi+1 _ (11410122
Pi+1 Pi+1
A A A 1/2
2[C :1: V 02 — 1]D a,- 2 2
(mm)2 p.-

where C and D are the solutions with minimum sum of the square error for C and
D. The :I: sign in Eqs. 4.68 and 4.69 is determined by the same rule as for Eqs. 4.53
and 4.54.

It can be seen that the velocity-density product and attenuation-density quotient
of the second medium can be found if the parameters of the first medium are known.
To generalize the approach, it is possible to evaluate material properties of an n-

layered structure by repeating the procedure successively.

54

4.3 Derivations under the Assumption that the

Attenuation Coefficient Is Proportional to

the l-th Power of Frequency

For some media, the attenuation constants are not proportional to frequency

squared. In this section, a more general case will be discussed. We assume the

attenuation constant a is proportional to the l-th power of frequency, and )6 is still

proportional to frequency.

The following equations can be obtained from Eqs. 4.19 and 4.27:

From 4.25, we have

[3104012 =
[(bi/Pi) "' (bin/P4012 + 1(aiwl"1/Pi) — (0‘:‘+1¢1’l"‘”-1/P:‘+1)l2
[(bi/Pi) + (b1+1/p.-+1)]2 + [(aiw"‘1/ 101) + (01+1w"+“‘/ P141)?

1+ IR.-(c«2)|2 :
1 - ӣ101012
(bi/[’02 + (bin/PHI)2 + (aiwl‘—l/Pi)2 + (“1+1‘4’l‘“'1/Pi+1)2
21(51/ pa)(b.-+1/ PM) + (aiw""‘/ p,)(a,-+1w‘1+1-1/ Pi+1)l

 

(4.70)

(4.71)

(4.72)

(4.73)

If we assume the attenuation of the i-th layer is zero, that is, a,- = 0, and a; = 0,

55

then Eq. 4.73 can be simplified to:

1 + I1‘Z1'(w)l2 _ (be/101)2 +(bz'+1/p.-+1)2 + (a;+1w"+1‘1/p;+1)2

 

 

 

 

 

1- |R1(w)|’ ‘ 21((bibi+1)/(P1Pi+1)l (4'74)
1 + IR:'(¢1’)|2
”“1 = 1411-01))?
2 A + Bw'" (4.75)
where
= (bi/Pilz + (bi+1/Pi+1)2
A ‘ 21(bibi+1)/(Pipi+1)l (4°76)
.__ (0:41 /Pi+1 )2
B ‘ 21(bibi+1)/(P£Pi+1)l “'77)
m E 2(l,-+1 — 1) (4.78)
The total sum of square error (SSE) of 11 data points will be
(SSE) E En:[F(w;) — (A + B(D;")]2 (4.79)

i=1

The minimum values of (SSE) will occur when 21%? = 0 , i3? = 0 , and 2%)) = 0

6(SSE)
a4

= :(—2)[F(w1)—(A+Bw:")] (4.80)

 

56

 

0 0(SSE)
BB
= ;(—2w.’")[F(wi)-(A+Bwf")] (4.81)
_ 6(58E)
° — 77m—
= Z(-2B(lnw1)wm)[F(w.-)—(A+Bw:")] (4.82)

8:]

Eqs. 4.80,4.81, and 4.82 can be written as

:(A + sz") = :Ffiug) (4.83)
:(A + Beam? = E; F(w.-)w:" (4.84)
ERA + Bwf")w}" 111w.- = in: F(w,-)wf" In (D; (4.85)

i=l 5:1
This provides three unknowns (A,B,m) and three equations (Eqs. 4.83, 4.84, and
4.85). This nonlinear regression problem can be solved for A, B, and m by iteration.
First, the two linear parameters A and B can be expressed in terms of F(w,~),w.- , and

m. From Eq. 4.83,
?=1(F(w1) — Bwtm)

 

 

A = n (4.86)
Inserting Eq. 4.86 into Eq. 4.84, allows B to be expressed as
B : comm.) — Z?=1(F(wj)/n)] (4,7)

211:1 wAm[wAm_ ;-‘=1(w;"/n)]

57

From Eq. 4.85 we can define a new function f (m),

f(m) '_=_ i[F(w;) — (A + sz")]w:" ln (D; (4.88)

i=1

with the solution determined by f (m) = 0. Therefore, the problem reduces to finding
zeroes of f (m) when A and B are known. First, an initial value of mo is assumed,
a new value m1 is found by Newton’s method to make f (m1) approach zero. This

method is continued according to the recursion relationship:

mk+1 = m. — (f(mk)/ [431?] m...) (4.89)

where £3131 can be obtained from Eq. 4.88 as follows:

’ df(m)

dm = En:[F(w.-) — (A + 2Bw?)]wf"[ln w;]2 (4.90)

i=1
Therefore, the values of A), and B), can be obtained from Eqs. 4.86-4.90.

The algorithm for the nonlinear regression can be described as follows:

58

Step 1: set k = 0, assume an initial value for ma,

B = 211:1wrlefwil‘ZLJFle/nll

L1 w?“ lw?‘ {3:1 (5311/ ")1

Step 2:

 

A = wed-331"")
Step 3: mk+1 2 ml: " (f(mk)/ [%2 m=mk)

Step 4: if If(mk+1)| > 6 and k < km“,
then 10 = 11: +1; go to Step 2;
else STOP.

Following the discussion in section 4.2, the iteration will continue until the abso-
lute value of f (mk) is less than a small quantity 6, at which point a good approximate
solution is obtained. To avoid solutions which are from local minimum squared errors,
other initial values should be tried to make sure the solutions with global minimum
are found. If the iterative number, k, is equal to the maximum allowable number
of iterations, km“, and the absolute value of f (mk) is greater than 6, the iteration
does not converge so that a good approximation is not obtained. If an unacceptable

solution results, other initial values should be tried.

CHAPTER 5

EXPERIMENTAL PROCEDURE
AND RESULTS

In order to analyze the characteristics of an unknown medium by acoustic waves,
we can use either the transmission method or the reflection method. For the transmis-
sion method, at least two transducers, one for transmitting and the other for receiving,
must be used. In addition, the two transducers must be aligned on the same axis to
maximize the detected power. On the other hand, the reflection method only re-
quires one transducer. This transducer acts as both the transmitter and the receiver.
Since the reflection mode measurements are easier to implement, this configuration
is chosen for experimental confirmation.

Figure 5.1 shows a computer simulation for reflected signals, which are Gaussian
pulses modulated with 2.25 MHz, from the interface of water and three different
materials. The three materials are aluminum, plexiglass, and wood. There are several
facts we can find from the simulation results; first, the amplitude of the reflected signal
depends on the impedance of materials. The phase of each signal is also related to the
impedance difference at the reflected interface. For example, the phases of reflected
signals from aluminum and plexiglass are the same as for the incident signals. This is

because the impedances of aluminum and plexiglass are greater than the impedance

59

60

of water. On the other hand, the impedance of wood is less than the impedance of
water, therefore the phase of the reflected signal from the interface of water and wood
will be different than the phase of the incident signal. Secondly, the reflected signal
may not be as symmetrical as the incident signal because of the dispersive properties
of the material. However, the shape of the reflected signal contains the characteristics

of the materials and provides information for material feature identification.

5.1 Experimental Procedure

An essential part of the measurement procedure is to record the detailed time and
frequency characteristics of the incident broadband signal, since the reflected wave
contains all frequency components modified by their individual frequency dependent
reflection coefficients. Experimentally, the incident signal can be calibrated by using
a water / air interface which gives almost complete reflection. The effects of diffraction
are reduced by using the water/ air interface as a reference.

Replacing the water/air interface with a water/object interface will allow inves-
tigation of the detailed reflection spectral characteristics for feature identification.

The block diagram of the experimental arrangement for reflected mode are shown in

Figure 5.2

The incident pulse is obtained from a 2.25 MHz transducer (Panametrics V306)
and a Panametrics 5050 pulser, which provides an approximately Gaussian shaped
pulse of .5 psec pulse width. The reflected signal is sampled by an A / D converter with
a sampling rate of 40 MHz and 8-bit resolution. The test materials used are aluminum,
composite material, plexiglass, and de-gassed wood. For each test material, 256

signals are transmitted and the average of those reflected signals is recorded with

61

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2(1)] ' - 1500 r
l
1
2:21 lwr «
5 1 . In”
Incident g 0&4] J\———— . .22:
' 3 i i 500-
§ .100» .
.2“) 1 0 1 1
4 5 6 7 8 0 1
2w I 15m 1m V
3' mi ( i
. < I
Alummum 1; 0 I . 3.“
5 i 2
11* .100» .
-200 ‘
4 5 6 7 8
200 e e
a. 100+ ‘
. 5 - /\ . ..
Plemglass 2 0' ‘ ‘2'
3- 500»
i 4001- '1
-200 * ‘
4 5 8 O0
zoor . 1 1500 e 1 f
a 100+ “
Wood g , 2
§ 400+ <
-2(X) * ‘
4 5 6 7 8
Timemsec) Frequency(MI-lz)

Figure 5.1. Computer simulation of the reflected Gaussian pulse from four different
materials

62

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Water Tank
T i R
V
Panametrics
IBM PC Pulsar/Receiver
Compatible
Microcomputer
Signal Processing
and
Display '
A/D Convertor
r—— 40 MHz

 

 

 

 

 

Figure 5.2. The block diagram of the experimental arrangement for reflected mode

63

1024 sampling points. From the average of 256 signals, the signal-to-noise ratio was
improved by a factor of m, or 16. Figure 5.3 shows the procedure of a computer
program which will control the sampling unit and take the average of the reflected
signals. This program is written in C language. The starting sampling points, the
number of sampling points, the trigger level, and the number of signals to be averaged,
are set at program initiation. The averaged signal is displayed on a monitor screen

and the value of each sampling point is stored in a computer memory data file.

Figure 5.4 shows the procedure for data processing. Initially, the incident signal
and the reflected signal are recorded in the time domain. After acquiring these signals,
the spectra of the incident signal and of the reflected signal can be obtained by a
1024-point Fast Fourier Transform done in software. From the spectra of the incident
and reflected signals, the reflection coefficient , 12(0)), of each frequency in the 3 dB
bandwidth can be obtained. Using a least-squares regression procedure, the velocity-

density product and attenuation-density ratio for each material can be determined.

5.2 Experimental Results

Using the experimental setup mentioned in the previous section, the incident signal
and reflected signal from the test materials can be measured. The signals in the time
domain and their frequency spectra are shown in Figures 5.5 and 5.6. The measured
incident signal in the time domain is close to the Gaussian pulse which is used in
computer simulation. As we predicted, the phase of the reflected signals from the
interfaces of water and aluminum, composite material, and plexiglass are the same
as for the incident signal, while the phase of the reflected signal from the interface

of water and wood is different. The spectrum of the incident signal although not a

64

 

set parameters for sampling

N: no. of signals to be averaged
L: no. of sampllng points
k = 0;

 

 

 

I'IO

 

 

 

 

yes
Trigger pulser
k = k + 1

 

1

Take samples, S(n)

 

I Sum(n)=Sum(n)+S(n)
I I1 = 1 .. L:

 

 

 

 

Find averaged signal

 

 

 

Figure 5.3. The procedure of the sampling program

65

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Reflected Incident
Signal Signal
Find Spectrum Find Spectrum
of reflected signal of incident signal
Find reflection

coefficient: R(w)

1

Regression

1

Find
attenuation const.
and velocity.

 

 

 

 

 

 

 

 

 

 

 

Figure 5.4. The procedure for data processing

66

perfect Gaussian is sufficiently close that the approximation of a Gaussian shaped

spectrum is quite acceptable.

From the spectra of the reflected signal and incident signal, the reflection coef-
ficients in the 3-dB bandwidth can be obtained. Figure 5.7 shows the reflection
coefficients, R(w), for the interfaces between water and four test materials. From
Figure 5.7 we can see that all the reflection coefficients of the test materials have the
tendency to increase as frequency increases. This matches the theoretical derivation

shown in Eq. 4.32.

From the reflection coefficients in the 3-dB bandwidth, the values of A and B
in Eq. 4.35 can be obtained for each material. The results of regression, with 90%
prediction confidence intervals, for the above materials can be found from Eqs. 4.51
and 4.52 as shown in Figure 5.8. In the experiments, the 3-dB bandwidth is about 1
MHz and the variation of the reflection coefficients of the test materials is about 0.03
from 1.7 MHz to 2.6 MHz. From another point of view, the reflection coefficient of
aluminum is the largest and the reflection coefficient of plexiglass is the smallest for

the test materials studied.

Based on the results of Figure 5.8, the velocity and attenuation of test materials
can be obtained. The standard deviation of the velocity can also be obtained from
Eqs. 4.51 and 4.52. Table 5.1 summarizes the experimental results for a number
of different materials using this approach. The values of R2 statistics of regression
are also calculated according to the constraints of the regression. By comparing the
experimental results with the published values, it is observed that the velocities of
aluminum and plexiglass are very close to the published values in [71, 4]. The variation

between the experimental mean values and the published values is less than 3%.

67

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

200
2‘- 1001' \ i
c 1
. < A \
InCtdcnt 1; 01- .4.“! 5. / \ mm.
:3: 'v ‘1/
i. 400- -1
-2(X)'
2 3 4 5 6
2(1) j
l
3 HIV 4
5 ,‘. .
Aluminum g 01- WIN I ‘ )l\\/~——\M
.- ' l .'
3 1.2] 1.]
a 'l“)" 4
-2(X)
2 3 4 5 6
200
3 1(1)1 .
E
_ < z\
Composxte 1,1 0+ ~/\ / 1, [\IA
3 ‘J v
2 4001 4
-2(X)
2 3 4 5 6
200
a 100+ .
. 5 ,1
Plcx1glass g 0t ~“V -f\—
3
$2 400-
-2(X)
2 3 4 5 6
200
a 100+ i
i
3
a? 4001' «
-200 *
2 3 4 5 6

Timemsec)

Figure 5.5. The incident signal and reflected signals from test materials

68

 

‘8‘

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

F1311
‘ 'F? i3
g 1000 - 11.1:-
wt
. U .1]? :1 ‘1
lnc1dcm E '1 7911‘ [i j.
5 500 ~ $3,111?
1 "”[Milll' ['Al':
11111191 » ..:Elilil["[ll E 1.1
:."li1Ԥ. . 111': 1.1111111.
0‘ - l l 'l 1. ' 11
0 l 2
1500 1 e
("1 l
'2’ 1°00 ' .11]
Aluminum .1): [1 [1;]
2 ll
3: 5°° * 111:2”)?
111‘ 'l’ 1.1
11- 1 1 1
0 .Ti: [aililll (1'1'111 I111
0 l 2
1500 f
g 1000 1- a
c ' ‘>’
cmposxte 3
U
a:
23‘ 1000 - «
2
Plexiglass ‘-'-' .11,
3 500 1' (:‘ili 11111111., ‘
:6
a
2
3
Wood ‘7;-
‘15
a:

 

 

I l i 1' 2‘ 3
Frequencyfllrfl-Iz)

Figure 5.6. The spectra of the incident signal and reflected signals from test materials

 

69

 

 

 

 

 

 

 

 

lit I! ill II .I iA
Iii ill- i Iii-A
i111 -11 - -1 1.1
iii 11 111. I. . 1.;
iii-11 i 11. i 1...;
a 1 - l -1 3 r

iiilli -1 1 .1...

11 1-l- 1- 1 11- 1A

1 1 - 1 - -1 1-1..

 

 

11H- ”111-HI.-.” 11H.#

2.5

 

 

 

 

'
' 2‘5

I'll 1-

[.1
6. 4
0

.tuoU 43:3—

Plexiglass

r
-- -r «J

0.2 1-

 

 

 

 

 

 

 

 

 

 

 

 

ii ii I iiil
iii- 11101! ii‘
I- .1 1 III-i I I;
i Iii-i I I5;

 

 

 

2.5

Frequency (MHZ)

Figure 5.7. The reflection coefficients of test materials within their 3-dB bandwidth

Reflection Coefficient

70

 

 

 

 

 

 

0-9
Aluminum

0.8 -

0.7 a ----------- ________________________________________
-------------------- Wood

0.6 : ------------------------------------- fl ______________________________
------—- -------------------------------- Composite

0.5 -

0'4 f::::‘:::2::-'----"“"‘ ---------------------------- Plexiglass

 

1.7 1.8 1.9 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7
Frequency (MHz)

Figure 5.8. The results of regression with 90% prediction confidence interval

71

 

 

 

Material Aluminum Composite Plexiglass Wood
Density 2695 1609 1182 525
(kg/m3)

Velocity-density product 1.695 x 107 4.884 x 106 3.075 x 106 3.687 x 105
(kg - sec"1 - m'z)

Attenuation-density ratio 6.575 4.238 3.866 25.92
(7712/ kg)
Velocity, Mean :1: S.D. 6287 j: 117 3036 :l: 44 2602 i 29 702 :1: 13
(m/ sec)
Velocity (published) 6350 [71] — 2670 [4] —
(m/ sec)
R2 Statistics 0.63 0.60 0.82 0.78

of regression

 

 

 

 

Table 5.1. Experimental results of measuring velocity and attenuation by broadband
signals

 

C)"

17

1..

72
5.3 Two Dimensional Imaging by Using Velocity-
Density Product

To improve on the measurement results of a one-dimensional experimental setup,
a scanning system was built to obtain two-dimensional images. Figure 5.9 shows the
block diagram for this scanning system. The system is controlled by an IBM PC / XT
compatible computer. Two stepper motors are used to locate the transducer in the
proper position. A Panametrics 5050 pulser/ receiver is used to excite the transducer
and receive reflected signals. Again, the received signal is sampled at 40 MHz and
the samples are stored in the computer. By using the technique of broadband signals
the velocity-density product of materials can be obtained and the two dimensional
image can be constructed by scanning. A computer program written in C language
is developed to control the scanning motor, get sampling data, and display color
images. The program list of the main program is in Appendix A. This program
can get data from the scanning system which is described in Figure 5.9, get data
from a previous recorded data file, or display a constructed image file. This program
will automatically send out control signals to operate the scanning system according
to the preset values. The data obtained from each (X ,Y) location in the scanning
system can be stored in a file for further analysis. Once the desired feature, such
as velocity or attenuation, has been obtained, a two dimensional color image can be
displayed on a monitor screen. The colors of images can be selected from 16 preset

colors to get the best display.

Figure 5.10 shows a reference sample to demonstrate the performance of the two
dimensional scanning system. The sample is made from two concentric cylinders of
Teflon and aluminum, respectively. This sample is placed inside a water tank and

the ultrasonic transducer scans the top of the sample. Therefore, this represents two

73

 

 

 

 

Panametrics AID converter Data processing
5050
Pulsar/Receiver 40 kHz and Display

 

 

 

 

 

   

 

 

 

landlo-
7 rs'm‘m W” Y IBM 9cm
compatible
miaocomputor

 

 

 

 

 

 

 

 

 

 

Water Tank

Figure 5.9. The block diagram of the scanning system

nan 'i

.y..~r

iii;

74

contiguous 2-layered sample. The first layer is water and the second layer is either
Teflon or aluminum. Using the techniques outlined above, an image of the velocity-
density product, as shown in Figure 5.11, is obtained. The blue area is the reflected
signal from Teflon and the red area is the reflected signal from aluminum. It is clear

that two different materials are present in the image.

A 3-layered model, as shown in Figure 5.12, consisting of a layer of plexiglass
placed over the first sample was also investigated. The first layer is water, the second
layer is plexiglass, and the third layer is Teflon or aluminum. This model is used
to demonstrate that the technique can be applied to retrieve image of deep-lying
targets. After evaluating the properties of the first layer, the properties of the second
layer can be obtained by using Eqs. 4.53, 4.54, and 4.68. Figure 5.13 shows the
image obtained from the velocity-density product for material in the second layer.
This image is not as clear as the image obtained from the 2-layer model because the
calculation error will accumulate from layer to layer. Although the error is larger than
the error in the previous model, the areas reflected from Teflon and aluminum still
can be distinguished. This result demonstrates that multilayer material identification

is feasible using the techniques described in this research.

75

 

Figure 5.10. A 2-layered 2-dimensional sample model

76

 

Figure 5.11. Imaging of the 2-layered 2-dimensional sample model using velocity-
density product data

77

Plexiglass

 

 

~ «A e
‘mlil‘ll‘I "

 

Figure 5.12. A 3-layered 2-dimensional sample model

78

1

I

1.
0.
O.
0.
0.
0.
0,
0.

Ft I. - rum-«:2. In.)

 

Figure 5.13. Imaging of the 3-layered 2-dimensional sample model using velocity-
density product data

CHAPTER 6

ULTRASONIC IMAGES AND
TISSUE CHARACTERIZATION
BY HIERARCHICAL
CLUSTERING TECHNIQUES

In the previous chapter the information of velocity and attenuation of test mate-
rials has been determined by broadband signals. In addition to the one-dimensional
processing, two-dimensional images using velocity-density product have been con-
structed. In order to investigate the properties of materials (or tissue), classification
techniques are used for automation. The objective of this chapter is to classify im-
ages constructed by ultrasonic broadband signals using hierarchical clustering meth-
ods [75, 76]. Section 6.1 introduces the hierarchical clustering techniques. In section
6.2 the image from the velocity-density product obtained from a 2-layered model
is used for classification. The clustering results for this known model provides the
confidence required to classify images using hierarchical clustering methods. Finally,
the clustering technique is applied to a section of human brain with a hemorrhaged

tumor. in section 6.3.

79

80
6.1 Hierarchical Clustering

A hierarchical clustering is a sequence of partitions in which each partition is
nested into the next partition in the sequence [76]. This method provides techniques
to Show the relationship among test patterns. In this research, the test pattern will
be obtained from the value of each pixel in an ultrasonic image. By this technique
the patterns with similar properties will be classified as in the same group. Therefore,
the patterns with different properties will be distinguished and material identification
or tissue characterization can be done. Before using a hierarchical clustering method

to classify images some matrices and preprocessing should be performed.

Pattern Matrix

The first step in clustering data is to construct the pattern matrix 1’. At each
(X ,Y) location of the test sample several features, such as total energy, peak fre-
quency, or velocity-density product, can be measured. Once the information of each
feature is obtained a pattern matrix can be construct for further analysis. A pattern
matrix is a matrix which contains the information of each feature in a column and

all information obtained from a certain location in one row.

Normalization

Since the scale and units of each feature are different, all features must be nor-
malized before processing. The normalization is based on the mean and variance of
each feature. If the number of patterns in the analysis is n, the j -th feature value

for the i—th pattern denoted by $3}, the data set can be normalized by:

1 l
.r
$18.:

‘ I
*1
‘r
‘11
...
q
(a,

81

 

where
21:1 CC:
m,- = —n——J (6.2)
and
2 ?=1(x:‘j _ TnJ)2
. = 6.3
s. n ( )

$.‘j is the normalized value of :53}, m,- is the sample mean of the j—th feature, and s}
is the sample variance of the j-th feature. The normalized pattern matrix is denoted

as ’PN = [mg].

Eigenvector Projection

Since the features may be correlated, the normalized pattern matrix should be
further transformed to uncorrelated features by eigenvector projection. This projec-
tion is performed by the eigenvectors of the covariance matrix C, which is defined

as:

705%

n

C:

 

(6.4)

where P}; is the transpose matrix of 'PN.
Finding the eigenvectors of the covariance matrix C, the pattern matrix with

uncorrelated features, 791;, can be found.

82

19;; = 191v ' (6.5)

 

 

IV}.

where V,- is the ith eigenvector of the covariance matrix, d is the total number of

features.

Proximity Matrix

A proximity matrix should be constructed by using the pattern matrix with un-
correlated features. The proximity matrix contains information on the Euclidean
distance of each pair of patterns. The (2', j)—entry in the proximity matrix provides
the distance between pattern 2' and pattern j. This results in a symmetrical proximity
matrix. The size of a proximity matrix of an image with n pixels is n by n, which is
usually very large. For example, a 50 by 50 image will contain 2601 pixels, and the

proximity matrix will be 2601 by 2601.

Hierarchical Clustering Methods

The most commonly referenced hierarchical clustering methods are [76]:

1. Single link method.
2. Complete link method.

3. Ward’s method.

4. Unweighted Pair Group Method using Arithmetic averages (UPGMA).

83

For the single link and the complete link methods, the clustering results depend on
the proximities through their rank order only. Ward’s method provides a minimum
square error and is therefore also called the minimum variance method. UPGMA
(Unweighted Pair Group Method using Arithmetic averages) computes the distance
between two groups by arithmetic averages. The effectiveness of each of the above
algorithms is discussed below:

Assume the proximity matrix is D = [d(i, j )], and there are n patterns. L(m) is
the distance of clusters after the mth merging operation. The algorithms each follow

the initial five step sequence:

84

Step 1: Set level L(0) = 0 and the sequence number m = 0.

Step 2: Find the cluster with the minimum distance between each pair.
Assume the minimum distance pair is cluster 1' and cluster 3.

The distance between cluster 7' and cluster 3 is

d[(1‘), (8)] = min{d[(i), (J')l}

Step 3: Increment the sequence number, m = m + 1.

Set the mth level L(m) = d[(r), (3)].

Step 4: Calculate the distance between the old cluster k and
the new cluster (1‘, s) , d[(k), (r, 3)], in

the proximity matrix, ’D.

Step 5: After all objects are in one cluster then STOP,

else go to Step 2.

The distance between the old cluster k and the new cluster (1', s) , d[(k), (r, 3)], is

calculated according to the different clustering methods as follows:

for the single link method:

d[(k), (133)] = min{d[(k), 0)], d[(k), (3)1}, (6-6)

for the complete link method:

(“(1“), (r9 3)] = max{d[(k)a (7')], duff), (3)”, (6'7)

85
for the Ward’s method:

(n. + nl. + m.) [(n' + n")d[(k)’ (r)]+

("a + nk)d[(k), (3)] - nkdlUC), (7‘)] 1, (6-8)

 

(40¢). (r, 8)l

and for the UPGMA:

n,d[(k), (r)] + n.d[(k). (8)]. (6.9)

duh), (r, 3)] = (n. + n.)

 

Where n,, n,, and n], are the number of patterns in cluster 1', s, and k, respectively.

Cluster Validation

After getting the results of clustering, a quantitative and objective method is
usually necessary to evaluate the results of clustering. The discussion of cluster
validation provides an index of the accuracy of clustering. Usually, the validity of a

clustering result can be expressed by the following three indices:

1. External index: a measurement by matching a clustering result to a priori

information.

2. Internal index: a measurement by using the original data only, no a priori

information is needed.

3- Relative index: comparison of two or more clustering results to decide which

one is more appropriate for the data.

Since there is no a priori information on the human brain sample in our exper-
lments, only the internal index and the relative index will be discussed. The most

COlrrlrnon internal index is called CPCC, or cophenetic correlation coefficient. Be-

86

fore defining the CPCC, the cophenetic proximity matrix must be constructed. The
cophenetic proximity matrix is a matrix which contains the clustering results from a
certain method. The (i, j ) entry of the cophenetic proximity matrix is the level of

object i and object j merged into the same group.

The definition of the CPCC is:

 

opcc = (W >221 d(i.j)d.(i.j) — (mama)
[(l/M) 2:21 0120.1) - mall/ml/M) ?=1Z;-‘=1d?:(i,j)- mall/2
(6.10)
where
_ n(n — 1)

M — .2 (6.11)

m = (fiéidus) (6.12)

me = (fiiédcw) (6.13)

d(i, j) is the given proximity between objects 2' and j, and dc(i, j ) is the cophenetic
proximity.

The value of CPCC is between -1 and 1. When CPCC is closer to 1, it means
the clustering result fits the original data better. If CPCC is equal to 1, it means the
clustering is perfectly matched to the original data, whereas, if CPCC is equal to -1,
the clustering result indicates the data is totally mismatched. Once the internal index
for each method is obtained, the values can be used as a relative index to compare

the clustering results between the different methods.

87
6.2 Material Identification by Ultrasonic Images

Material identification is an important topic in ultrasound. By using broadband
Signals, more information can be obtained in one transmitted signal than by con—
ventional detection techniques. In previous chapters, it has been proved that the
velocity-density product and attenuation-density ratio can be obtained from broad-
band signals. By using ultrasonic scanning systems two-dimensional images can be
constructed. In this section, the images obtained from a known sample will be used
to identify test materials. The sample used in the research is the same as the sample
used in section 5.3, which is shown in Figure 5.10. Using the image of velocity—density
product as shown in Figure 5.11, material identification can be done by hierarchical
clustering methods. In the image of Figure 5.10, only the velocity-density product
feature of the material is shown. Because the size of the image is 51 by 51 and
there are only two different materials in this model, we can consider this problem as
classifying 2601 patterns into two groups. A pattern matrix of 2601 by 1 and a prox-
imity matrix of 2601 by 2601 are constructed before clustering. Since there is only
one feature, it is not necessary to do normalizations and eigenvector transformations.
By the single link, complete link, Ward’s, and UPGMA methods, which have been
mentioned in section 6.1, all the patterns will be classified into proper groups.

Once the results of clustering have been obtained, color images can be constructed
and show the areas of different materials. The results of the single link, complete link,
Ward’s, and UPGMA methods are shown in Figures 6.1, 6.2, 6.3, and 6.4, respectively.
The results obtained from the complete link, Ward’s, and UPGMA methods are fairly
similar except for the boundary between the Teflon and the aluminum. All three
methods provide clear and good results to show two different materials from the
images. The single link method classifies the area reflected from aluminum by two

subgroups while the area reflected from Teflon is still in one cluster. Therefore, a

88

three-cluster image is shown instead of a two-cluster image.

In order to test the validation of clustering by different methods the values of
CPCC are calculated according Eq. 6.10. The values of CPCC by different clustering

methods for the 2-layered model are shown in Table 6.1 and Figure 6.5.

 

Method Value

 

Single Link 0.50
Complete Link 0.50
Ward’s method 0.52

UPGMA 0.81

 

 

 

Table 6.1. The values of CPCC by different clustering methods for the 2-layered

model

From Table 6.1 the values of CPCC by the single link, complete link, and Ward’s
method are very close, which imply that the results of these three methods provide
about the same degree of accuracy. According to the values of CPCC, the best
clustering result should be the image constructed by UPGMA because its CPCC is
the highest. On the other hand, comparing the original model, which is shown in
Figure 5.10, and the clustering results, we also find that all the images constructed

by the above clustering methods provide good classification.

89

 

Figure 6.1. Clustering result of a simple 2-layered model by single link clustering
method with 3 clusters

 

Figure 6.2. Clustering result of a simple 2-layered model by complete link clustering
method with 2 clusters

90

 

Figure 6.3. Clustering result of a simple 2-layered model by Ward’s clustering method
with 2 clusters

’{i‘hfiJ—‘I :Efl‘fifi
.; ‘ 1 15‘ “ '

 

Figure 6.4. Clustering result of a simple 2-layered model by UPGMA clustering
method with 2 clusters

91

Internal Index
CPCC

 

 

 

OO‘UO

 

 

SL CL Ward’s UPGMA

Figure 6.5. Comparison of CPCC by the different methods for the 2-layered model

92
6.3 Tissue Characterization by Ultrasonic Broad-
band Signals

Presently, x-rays are one of the most widely used methods of flaw detection. Be-
cause of the potential cumulative risk from body cell ionization use of radiation must
be limited and thus ultrasonic diagnosis is preferable. Ultrasonic imaging provides
anatomic information based mostly on the differences in impedance at the interface.
The echo amplitude of ultrasound is directly proportional to the difference in acous-
tic impedance across the boundary. Unfortunately, it has been shown that there is
no significant variation of acoustic impedance between normal and cancerous tissues
[77]. Thus it is extremely difficult to differentiate between cancerous tissue and nor-
mal tissue. On the other hand, tumors have quite different attenuation characteristics
and propagation velocity.

In order to identify cancerous tissue, additional features must be extracted from
the reflected signals, such as velocity, attenuation, peak frequency, bandwidth, etc.
Once these features have been extracted, the problem is to process the information to
provide tissue classification. Pattern recognition and segmentation are two techniques

used for this classification.

6.3.1 Ultrasonic Imaging of Tissues

Since most ultrasonic imaging systems utilize echo returns from boundaries of dif-
ferent acoustic impedance, the image shows only the outline of the boundaries rather
than the material properties. However, since the echo return contains not only bound-
ary information but also the acoustic wave propagation dispersion and attenuation
which are both frequency dependent, techniques for evaluation of acoustic velocity

and attenuation of a single medium have been investigated [10]. Typically, a sequence

93

of measurements are made at discrete frequencies over the range of interest [16]. The
procedure is quite tedious and time consuming, so instead of using discrete frequencies
for each echo measurement, a narrow broad bandwidth pulse is used to investigate the
properties of tissue. Each frequency component of the transmitting signal spectrum
propagates through the material with a different velocity and attenuation depending
on the properties of the medium. Consequently, comparing the pulse shape of the
echo return to that of the known incident pulse will allow us to identify the target
material having different properties such as density, velocity and attenuation.

The sample used for experimental demonstration consists of a section of human
brain with a hemorrhaged tumor [75], as shown in Figure 6.6. This sample is fixed
in Formalin and packed inside a plastic bag. The hemorrhaged part of the brain is in

the upper left corner of the sample.

In order to evaluate the sample the ultrasonic scanning system previously de-
scribed, and shown in Figure 5.9, is used. A Panametrics 5050 pulser/receiver is
used to excite the transducer and receive reflected signals. The received signal is
sampled at 40 MHz with 8-bit resolution. The scanning area is 100 mm by 100 mm,
and the scanning separation is 2 mm in both the X axis and the Y axis. For each
position, 400 sampling points are recorded and the sampling data is stored in the
computer for further analysis. Since the original images have some noise and low
resolution, some image processing techniques are used to enhance the images. The

details of the feature selection and image enhancement will be discussed next.

6.3.2 Feature Selection

In order to obtain maximum information from a sample, five features from the
ultrasonic echo return are extracted; they are the total energy, central frequency, fre-

quency at which the peak amplitude occurs, 3-dB bandwidth of the echo spectrum,

94

and the correlation coefficient between the incident and reflected signals. These fea-
tures contain not only the information on acoustic impedance variation but also the

attenuation and velocity characteristics.

Total energy
The first feature used is the total energy which can be obtained from a time
domain signal 3(n). We can define the total energy as the summation of the square
of each sampling value.
N
Total Energy = Z [s(i)|2 (6.14)
i=1
The total energy of the reflected signal is related to the reflection coefficient, which
contains information on the impedance of the medium. Therefore, the total energy
will provide us with information on the impedance of the tissue. An image constructed
by utilizing total energy is shown in Figure 6.7. From Figure 6.7 we can see the shape
of the brain is very close to the sample. The gray matter of the brain is displayed as

dark blue in the image.

Central frequency, frequency with peak amplitude, and bandwidth

It has been shown that the central frequency and frequency with peak amplitude
of reflected signals will shift according to the attenuation coefficient of the material
[22]. The bandwidth of the reflected signal will also change. Figure 6.8 shows
these three features in a typical ultrasonic spectrum, obtained from the Fast Fourier
Transform. The central frequency is determined by the center of the 3-dB bandwidth.
These three features will provide information on the attenuation properties inside the
tissue. The images constructed by these three features are shown in Figures 6.9, 6.10,

and 6.11. For a symmetrical spectrum the central frequency will be close to the peak

95

 

Figure 6.6. A section of human brain with hemorrhaged tumor

 

Figure 6.7. Ultrasonic imaging constructed using the total energy of the reflected
signals

96

frequency, but for a typical ultrasonic signal, the spectrum is not symmetrical, so that
the central frequency and the peak frequency are not necessarily equal. The three
features are fairly closely correlated since the central frequency is determined by the
center of the 3-dB bandwidth. From the images obtained from the peak frequency
and 3-dB bandwidth, the gray matter is not as clear as the image obtained from total
energy, but the upper part and upper left corner did show a difference with respect to
other areas. This is due to the difference in the attenuation coefficient of tissue close
to the tumor versus the attenuation coefficient in other areas. In the image obtained
from the central frequency, the area of gray matter is still clear, but other areas are

not as well defined.

Correlation coefficient between the incident and reflected signal
The correlation feature shows the difference between the incident signal and the
reflected signal in the time domain. From the shape of the reflected signal the prop-
erties of the medium, such as elasticity, stiffness constant, velocity, and attenuation
can be determined. The image constructed by the correlation feature is shown in Fig-
ure 6.12. The image obtained from the correlation coefficient is close to the image
obtained from the total energy. The shape and the gray matter is very clear in this

image.

6.3.3 Image Enhancement

Noise represents a major problem for generating ultrasonic images. In order to
reduce the effect of noise, a modified median filter has been used to reduce the noise in
images. This two-dimensional modified median filter has been found to be very pow-

erful in removing noise from two-dimensional signals without blurring the edges [62].

Magnitude

97

x10'5 Frequency spectrum
1 .6 l I I fl I I I I I

 

1.4 -

    

I

1.2
....... ,-_---- 34113 Bandwidth

 

0.8

I

0.6

0.4

0.2

 

 

 

 

h------------------------.
-----------p---------------

l

0 0.5 1 1.5 2

)-
-

 

.N
u.

3 3.5 4 4.5 5
Frequency (Hz) 11(10‘5

Figure 6.8. Central frequency, peak frequency, and bandwidth of a typical ultrasonic
spectrum

98

 

Figure 6.9. Ultrasonic image constructed using the peak frequency of the reflected
signals

 

Figure 6.10. Ultrasonic image constructed using the central frequency of the reflected
signals

99

 

Figure 6.11. Ultrasonic image constructed using the 3-dB bandwidth of the reflected
signals

 

Figure 6.12. Ultrasonic image constructed using the cross correlation of the incident
signal and the reflected signals

100

In computing the median, the set of pixels is restricted to those with a difference in
gray level no greater than some threshold.

For a high resolution of the images, a two-dimensional interpolation program has
been used to enhance the images. This technique will keep the smooth edges while
magnifying images. Assume we have an R by C image, and the intensity of pixel
(2°, j ) is AM. If the size of the magnified image is (RK,) by (CKC), the intensity of

the new image,B;.K,+m,J-.Kc+n, will be:

1
BisKr+mJ¢Kc+n = fi— { (K, — m)(Kc — ")Ai+1.j+1 + m(Kc " n)A.-+1,,-+1 +

"(Kr - m)Ai+l.j+1 + nmAi+1.j+1 } (6-15)

where 0§i<R
OSj<C
OSm<Kr

OSn<KC

Figure 6.13 shows an original image without using enhancement techniques. This
image was obtained from a composite material with a hole in the center. Each pixel
in the image is displayed as a box, therefore the resolution is not good. For an image
using the enhancement technique as shown in Figure 6.14. The edges become smooth

and an overall image quality improvement is obtained.

6.3.4 Tissue Characterization by Hierarchical Clustering

Once the information of each feature is obtained a pattern matrix can be con-
structed for further analysis. A pattern matrix is a matrix which contains the in-

formation of each feature in a column and all information obtained from a certain

101

 

Figure 6.13. The original image obtained from a composite material with a hole in
the center

 

Figure 6.14. The enhanced image by interpolation

102

location in one row. In our experiments, we have 2601 rows, (51 by 51 locations),
and 5 columns (5 different features). We can use an eigenvector projection to see the
relationship between these five features. In the eigenvector projection, the two eigen-
vectors with the two most significant eigenvalues are chosen. The projection of five
features in a two-dimensional space is shown in Figure 6.15. From Figure 6.15 we
can observe that the power and cross correlation provide similar information, while

the peak frequency and central frequency are close.

6.3.5 Data Reduction

Since the number of pixels in the images is large, it will take a long processing time
and large memory space. The first step is to remove the informationless background
of the images. This process will reduce the number of patterns necessary and thus
simplify the problem with decreased processing time. The background of the images
can easily be removed because it is very uniform. For the images obtained from
the sample of human brain, the original image had 2601 pixels, while the image
without background has only 1905 pixels. By these 1905 pixels a proximity matrix
which is 1905 by 1905 will be constructed using the normalized pattern matrix with
uncorrelated features. The proximity matrix will be symmetrical and the diagonal
entries are all zero. Due to a 17% decrease in the number of pixels a 47% decrease
in the size of the proximity matrix will result. Therefore the processing time will be

greatly improved.

6.3.6 Image Reconstruction

Once the sequence of partitions is obtained, the pixels can be classified as numbers
of groups desired. By giving each group of pixels a new label, new images can be

constructed. The results of clustering by the single link, the complete link, Ward’s

103

Projection of 5 Different Features

 

 

 

 

10 I I l
8 — _
. Correlation
6 - . Total Energy ‘
4 " th “
2 — -
N
'3 0“ ‘
-2 _ -
.4 .. -
‘ Peak Freq.
‘ Central Freq.
-3 - A
_10 l l l
-10 -5 O 5 10

Figure 6.15. Projection of five features in a two-dimensional space

104

method, and the UPGMA are shown in Figures 6.16, 6.17, 6.18, and 6.19, respec-
tively. In Figure 6.16, which is constructed by the single link method, the Size of
most clusters are small and separated, so it does not provide good information about
the properties of tissue. We can only see that the upper part of the image is different

' from the other part of the tissue.

Figure 6.17 is constructed by the complete link method. This image is classified
as 4 clusters including two clusters are mixed together along the edge of the sample.
The distribution of these two clusters are close to the gray matter of the brain. This

image looks better than the image obtained from the single link method.

The result of the Ward’s method is shown in Figure 6.18. The gray matter and the
area of tumor is clearly shown in the image, but some small clusters are distributed

over the samples.

Finally the result of the UPGMA method is shown in Figure 6.19. The clustering
result by this method looks pretty good. We can clearly see the gray matter, the
white matter, and the hemorrhaged tumor which are all verses close to the original

sample.

6.3.7 Cluster Validation

The discussion above focussed on direct observation for image validation. This
may not be objective and is certainly not a quantitative way to evaluate the clustering
results. In order to determine which method best fits the original data, the internal
index, CPCC, of each clustering result is calculated. The value of the CPCC by
different clustering methods for the section of human brain with a hemorrhaged tumor

are shown in Table 6.2 and Figure 6.20.

From Table 6.2, the value of the CPCC by the single link method is only 0.07,

105

 

Figure 6.16. Clustering result of a section of human brain with a hemorrhaged tumor
by single link clustering method with 8 clusters

 

Figure 6.17. Clustering result of a section of human brain with a hemorrhaged tumor
by complete link clustering method with 4 clusters

106

 

Figure 6.18. Clustering result of a section of human brain with a hemorrhaged tumor
by Ward’s clustering method with 9 clusters

 

Figure 6.19. Clustering result of a section of human brain with a hemorrhaged tumor
by UPGMA clustering method with 3 clusters

107

which is very low. This means that the results obtained from the signal link in this
experiment is not good but instead is close to random noise. Thus, we should not
trust this clustering result.

The results obtained from the complete link and the Ward’s method have higher
values of CPCC. This implies that the result obtained by the complete link and by
Ward’s method are better than the result obtained by the single link method.

Comparing with the above three methods with the UPGMA method indicates that
the latter method has the highest value of CPCC. Therefore, the results obtained from
the UPGMA method should provide the best information about the original data.
This is confirmed by direct observation of the images. We can see from Table 6.2

and the constructed image that the best result for our experiments is obtained by the

UPGMA method.

108

 

Method Value

 

Single Link 0.07
Complete Link 0.18
Ward’s method 0.27

UPGMA 0.71

 

 

 

Table 6.2. The values of CPCC by different clustering methods for a section of human
brain with a hemorrhaged tumor

 

00500

109

Internal Index

CPCC
0.8 _--....

 

 

 

 

 

 

 

 

 

 

 

 

SL CL Ward’s UPGMA

Figure 6. 20. Comparison of CPCC by the different methods for a section of human
brain with a hemorrhaged tumor

CHAPTER 7

CONCLUSIONS AND
DIRECTIONS FOR FUTURE
WORK

7 .1 Conclusions

In this dissertation, some basic principles of acoustic waves have been reviewed and
the techniques to find velocity and attenuation information by broadband signals have
been derived and experimentally verified. Based on the unique spectral reflectance
features of different materials, the velocity-density product and attenuation-density
ratio can be obtained for multilayered models. From this information, material iden—
tification using the experimentally measured signal returns was achieved.

Acoustical images were obtained by broadband signals derived from narrow time
video pulses. The images were characterized by hierarchical clustering techniques.
As a result, the technique developed is promising for applications to nondestructive

identification, evaluation of materials, and biological applications in general.

110

111

7 .2 Directions for Future Work

Acoustical images are very important in nondestructive evaluation and biological
applications. There are many topics to pursue in future research. For example, a
three-dimensional acoustic image could be constructed with a processing computer
which is sufficiently large and fast. Also a real time acoustical image system could be
designed with digital signal processing chips. The accuracy of tissue characterization
could be improved by building a database of the response from different biological
tissue. Finally, using the technologies of artificial neural networks and artificial intel-

ligence, an automatic diagnosis system could be constructed.

APPENDIX

 

APPENDIX A

PROGRAM LIST - SCANNING
AND COLOR DISPLAY

/*****************************************************************¥******

* a
* Scanning and Color Display: *
* at
* Control motors X t Y to take data, find an echo in a certain *
* range, and Show the image on screen. *
* *
* By Nathan N. H. Hang at:
* *

*ttttttt*********#************************#******t*****t****tt******#*#*/

#include <3tdio.h>
tinclude <3tdlib.h>
tinclude <math.h>
tinclude <con10.h>
tinclude <ctype.h>
tinclude <graphice.h>

#define VERSION "Color Scan V 9.3"
Idefine AUTHOR "Kai-Helen Hang"
tdefine H_DATE "Apr. 29, 1991"

/* tdefine DEBUGECHO 1*/
It tdefine DEBUGNOINP 1*/
It tdefine DEBUGPLOT 1*/
tdefine SIZEX 61
idefine SIZEY 61

112

idefine

tdefine
#define
tdefine
#define
tdefine

#define
#define
#define
#define
#define
#define

#derine
#define

113

idefine BEEP
idefine CTRLC
idefine ESC

unsigned int
signed int
signed int
float
unsigned char

unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned int
unsigned char
unsigned int
signed int
unsigned int
unsigned int
float

signed int
long

HAXSIZE 6600

MAIN 5

HAXERR 300

SIGTHRE 10

CHECKHT 0.2

NTRY 2

HOVE-PDS 1

MOVE-NEG -1

MOVEX 0

HOVEY 1

XSTEP 0.033333

YSTEP 0.033333

STANDARD "INPUTX.DTA"

ECHOTEST "ECHOTEST.DTA"
’\OO7’
’\003’
27

sumEMAXSIZE];
sigEMAXSIZE];

CESIZEXJESIZEY];
chtESIZEXJESIZEY];

isumEHAXSIZE];

Israte;
Start;
Length;
Trigger;
Neignal;
Xmax;
Ymax;
Ster;
StepY;
Xdot;
Ydot;
Xdot2;
Ydot2;
thabs;
Delay;
inpxCSOO];
inpx-length;
pmax_inpx;
smean;
ismean;
xt2sum;

ls
It
It
Is
[a
It
/*
It
It
/s
It
la
It
It
It
It

index of sampling rate */
starting point of sampling */
length of sampling interval */
trigger level of sampling */

No. of signals to be averaged*/
max. step of motor X *I
max. step of motor Y */
scanning step of motor X */
scanning step of motor Y */
display points in X axis */
display points in Y axis */
distance between two ponits X*/
distance between two points Y*/
flag to take abs or not */
delay time of stepper motors */
input xt */

FILE *fopdata;
FILE *fipdata;
char fnamef40];

int graphdriver;
int graphmode;

int th,wOb,wOr,wOl,w0clip;
int w1t,w1b,w1r,w11,w1c1ip;
int w2t,w2b,w2r,w21,w2clip;

114

/* pointer of output data file
/* pointer of input data file

/* graphics card, eg. Hercules.

/* graphics mode

/* window 0: color bar
/* window 1: image
/* window 2: command

int pcolorf9] = {RED,LIGHTRED,YELLOH,LIGHTGREEN,LIGHTCYAN,
CYAN,LIGHTBLUE,BLUE,BLACK};

int threfii];
float htmin;
float htmax;
float htscale;
float htscalez;
int x0,y0;

int 1x,1y;
float fthreEii];

void main()

{

clrsch);

/* normalized threshold
/* real value of threEO]

/* real value of threEiO]

/* scale of threE]

It scale of threE] for display

/* origin of the image
/* size of the image

printf("Xs :\n%8. %8\n", VERSION. AUTHOR. H-0ATE);

mainfuncC);
fcloseallC);
}

mainfuncC)

{
unsigned int ansflag;
unsigned char ans;
unsigned char saveall;
unsigned char datafile;
int i;

ansflagso;
thre[10]=0.0;
htmin=0.0;

while( ansflag == 0 )
{

/* end of main()

printf("\n\n D: Display image from a file.\n");
printf(" T: Take data from the scanning system.\n");
printf(" F: Take data from a data file.\n");
printf(" M: Get monochrome image by threshold.\n");

*/

*/
*/

.*/

*l

*l
*l
*/

*l
*l
*l
*I
*l
*/
*l

115

printf(" Which one ? ");

ans=toupper(getch());

if ( ans == ’D’ ll ans == ’T’ I] ans == ’F’ I] ans ==’M’) ansflag=1;
else printf("Xc",BEEP);

}

printf("\n\n");

switch (ans)

{

case ’D’:

Xdot=4;
Ydot=4;
Xdot2-4;
Ydot2-4;
loadimgC);
getcommandC);
break;

case ’T’:
loadxtC);
for(i=0;i<10;i++)

thre[i]=(int)((9-i)*256/9.0+.5):

initmotor();
resetmotorC);
inputdata(&saveall,tdatafile);
scandata(saveall,datafile);
getcommandC);
break;

case ’F’:
readdata();
getcommand();
case ’H’:
setthre();
getcommand();
} /* end of switch */
} /* end of mainfunc() */

readdataf)

{
int 1;
char fname1[40];
char tmp;
unsigned char datafile;

loadxt();

for(i=0;i<10;i++) thre[i]=(int)((9-i)*256/9.0+.5);
printf("\nData filename =?");

scanf("%s", fname1);

116

fipdata = fopen(fname1, "rb");
if (fipdata <= 0)
{
printf("\n 2c Cannot open file: Zs.\n",BEEP,fname1);
return (-1); /* cannot open file, return */
}
fscanf(fipdata,"%d %d %d 2d 2d 2d Xd 1d ",tIsrate,&Xmax,&Ymax
,tSter,&StepY,&Start,lLength,&Nsigna1);
fscanf(fipdata,"%d 2d 1d 2d 1d 2d Zg",&Xdot,&Ydot,thot2
,tYdot2,&Trigger,&thabs,thtmax);
fscanf(fipdata,"%c",ttmp); It delete LF */
if (thabs >90) datafile=’S’;
thabsstoupper(thabs);

scandata2(datafile);
} /* end of readdata() */
setthre()
{
int 1;
char fname1[40];
char tmp;

unsigned char datafile;

for(i=0;i<10;i++) threfi]'(int)((9-1)*2$6/9.0+.5);
printf("\nData filename .?n);
scanf("%s", fname1);
fipdata I fopen(fname1, "rb");
if (fipdata <= 0)
{
printf("\n Xc Cannot open file: Zs.\n",BEEP,fname1);
return (-1); /* cannot open file, return *I
}
fscanf(fipdata,"%d %d %d 1d 2d 2d %d 2d ",tIsrate,&Xmax,&Ymax
,tSter,&StepY,&Start,lLength,&Nsignal);
fscanf(fipdata,"%d Zd 2d Zd %d %d Ze",&Xdot,&Ydot,&Xdot2
,&Ydot2,&Trigger,&thabs,thtmax);
fscanf(fipdata,"%c",ttmp); It delete LF */
if (thabs >90) datafile=’S’;
thabs=toupper(thabs);

scandata3(datafile);
} /* end of setthre() */
loadimgC) /* load image from a file */
{
FILE *fipimg; /* pointer of input image file */
int 19:];
char fname2[40];

 

117

closegraph();
printf("Input image data filename = ? ");
scanf("Xs",fname2);
fipimgsfopen(fname2,"r");
if (fipimg <= 0)
{
printf("\n 2c Cannot open file: Xs.\n",BEEP,fnam62);
printf("\n Press any key to continue.\n");
getch();
}
else
{
strcpnyname,fname2);
fscanf(fipimg,"%d Zd Xd 2d %d If If",tXmax,lYmax,&Ster,&StepY,
chtabs,&htmax.&htmin);
focanfoipimg,"Xd 2d 1d 1d 1d",lIsrate,&Start,&Length,

tTrigger,&Nsigna1);
for (i-O;i<9;i++)
{
fscanf(fipimg,"%d ",tpcolor[i]);
}
for (i-O;i<10;i++)
{
fscanf(fipimg,"%d ",tthreEiJ);
}
for (j-o;j<-Ymax;j++)
{
for (i=0;i<=xmax;i++)
{
fscanfoipimg,"%d Xe ",&c[i][j].&cht[i][j]);
}
}
fcloseCfipimg);
}
showall();
} /* end of loadimg() */
saveimg() /* save image to a file */
{
FILE *fopimg; /* pointer of output image file */
int i,j;
char ans;
char fname2[40];
closegraphC);

printf("Output image data filename = ? ");
scanf("%8",fname2);
if (access(fname2,0) == 0)

118

{
printf("\nFile [23] exists, Replace (Y or N) ? ",fname2);
do {
ans=toupper(getch());
if (ans 8' ’Y’) printf("\nFile will be replaced.\n");
if (ans 8- ’N’)
{
printf("\nFile not saved, press any key to return.\n");
getch();
showall();
return(-1);
}
} while(ans != ’Y’ at ans != ’N’);
}
fopimg-fopeannameZ,"w");
if (fopimg <= 0)
{
printf("\n 2c Cannot open file: Zs.\n",BEEP,fnam92);
printf("\n Press any key to continue.\n");
getch();

else

strcpy(fname,fname2);
fprinthfopimg,"%5d 25d 16d 15d 15d 2g Xg\n",Xmax,Ymax,Ster,StepY,
thabs,htmax,htmin);
fprintf(fopimg,"%5d 25d 15d 15d %5d\n",Israte,Start,Length,
Trigger,Nsignal);
for (i-O;i<9;i++)
{
fprintf(fopimg,"%5d ",pcolorEiJ);
}
fprinthfopimg,"\n");
for (i=0;i<10;i++)
{
fprinthfopimg,“%5d ",threEiJ);
}
fprinthfopimg."\n"):
for (j-O;j<=Ymax;j++)
I
for (i=0;i<=Xmax;i++)
{
fprintf(fopimg,"26d 26.3f \n",c[i][j],cht[i][j]);
}
}
fclose(fopimg);
}
showall();

119

} /* end of saveimg() */
loadxt() /* load x(t) from INPUTX.DTA */
{
FILE *fip;
unsigned int i;
float inpr;
float samplerate; /* samplerate,israte,startpoint,*/
unsigned int israte; /* ,and j1 are not used in */
unsigned int startpoint; It this subroutine. They are set*/
unsigned int j1; /* to be compatible with the *l
/* data file which are generated*/
/* by AD?.c. I"I
float stest;
FILE *ftest;
fip B fopenCSTANDARD,"r");
if (fip<-O) /* Cannot open file *I

{

Printf("\n\n%cCannot find standard input file Zs !\n\n",

BEEP,STANDARD);

return(-1);

}

printf("\nRead standard signal from Zs.",STANDARD);

fscanfoip,"%f 1d Xd Xd If ".
tsamplerate,kisrate,&startpoint,&j1,tsmean);

for (i I O; lfeof(fip); i++)

{
fscanf(fip, "Xf ", linpr);
inpri]=inpr-smean+0.5;
}
fclose(fip);

inpx_1ength=i-1;
pmax-inpx-O;
xt23um=0;
for (i=0; i<inpx_length; i++)
{
if (abs(inpx[i]) > abs(inpx[pmax_inpx])) pmax_inpx=i;
xt2sum=xt2sum+inpx[i]*inpri];
} /* end of for */
#ifdef DEBUGECHO
ftestsfopenCECHOTEST,"r");
printf("\nRead data from ECHOTEST.DTA. ");
for (i = O; lfeof(ftest); ++i)

120

{
fscanfotest, "1f", tstest);
sig[i]=(unsigned int)(-(stest-smean)*0.75+smean);
printf("xd ",sigEiJ);
}
fclose(ftest);
Length-i-i;
Start=0;
Cendif
ismean=(int)(smean+0.5);
} /* end of loadxt() */

inputdata(unsigned char *psave,unsigned char *pdatafile)

/* get setup parameters */
{
unsigned char saveall;
unsigned int loopflag;
char fname1[40]; /* filenmae of output file */
float xmaxl;
float ymaxl;
float xstepl;
float ystepl;
char ans;
do {
printf("\nChoose the sampling rate:\n");
printf("\n 0: 2OHHz 8: 4OHHz");
printf("\n 1: ZHHz 9: 4HHz");
printf("\n 2: 2OOKHz 10: 4OOKHz");
printf("\n 3: ZOKHz 11: 4OKHz");
printf("\n 4: 2KHz 12: 4KHz");
printf("\n\n Sampling rate I ?");

scanf("Xd",&Israte);
} while (Israte<0 ll Israte>12 ll ( Israte>4 kt Israte<8));

Cifndef DEBUGECHO

printf("\n\nStarting point of sampling - ?");

scanf("%d",&8tart);

printf("\nLength of sampling points (max no. a 1d) = ?",HAXSIZE);

scanf("%d",&Length);

printf("\nTrigger signal offset (0-255, 128 = 0 volt) = ?");

scanf("%d",&Trigger);
#endif

thabs=’Y’;

printf("No. of signals to be averaged = ?");

scanf("%d",&Nsignal);

do

{

loopflag=0;

121

printf("\nHax scanning range in X direction (in mm) I ?");
scanf("%f",&xmaxl);
printf("Hax scanning range in Y direction (in mm) = ?");
scanf("%f",&ymaxl);
printf("\nScanning step in X direction (in mm) = ?");
scanf("%f",&xstepl);
printf("Scanning step in Y direction (in mm) = ?");
scanf("%f",&ystep1);
Xmax=(int)(xmax1/xstepl);
YmaxI(int)(ymaxl/ystepl);
Ster=(int)(xstepl/XSTEP);
StepYI(int)(ystepl/YSTEP); I
printf("\nXmax I 16d, xstep I 16d, i
XSTEP I Zf",Xmax,Ster,Ster*XSTEP);
printf("\nYmax I 16d, ystep I 26d,
YSTEP = Zf\n",Ymax,StepY,StepY*YSTEP);
if (Xmax>SIZEX)

 

{
printf("\nToo many points in X, ", h
"max. points in X is 2d, try again\n",
SIZEX);
loopflag=1;
}
if (Ymax>SIZEY)
{

printf("\nToo many points in Y, ",
"max. points in Y is Xd, try again\n",
SIZEY);
loopflain;
}
} while (loopflag II 1);

printf("\nDisplay size of each pixel in X axis I ?");
scanf("2d",&Xdot);
printf("Display size of each pixel in Y axis I ?");
scanf("%d",&Ydot);
Xdot2IXdot;
Ydot2IYdot;
printf("\nflax. ht I ?");
scanf("%f",&htmax);
printf("Delay time (msec) of stepper motors (>= 3) = ?");
scanf("%d",&Delay);
do {
printf("\nSave all sampling points (Y or N) = ?");
saveall=getch();
printf("2c\n",saveall);
*psave=toupper(saveall);
ansI’Y’;

if (Ipsave II ’Y’)
{
printf("\nSmall size of Large size data file (8 or L) I ?");
*pdatafi1e=toupper(getch());
if (*pdatafile II ’S’) printf("Small data file.\n");
else printf("Large data file.\n");
printf("\nOutput filename I?");
scanf("%s", fname1);
if (access(fname1,0) II 0)
{
printf("\nFile [Is] exists, Replace (Y or N) ? ",fnamei);
ans=toupper(getch());
if (ans == ’Y’) printf("\nFile will be replaced.\n");
if (ans II ’N’) printf("\nNot replaced\n");
} while(ans !I ’Y’ n ans !I ’N’);
if (ans II ’Y’)
{
fopdata I fopen(fname1, "wb");
if (fopdata (I O)
{
printf("\n 2c Cannot open file: Xs.\n",BEEP,fname1);
return (-1); [I cannot open file, return I/
I
if (Ipdatafile II ’S’) thabsItolower(thabs);
[I if thabs is captial, then the data I"I
[I file is large, otherwise is small I/

 

fprintf(fopdata,"%6d XSd 15d 25d 25d 25d 15d 25d\n",
Israte,Xmax,Ymax,Ster,StepY,Start,Length,Nsignal);
fprintf(fopdata,"25d 25d 25d 25d 25d 15d 211.4e\n",
Xdot,Ydot,Xdot2,Ydot2,Trigger,thabs,htmax);
thabs=toupper(thabs);

}
}
} while(ans == ’N’);
} /* end of inputdata() */
scandata(saveall,datafile) /I control motors and get data */

unsigned char saveall;

unsigned char datafile;

{
signed int xp;
signed int yp;
unsigned int i,j;
float err;
unsigned int goodsig;
unsigned int ntotal;
unsigned int echo;

123

float echoht;
int try;
int tryflag;

detectgraph(&graphdriver,agraphmode);
initplotC);

for (ypIO; yp <I Ymax; yp+I2) {
Cifndef DEBUGECHO
tryIO;
do
{
try++;
tryflagIO;
takesample(&echo,&echoht,lgoodsig,&ntota1);
II sampling at (O,yp) */
if ( (yp>o) it (try <uraY) )
{
if ((fab8((doub16)(echoht-cht[0][yp-1])) > CHECKHT) it
(fabs((double)(echoht-cht[1][yp])) > CHECKHT))
tryflain;
}
} while ( (tryflagIIi) Ct (try<NTRY));
telse
findiecho(sig,&echo,kechoht);
tendif

movemotor(HOVEX,NOVE_POS,Ster,Delay);

if (saveall II ’Y’) saveallp(0,yp,sum,goodsig,ntotal,datafile);
C[0][yp]=echo;

cht[0][yp]=echoht;

display(0.yp);
for (pr1; xp <I Xmax; xp++)
{
Cifndef DEBUGECHO
try-O;
do
{
try++;
tryflagIO;

takesample(&echo,&echoht,lgoodsig,&ntotal);
II sampling at (xp,yp) */
if ( try (NTRY )
{
if (yp >0)
{
if ((fabs((double)(echoht-chtfxp][yp-1])) > CHECKHT) &&

 

124

(fabsCCdouble)(echoht-chtExp-i][yp])) > CHECKHT))
tryflagI1;
}
else
if (fabs((doub1e)(echoht-chtfxp-il[yp])) > CHECKHT)
tryflag=1;
}
} while ( (tryflagIII) ll (try<NTRY));
telse
findiechoCsig,techo,&echoht);
Cendif
if ( xp < Xmax ) movemotor(HOVEX,NOVE_POS,Ster,Delay);
else movemotor(HOVEY,NOVE_POS,StepY,Delay);
if (saveall II ’Y’) saveallp(xp,yp,sum,goodsig,ntotal,datafile);
CEprEypl-echO;
cht[xp][yp]=echoht;
display(xp.yp);

( yp < Ymax )

Sifndef DEBUGECHO
try-O;
do
I
try++;
tryflagIO;
takesample(&echo,&echoht,tgoodsig,&ntotal);
II sampling at (Xmax,yp+1) */
if ( try (NTRY )
{
if (fabs((double)(echoht-chtEXmax][yp])) > CHECKHT)
tryflag=1;
}
} while ( (tryflag==1) a: (try<NTRY));
Celse
findiecho(sig,kecho,&echoht);
Iendif
movemotor(NOVEX,HOVE_NEG,Stepx,Delay);
if (saveall II ’Y’)
saveallp(Xmax,yp+1,sum,goodsig,ntotal,datafile);
c[Xmax][yp+1]Iecho;
cht[Xmax][yp+1]Iechoht;
display(Xmax,yp+1);

for (xp=Xmax-1; xp >= 0; xp--)
{
#ifndef DEBUGECHO

125

try=0;
do
{
try++;
tryflag=0;
takesample(techo,techoht,tgoodsig,&ntota1);
II sampling at (xp.yp+1) */
if ( try (NTRY )
{
if ((fabs((double)(echoht-chtfxp][yp])) > cascxnr) a:
(fabs((double)(echoht-chtfxp+1][yp+1])) > CHECKHT))
tryflag=1;
}
} while ( (tryflag==1) Ia (try<NTRY));
telse
findiecho(sig,techo,&echoht);
Cendif

 

if ( xp > O ) movemotor(NOVEX,HDVE-NEG,Ster,Delay);
else movemotor(MOVEY,MDVE_POS,StepY,De1ay);

if (saveall == ’Y’) .
saveallp(xp,yp+1,sum,goodsig,ntota1,datafile);

cExp][yp+1]Iecho;

chtExp][yp+1]Iechoht;

display(xp,yp+1);
} [I end of for xp */
} [I end of if (yp < Ymax) */
} /* end of for loop (yp) */
if (saveall I: ’Y’) fclose(fopdata);
} [I end of scandata() */
scandata2(datafile) [I readdata from file II
unsigned char datafile;

{
signed int xp;
signed int yp;
unsigned int i,j;
float err;
unsigned int goodsig;
unsigned int ntotal;
unsigned int echo;
float echoht;

detectgraph(&graphdriver,tgraphmode);
initplotC);

for (ypIO; yp <2 Ymax; yp+I2) {
loadallp(datafile);

126

findiecho(sig,techo,techoht); II echo at (O,yp)
c[0][yp]Iecho;

cht[0][yp]Iechoht;

displayCO.YP);

for (pri; xp (I Xmax; xp++)

{
loadallp(datafile);
findlecho(sig,&echo,aechoht); /* echo at (xp,yp)
c [xp] [yp] =echo ;
chtfxp][yp]=echoht;
display(xp,yp);

if ( yp < Ymax )

loadallp(datafile);

findiecho(sig,&echo,aechoht); /* echo at (Xmax,yp+1)
c[Xmax][yp+1]Iecho;

cht[Xmax][yp+1]=echoht;

display(Xmax,yp+1);

for (prXmax-i; xp >I O; xp--)
{
loadallpCdatafile);

find1echo(sig,&echo,techoht); II echo at (xp,yp+1)

c[xp][yp+1]Iecho;
chtExp][yp+1]Iechoht;

display(xp,yp+1);
} [I end of for xp
} [I end of if (yp < Ymax)
} [I end of for loop (yp)
} II end of scandata2()
scandata3(datafile) /* readdata from file

unsigned char datafile;

{

signed int xp;
signed int yp;
unsigned int i,j;
float err;
unsigned int goodsig;
unsigned int ntotal;
unsigned int echo;
float echoht;
unsigned int flagi;
unsigned int flag2;
int stemp;

*/
*l
*/
*l

*I

*/

 

127

float upthre;

float lowthre;
int iupthre;
int ilowthre;

printf("\nInput the range of threshold.\n");
printf("Upper bound of the signal (0.0 - 1.0) = ?");
scanf("%f",&upthre);

printf("Lower bound of the signal (0.0 - Zf) = ?",upthre);
scanf("%f",&lowthre);

iupthreaupthreI128;

ilowthre=lowthre*128;

detectgrapthgraphdriver,tgraphmode);
initplot();

for (yp=0; yp <I Ymax; yp+=2) {
loadallpCdatafile);
flagiIO;
flag2I0;
for (i=0;i<Length;i++)
{
stemp=abs(sig[i]-128);
if (stemp >I ilowthre) flagiIi;
if (stemp >I iupthre) flag2=1;
}
if ((flagl == 1) a: (flag2 == 0))
{
CEO] [yp1=1;
cht[0][yp]=lowthre*htmax;

CEOJEyplIO;
chtEOJEprIO;
}
display(0,yp);

for (XPII; xp <I Xmax; XP++)
{
loadallp(datafile);
flagiIO;
flag2=0;
for (i=0;i<Length;i++)
{
stemp=abs(sig[i]-128);
if (stemp >= ilowthre) flag1=1;
if (stemp >= iupthre) flag2=1;

 

128

}
if ((flagl as 1) && (flag2 == 0))
{
CExplfypl-iz
chtExp][yp]=lowthre*htmax;

CEprEpr=0;
chtExplfyp1=os
}
displanyp.yp);

( yp < Ymax )

loadallp(datafile);
flagiIo;
flag2I0;
for (iI0;i<Length;i++)
{
stemp=abs(sig[iJ-128);
if (stemp >= ilowthre) flag1I1;
if (stemp >I iupthre) flag2I1;
}
if ((flag1 == 1) as (f1ag2 == 0))
{
CEXmaxltyp+11I1;
chtEXmax][yp+1]Ilowthre*htmax;

CEXmaXJEyp+1]-0;
chthmax][yp+1]=O;
}
displayCXmax,yp+1);

for (xp=Xmax-1; xp >= 0; xp--)
{
loadallp(datafi1e);
flag1=0;
flag2I0;
for (iI0;i<Length;i++)
{
stemp=abs(sig[iJ-128);
if (stemp >= ilowthre) flag1=1;
if (stemp >= iupthre) flag2=1;

129

if ((flagi == 1) && (flag2 == 0))

{
CExplfyp+13I13
chtExp][yp+1]=1owthre*htmax;
}
else
{
c[xp][yp+11=0;
chtfxplfyp+1JI0;
}
display(1p.yp+1);
} [I end of for xp */
} [I end of if (yp < Ymax) */
} [I end of for loop (yp) */
} /* end of scandata3() */
saveallp(xp,yp,sum,goodsig,ntotal,datafile)
signed int xp;
signed int yp;
unsigned int sum[];

unsigned int goodsig;
unsigned int ntotal;
unsigned char datafile;

{

}

unsigned int i;

fprintf(fopdata,"ZcXchZc",xp,yp,goodsig,ntotal);
for (i=0;i<Length;i++) sigEiJ+Iismean;

if (datafile == ’3’)
{
if (goodsig I=1)
for (i=0;i<Length;i++) isum[i]=sig[i];
else
for (i=0;i<Length;i++) isumEi]=(sum[i]+goodsig/2.0)/goodsig;
fwrite(isum,sizeof(char),Length,fopdata);

else
{
if (goodsig I=1)
fwrite(sig,sizeof(int),Length,fopdata);
else
fwrite(sum,sizeof(int),Length,fopdata);
}

loadallp(datafile)
unsigned char datafile;

{

 

130

unsigned int i;
unsigned int xp;
unsigned int yp;
unsigned int goodsig;
unsigned int ntotal;
unsigned char tmp;

fscanf(fipdata,"%c",ltmp); /* read xp */
fscanf(fipdata,"%c",ttmp); /* read yp II
fscanf(fipdata,"%c",ltmp); I* read goosig II
goodsig=tmp;

fscanf(fipdata,"%c",ttmp); II read ntotal *I

if (datafile == ’8’)
{
fread(isum,sizeof(char),Length,fipdata);
for (iI0;i<Length;i++) sigEiJIisumEi];

else
{
fread(sig,sizeof(int),Length,fipdata);
for (i=0;i<Length;i++)
sigEi]I(sig[i]+goodsig/2)/goodsig;

}
}
initplotC)
{

unsigned int maxx;
unsigned int maxy;

int i,j,k;
int x.y;

float xmaxl;
float ymaxl;
float xstepl;
float ystepl;
int iscale;

initgraph(&graphdriver,tgraphmode,"\\tc\\");
setpalette(5,36);

setpalette(13,32);

maxxsgetmaxxC);

maxngetmaxy();

wOl I maxxI3/4; It parameters for window 0 II
th = 50; /* for color bar. */
wOr I maxx;

wOb I (maxy*3)/4;

 

131

w0clip =1;

w11 I 10; II parameters for window 1
wit I 10; II for image.

w1r w01-1;

wib I wOb;

w1clip =1;

w21 0; II parameters for window 2
w2t w0b+20; II for command line.

w2r maxx;

w2b = maxy;

w2clip =1;

if(thabs == ’Y’)
htscale=(htmax-htmin)/256.0;
else

I

htminI-htmax;
htscale=htmaxI128.0;

}
setviewport(w01,w0t,wOr,w0b,w0clip);
for ( i=0;i<9;i++)

I
k860+i*20;
setfillsty1e(SOLID_FILL,pcolor[i])3
bar(3,k,52,k+19);
}
setcolor(BROHN);
rectangle(1,59,54,241);

x=60;
setcolor(HHITE);
iscaleIO;
htscale2I1.0;

while (htmax > 10.0Ihtscale2)
I
htscale2I=0.1;
iscale++;
}
while (htmax < 1.0Ihtscale2)

{
htscale2*=10.0;

iscale--;
}
if (iscale !=0)
I

*I
*I

*I
*I

132

yI20;
x=60;
gprintf(&x,ty," 21d".iscale);
gprintf(&x,&y," x 10");
}

for (i=0;i<10;i++)
I
y=56+iI20;
fthreEiJIhtscaleIthreEi]+htmin;
gprintf(&x,&y,"%6.3f",fthreEiJIhtscale2);
}
fthre[10]=htscaleIthre[10]+htmin;

xstepl=SterIXSTEP;
ysteplIStepYIYSTEP;
xmaxl=XmaxIxstepl;
ymaxlIYmanystepl;

 

x=1;

y=280;

gprintf(&x,&y,"Size:%4.1fmm x I4.1fmm",xmaxl,ymaxl);
gprintf(&x,&y,"Res.:X4.1fmm x Z4.1fmm",xstepl,ystep1);

setviewport(w11,w1t,w1r,w1b,w1clip);
lxIXdot2I(Xmax+1);
ly=Ydot2I(Ymax+1);
xOIwir-(Cwir-wil-lx)/2)3
yOIwib-((w1b-w1t-ly)/2):

setcolor(BRDNN);

rectangle(x0+2,y0+2,x0-lx-2,yO-ly-2);
rectangle(x0+6,y0+6,x0-1x-6,y0-ly-6);
} II end of initplot() II

getcommand()
{
Idefine HESi "Use arrow keys to change threshold."
#define NESZ " I: toggle instant update."
tdefine HESS "<CR>: return to main menu."
int x;
int y;
int flag;
char ans;
int 1;
char ans2;
char ansS;

133

int oldth;
int ith;
int flag2;

unsigned char saveall;
unsigned char datafile;
int dispflag;

clearcursor();
dispcdeNHITE);
flag=1;

while (flag ==1)
I
do {
} while (kbhitC) I=0 );
ans=toupper(getch());
switch (ans)
I
case ’L’:
loadimgC);
break;
case ’S’:
saveimgC);
break;
case ’T’:
dispcmd(BLACK) ;
xI150;
yI10;
setviewport(w2l,w2t,w2r,w2b,w2clip);
setcolor(NHITE);
gprintf(&x,&y,"%s",HES1)3
gprinthtx,&y,"Xs",NES2);
gprintf(&x,&y,"%s",HE83);
xI60;
ithIi;
setviewport(w01,w0t,wOr,w0b,w0c1ip);
setcolor(GREEN);
yI56+ithI20;
gprintf(&x,&y,"%6.31",(htscaleIthre[ith]+htmin)Ihtscale2);
flag2=1;
dispflag=1;
while(flag2==1)
I
setviewport(w01,w0t,w0r,w0b,w0clip);
do I
} while (kbhit() ==0 )3
ans2=getch();
if (ans2 == ’\r’)

134

flag2=0;
setcolor(NHITE);
y=56+ithI20;
gprintf(&x,&y,"%6.3f",(htscaleIthreEithJ+htmin)Ihtscale2);
x=150;
y=10;
setviewport(w2l,w2t,w2r,w2b,w2clip);
setcolor(BLACK);
gprintf(&x,&y,"%s",MESI);
gprintr(&x,&y,"zs",HEsz);
gprintf(&x,&y,"Zs",HES3);
dispimg();
dispcmd(HHITE);
}
else
I
if (ans2 =8 ’\0’)
I
ansB=getch();
switch(an83)
I
case ’P’: I* down II
setcolorCNHITE);
y=56+ithI20;
gprinthkx,&y,"%6.3f",
(htscaleIthre[ithJ+htmin)Ihtsca1e2);
if (ith < 8 ) ith++3
setcolor(GREEN);
y=56+ithI20;
gprintf(tx,&y."26.3f".
(htscaleIthre[ith]+htmin)Ihtscale2);
break;
case ’3’: /* Up */
setcolor(HHITE);
y=56+ithI20;
gprintf(&x,&y."%6o3f".
(htscaleIthreEithJ+htmin)Ihtscale2);
if (ith > 1 ) ith--;
setcolor(GREEN);
y=56+ithI20;
gprintf(&x,&y,"%6-3f".
(htscaleIthreEithJ+htmin)Ihtscale2);
break;
case ’K’: [I left II
oldth=threfith];
if (threEithJ > threEith+11> threEith]--;
setcolorCBLACK);

135

y=56+ithI20;

gprintf(&x,&y,"%6.3f",
(htscaleIoldth+htmin)Ihtscale2);

setcolorCGREEN);

y=56+ithI20;

gprintf(&x,&y,"26.3f",
(htscaleIthreEithJ+htmin)Ihtscale2);

if (dispflag II 1 ) dispimgC);

break;

case ’H’: II right II

oldthIthreIith];

if (threEithJ < thre[ith-1]) threEith]++;

setcolorCBLACK);

yI56+ithI20;

gprintf(&x,&y,"%6.3f",(htscaleIoldth+htmin)Ihtscale2);

setcolor(GREEN);

y=56+ithI20;

gprintf(&x,&y,"%6.3f",
(htscaleIthreEithJ+htmin)Ihtscale2);

if (dispflag II 1) dispimgC);

break;
} II end of switch ansS II
} II end of if ansz I/
else
if (ans2 II ’I’ ll ans2 II ’1’) dispflagI-dispflag;
}

} /* end of while flag2 */

break;

case ’2’:

closegraph();

printf("\nDisplay points in X axis (currentI Id) = ?",Xdot);
scanf("Xd",&Xdot);
printf("Display points in Y axis (currentI Xd) I ?",Ydot);
scanf("%d",&Ydot);
Xdot2IXdot;
Ydot2=Ydot;
showall();
break;
case ’A’:
if (thabs==’Y’) thabs=’N’;
else thabSI’Y’;
showall();
break;
case ’R’:
for(i=0;i<10;i++) thre[i]=(int)((9-i)I256I9.0+.5);
showall();
break;
case ’H’:

}
}
}

136

closegraph();

printf("\nHax. ht (current I 'I.f) = ? ",htmax);
scanf("1f",&htmax);

printf("\nnin. ht (current I If) I ? ",htmin);
scanf("%f",&htmin);

showall();

break;

case ’ ’:

dispcmd(BLACK);
do {
} while (kbhit() ==0 )3
getch();
dispcdeNHITE);
break;

case ’X’:

changeoolorC);
break;

case ’N’:

closegraph();

printf("Xs :\n%s, %s\n", VERSION, AUTHOR, N-DATE);
loadxt();

resetmotor();

inputdata(&saveall,tdatafile);
scandata(savea11,datafile);

clearcursor();

dispcdeNHITE);

break;

case ESC:

flagIO;
closegraph();
break;

default:

flag=1;
II end of switch() *I
II end of while *I
II end of getcommand() */

display(signed int xp, signed int yp)

{
int
int
int
float
int

x1;
yl;
1;
ht;
ith;

x1=xO-Xdot2pr;
y1=y0-Ydot2Iyp;
setfillstyle(SOLID_FILL,BLACK);

if ((xp II 0) ll (xp II Xmax))

I

bar(x0+3,y0-ly-1,x0+5,y0);

137

II clear cursor in Y axis

bar(x0-lx-5,y0-ly-1,xO-lx-3,y0);

bar(x0-1x,y0+3,x0,y0+4);

II clear cursor in X axis *I

bar(x0-lx,y0-1y-4,x0,yO-ly-3);
BetfillstleCSOLID_FILL,HHITE);

bar(x0+3,y1-Ydot+1,x0+5,y1);

II set cursor in Y axis

bar(x0-1x-5,y1-Ydot+1,xO-lx-3,y1);

}

else

I

bar(x1-Xdot-Xdot+1,y0+3,x1+Xdot,y0+4); II clear cursor in X axis
bar(x1-Xdot-Xdot+1,yO-ly-4,x1+Xdot,yO-ly-3);

}

setfillstyle(SOLID-FILL,NHITE);

bar(x1-Xdot+1,y0+3,x1,y0+4);

II set cursor in X axis

bar(x1-Xdot+1,yO-ly-4,x1,y0-ly-3);

#ifdef DEBUGPLOT
chtExp][yp]I(float)(xp+yp)I(Xmax+Ymax)I3.0-1.5;

#endif
ht-chtExplfyp];

if(thabs == ’Y’) ht=fabs((double)ht);

ithI-i;
do

{
ith++;

} while (( (ht < fthreEithJ) && (ht < fthrefith+1]) ) as (ith < 8) );

setfillstyle(SOLID-FILL,pcolorEith]);
bar(x1-Xdot+1,yi-Ydot+1,x1,y1);
if (kbhit() !I 0)

I
if (getchC)
I

== CTRLC)

closegraph();
fclosea11();
exitC-i);

}
}
}

find1echo(sig,pecho,pechoht)

signed int
unsigned int
float

sigf];
Ipecho;
Ipechoht;

II

I/

}

138

long ytsum;
long ytsummax;
int i,j;
unsigned int echo;
float ht;

for (i=0;i<Length;i++) sigEiJ-Iismean;
echOIO;
htIO.O;
if (Length>inpx-length)
I
ytsummaxIO;
for (iIO;i<Length-inpx_1ength;i++)
I
if (abs(sig[i+pmax_inpx]) > SIGTHRE )
I
ytsumIO;
for (jI0;j<inpx_length;j++) ytsumetsumIsigEin]Iinprj];
if (labs(ytsum) > labs(ytsummax))
I
ytsummax=ytsum;
echOIi+pmax_inpx+Start;
}
}
}
htI(float)ytsummax/xt2sum;
}
else
I
echOIO;
htI0.0;
}
Ipecho=echo;
IpechohtIht;
Cifdef DEBUGECHO
for (iI0;i<Length;i++) sigEiJIIIismean;
Iendif

dispimgC)

I

signed int xp;
signed int yp;
signed int ith;
float ht;
signed int x1;
signed int y1;

139

signed int i;
int x;
int y;

setviewport<w11,w1t,wir,w1b,w1clip);

setcolor(HHITE);

x=xO-1x;

y=w1b-w1t-20;

gprintf(&x,&y,"File . [Zs]",fname);

for (i=0;i<11;i++) fthreEi]=htsca1e*thre[i]+htmin;

for (yp=0;yp<= Ymax;yp++)

I

y1=y0-Ydot2*yp;
setfillstyleCSOLID_FILL,BLACK);
bar(x0+3,yO-ly-1,xO+5,yO);
baerO-lx-5,yO-ly-1,xO-lx-3,y0);

setfillstyle(SOLID_FILL,WHITE);
bar(xO+3,y1-Ydot-1,xO+5,y1);
bar(x0-lx-5,y1-Ydot-1,xO-lx-3,y1);

for (xp-O;xp<-Xmex;xp++)
I
xlaxO-Xdot2txp;
ht=cht [pr [yp] ;
if(thabs II ’Y’) ht-fabs((double)ht);
ithI-1;
do
I
ith++;
} while(((ht<rthre[1th]) u (htdthre [1th+1])) n (ith<8) );

setfillstyle(SOLID_FILL,pcolorfith]);
bar(xl-Xdot+1,y1-Ydot+1,xl,y1);
} /* end of for xp */
} /* end of for yp */
settillstyleISOLID-FILL,BLACK);
bar(x0+3,yO-1y-1,xO+5,yO);
barIxO-lx-5,yO-1y-1,xO-lx-3.y0);

} /* end of dispimg() */
showalII)
I

int x;

int y;

initplot();
dispimgC);

140

dispcmd(HHITE);
} /* end of showall()

takesamp1e(pecho,pht,pgood,pntota1)
unsigned int *pecho;

float *pht;

unsigned int *pgood;

unsigned int *pntotal;

I
signed int i,j;
signed int err;
unsigned int goodsig;
unsigned int ntotal;
unsigned int echo;
float echoht;

if (Nsignal =2 1 )
I
nsample(Nsigna1.Length.Israte,Start,Trigger,
sig,&goodsig,tntotal);

else
I
errBHAXERR*2;
for (i=0;( i<HAXN ) It I err >HAXERR );i++)
I

nsample(Nsignal,Length,Israte,Start,Trigger,

sum,&goodsig,&ntotal);
err=(100*(ntotal-goodsig))lntotal;
}

for (j-0;j<Length;j++) sigEj]-(sum[j]+goodsig/2)/goodsig;

}
findlecho(sig,&echo,lechoht);
*pecho=echo;

*phtsechoht;
*pgood=goodsig;
*pntotal-ntotal;

} /* end of takesample()

clearcursor()

I

int i;

setviewport(w11,w1t,wir,w1b,w1clip);
setfillstyle(SOLID_FILL,BLACK);
bar(xO+3, yO-ly-1,x0+5, yO);

*/

 

 

141
bar(xO-1x-5,yO-1y-1,xO-lx-3,y0)3
bar(xO-lx, y0+3, xO, y0+4);
bar(x0-lx, yO-ly-4,x0, yO-ly-3);
} /* end of c1earcursor() */
dispcmd(int color)
#define CHDI "Lzload stave thhreshold Zzsize A:abs "
#define CHDZ "Hzchg htmax R:reset threshold"
#define CHD3 "Spacezhide command line chhg color "
#define CHD4 "Nzrun again ESC:exit"
I
int x,y;
setviewport(w2l,w2t,w2r,w2b,w2clip);
1-10;
Y810;
setcolor(color);
gprintf(&x,&y,"%s%s",CHD1,CMD2);
gprintf(&x,&y,"%s%s",CHD3,CMD4);
} /* end of dispcde) */
changecolorC)
I

int texth,textw;
int i,jak;

int x,y;
int ith;
int flag2;

char ans2,ans3;

setviewport(w2l,w2t,w2r,w2b,w2c1ip);
clearviewportI);
textw=textwidth("H");
texthctextheight("H");
x=10;
y-10;
setcolor(HHITE);
gprintf(&x,&y," O: 1: 2: 3: 4: 5: 6: 7:");
gprintf(&X,&y,"");
gprintf(&x,&y," 8: 9: A: B: C: D: E: F:");
for (i=0;i<8;i++)
I

k=textw*(i+1)*7;

setfillstyleCSOLID_FILL,i);

bar(k+1,x+1,k+3*textw-1,x+texth-1);

setcolorIBROHN);
rectangleIk,x,k+3*textw,x+texth);

142

k-textw*(i+1)*7;
setfillstyleCSOLID_FILL,1+8);
bar(k+1,x+2*texth+5,k+3*textw-1,x+3*texth+4);

setcolor(BROWN);
rectangle(k,x+2*texth+4,k+3*textw,x+3*texth+4);
}

setviewport(wOl-20,w0t,wOr,wOb,chlip);
setcolor(GREEN);
ith-0;
y=66;
x=1;
gprintf(&x,&y,"=>");
flag2=1;
while(flag2881)
I
do I
} while (kbhitC) =-o ):
ans2=toupper(getch());
ssitch(an82)
I
case ’\r’:
flag2=O;
setcolorIBLACK);
y=66+ith*20;
gprinthtx,&y,"=>");
setviewport(le,w2t,w2r,w2b,w2clip);
clearviewportI);
dispcmd(HHITE);
break;
case ’0’:
case ’1’:
case ’2’:
case ’3’:
case ’4’:
case ’5’:
case ’6’:
case ’7’:
case ’8’:
case ’9’:
pcolorEith]=atoi(&ans2);
k860+ith*20;
setfillsty1e(SOLID-FILL,pcolorEith]);
bar(23,k,72,k+19);
break;
case ’A’:
case ’B’:

C889
C386
case
C880

’C’:
’D’:
’E’:
’F’:

pcolorEith]=ans2-55;
k-60+ith*20;
setfillstyleISOLID-FILL,pcolorEith]);
bar(23,k,72,k+19);
break;

case

’\O’:

an83=getch();
switch(an83)

I

}

case ’P’:
setcolor(BLACK);
y=66+itht20;
gprintf(&x,&y,"=>");
it (ith < 8 ) ith++3
setcolorCGREEN);
y=66+ith*20;
gprintf(&x,&y,"=>");
break;

case ’H’:
setcolor(BLACK);
y=66+ith*20;
gprintf(&x,&y."=>");
it (ith > o ) ith--;
setcolorCGREEN);
y-66+ith*20;
gprintf(&x,&y,"=>");
break;

break;

dispimgC):

143

/* down *I

/* up *I

I! end of switch ans3

/* end of switch ans2
/* end of while f1ag2

/* end of changecolor()

*/
*I

*/

 

 

BIBLIOGRAPHY

BIBLIOGRAPHY

[1] R. Kuc, M. Schwartz, and L. V. Micsky, “Parametric estimation of the acoustic atten-
uation coefficient slope for soft tissue,” in IEEE' Ultrasonics Symposium Proceedings,
pp. 44 — 47, 1976.

[2] K. R. Erikson, F. J. Fry, and J. P. Jones, “Ultrasound in medicine - a review,” IEEE
Transactions on Sonics and Ultrasonics, vol. SU-21, pp. 144 - 164, July 1974.

[3] T. E. Preuss and G. Clark, “Use of time-of-flight c-scanning for assessment of impact
damage in composites,” Composites, vol. 19, pp. 145 — 148, Mar. 1988.

[4] D. Ensminger, Ultrasonics-fnmdamentals, technology, applications. New York: Marcel
Dekker, Inc., 2nd ed., 1988.

[5] F. Dunn, “Attenuation and speed of ultrasound in lung: Dependence upon frequency
and inflation,” Journal of Acoustical Society of America, vol. 80(4), pp. 1248 — 1250,
Oct. 1986.

[6] C. Barnes, J. A. Evans, and T. J. Lewis, “A broad-band single pulse technique for
ultrasound absorption studies of aqueous solutions in the frequency range 200 kHz-1
MHz,” Ultrasonics, vol. 24, pp. 267 - 272, Sept. 1986.

[7] R. W. P. King and J. C. W. Harrison, “The transmission of electromagnetic waves and
pulses into the earth,” Journal of Applied Physics, vol. 39, pp. 4444 — 4452, Aug. 1968.

[8] L. C. Chan, J. L. Peters, and D. L. Moffatt, “Improved performance of a subsurface
radar target identification system through antenna design,” IEEE Transactions on
Antennas and Propagation, vol. AP-29, pp. 307 - 311, Mar. 1981.

[9] H. J. McSkimin, “Pulse superposition method for measuring ultrasonic wave velocities
in solids,” Journal of Acoustical Society of America, vol. 33, pp. 12 — 16, Jan. 1961.

[10] E. P. Papadakis, “Ultrasonic velocity and attenuation: Measurement methods with
scientific and industrial applications,” in Physical Acoustics, Vol. XII (W. P. Mason
and R. N. Thurston, eds.), pp. 277—374, Academic, 1976.

[11] N. P. Cedrone and D. R. Curran, “Electronic pulse method for measuring the velocity of
sound in liquids and solids,” Journal of Acoustical Society of America, vol. 26, pp. 963
— 966, Nov. 1954.

144

[1'2]

[13]

[14]

[15]

[16]

[17]

[18]

[19]

[20]

[‘21]

[22]

[23]

145

D. H. Dameron, “Determination of the acoustic velocity in tissues using an inhomoge-
neous media model,” IEEE Transactions on Sonics and Ultrasonics, vol. SU-26, pp. 69
- 74, Mar. 1979.

M. V. Zummeren, D. Young, C. H. G. Baum, and R. Treleven, “Automatic determina-
tion of ultrasound velocities in planar materials,” Ultrasonics, vol. 25, pp. 188 - 294,
Sept.1987.

J. Toulouse, “A modified version of the phase sensitive technique for measurements of
absolute sound velocity in solids,” in IEEE Ultrasonics Symposium Proceedings, pp. 407
— 410, 1987.

J. P. Weight, “A model for the propagation of short pulses of ultrasound in a. solid,”
Journal of Acoustical Society of America, vol. 81(4), pp. 815 - 826, Apr. 1987.

G. C. Steyer, R. Singh, and D. R. Houser, “Alternative spectral formulations for
acoustic velocity measurement,” Journal of Acoustical Society of America, vol. 81(4),
pp. 1955 — 1961, Apr. 1987.

H. J. McSkimin, “Measurement of ultrasonic wave velocities and elastic moduli for
small solid specimens at high temperatures,” Journal of Acoustical Society of America,
vol. 31, pp. 287 — 295, Mar. 1959.

J. Ophir, “Estimation of the speed of ultrasound propagation in biological tissues: a
beam-tracking method,” IEEE Transactions on Ultrasonics, Ferroelectrics, and Fre-
quency Control, vol. UFFC-33, pp. 359 - 368, July 1986.

T. Kontonassios and J. Ophir, “Variance reduction of speed of sound estimation in
tissues using the beam tracking method,” IEEE Transactions on Ultrasonics, Ferro-
electrics, and Frequency Control, vol. UF F C-34, pp. 524 — 530, Jan. 1987.

J. Ophir and T. Lin, “A calibration-free method for measurement of sound speed
in biological tissue samples,” IEEE Transactions on Ultrasonics, Ferroelectrics, and
Frequency Control, vol. UFFC-35, pp. 573 — 577, Jan. 1988.

J. Ophir, Y. Yazdi, T. S. Lin, and D. P. Shattuck, “Optimization of speed-of-sound
estimation from noisy ultrasonic signals,” IEEE Transactions on Ultrasonics, Ferro-
electrics, and Frequency Control, vol. 36, pp. 16 — 24, Jan. 1989.

R. Kuc, “Estimating acoustic attenuation from reflected ultrasound signals: Com-
parison of spectral-shift and spectral-difference approaches,” IEEE Transactions on
Acoustics, Speech, and Signal Processing, vol. ASSP-32, pp. 1 — 6, Feb. 1984.

R. Kuc, “Bounds on estimating the acoustic attenuation of small tissue regions from
reflected ultrasound,” Proceedings of the IEEE, vol. 73, pp. 1159 - 1168, July 1985.

['24]

[25]

[26]

[27]

[28]

[29]

[30]

[31]

[32]

[33]

[34]

146

R. Kuc, “Estimating reflected ultrasound spectra from quantized signals,” IEEE Trans-
actions on Biomedical Engineering, vol. BME-32, pp. 105 — 112, Feb. 1985.

K. Matsuzawa, N. Inoue, and T. Hasegawa, “A new simple method of ultrasonic veloc-

ity and attenuation measurement in a high absorption liquid,” Journal of Acoustical
Society of America, vol. 81(4), pp. 947 — 951, Apr. 1987.

R. L. Roderick and R. Truell, “The measurement of ultrasonic attenuation in solids
by the pulse technique and some results in steel,” Journal of Applied Physics, vol. 23,
pp. 267 — 279, Feb. 1952.

A. C. Kak and K. A. Dines, “Signal processing of broadband pulsed ultrasound: Mea-
surement of attenuation of soft biological tissues,” IEEE Transactions on Biomedical
Engineering, pp. 321 - 344, July 1978.

Y. Hayakawa, T. Wagai, K. Yosioka, T. Inada, T. Suzuki, H. Yagami, and T. Fu-
jii, “Measurement of ultrasound attenuation coefficient by a multifrequency echo
technique- theory and basic experiments,” IEEE Transactions on Ultrasonics, Fer-
roelectrics, and Frequency Control, vol. UFFC-33, pp. 759 — 764, Nov. 1986.

P. He, “Acoustic attenuation estimation for soft tissue from ultrasound echo enve-
lope peaks,” IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control,
vol. UFFC-36, pp. 197 - 203, Jan. 1989.

G. A. McDaniel, “Ultrasonic attenuation measurements on excised breast carcinoma at
frequencies from 6 to 10 MHz,” in IEEE Ultrasonics Symposium Proceedings, pp. 234
— 236, 1977.

K. J. Parker and R. C. Waag, “Measurement of ultrasonic attenuation within re-
gions selected from b-scan images,” IEEE Transactions on Biomedical Engineering,
vol. BME-30, pp. 431 — 437, Aug. 1983.

W. A. Verhoef, M. J. T. M. Cloostermans, and J. M. Thussen, “Diffraction and dis-
persion effects on the estimation of ultrasound attenuation and velocity in biological
tissues,” IEEE Transactions on Biomedical Engineering, vol. BME—32, pp. 521 — 529,
July 1985.

K. J. Parker, M. S. Asztely, R. M. L. E. A. Schenk, and R. C. Waag, “Frequency
dependent attenuation in normal and diseased livers,” in IEEE Ultrasonics Symposium
Proceedings, pp. 793 — 795, 1986.

E. Walach, A. Shmulewitz, Y. Itzchak, and Z. Heyman, “Local tissue attenuation
images based on pulsed-echo ultrasound scans,” IEEE Transactions on Biomedical
Engineering, vol. 36, pp. 211 — 221, Feb. 1989.

147

[35] N. H. Wang, “Remote sensing by acoustic video pulse techniques,” Master’s thesis,
Michigan State University, Department of Electrical Engineering, 1988.

[36] N. H. Wang, J. Nodar, B. Ho, and R. Zapp, “Video pulse techniques for ultrasonic sens-
ing,” in International Symposium on Ultrasonic Imaging and Tissue Characterization,
1988.

[37] B. Ho, D. Ye, N. H. Wang, J. Nodar, and R. Zapp, “High range resolution ultrasonic
imaging for evaluation of layered composite materials,” in Materials Research Society
Conference, 1988.

[38] N. H. Wang, B. Ho, and R. Zapp, “Attenuation and velocity imagings of biological
tissues by broadband ultrasonic signals,” in International Symposium on Ultrasonic
Imaging and Tissue Characterization, 1990.

[39] N. H. Wang, B. Ho, and R. Zapp, “Nondestructive evaluation of material properties
by broadband signal technique,” in 5th Technical Conference of the American Society
for Composites, 1990.

[40] N. H. Wang, B. Ho, and R. Zapp, “Velocity-density product and attenuation-density
ratio measurements using broadband signals,” IEEE Transactions on Instrumentation
and measurement, Dec. 1991. (accepted for publication).

[41] B. Ho, D. Ye, R. Zapp, and N. H. Wang, “Three-dimensional damage assessment in
composites by ultrasonic imaging techniques,” in 43rd annual conference of Reinforced
Plastics and Composites, 1988.

[42] W. Sachse, B. Castagnede, I. Grabec, K. Y. Kim, and R. L. Weaver, “Recent develop-
ments in quantitative ultrasonic nde of composites,” Ultrasonics, vol. 28, pp. 97 — 104,
Mar. 1990.

[43] J. A. Johnson, B. A. Barna, L. S. Beller, S. C. Taylor, and B. Walter, “A camac-based
ultrasonic dataacquistion workstation,” Materials Evaluation, pp. 934 - 938, Aug.
1987.

[44] S. Lees, “Ultrasonic measurement of thin layers,” IEEE Transactions on Sonics and
Ultrasonics, vol. su-18, pp. 81 - 86, Apr. 1971.

[45] I. Beretsky and G. A. Farrell, “Improvement of ultrasonic imaging and media charac-
terization by frequency domain deconvolution, experimental study with non-biological
models,” Ultrasound in Medicine, vol. 38, pp. 1645 - 1665, 1977.

[46] T. Sato, Y. N akamura, O. Ikeda, and M. Hirama, “Resolution enhancement of a sector-
scan image using a homomorphic transform and deconvolution,” Journal of Acoustical
Society of America, vol. 75(1), pp. 265 - 267, Jan. 1984.

 

148

[47] J. P. Steiner, E. S. Furgason, and W. L. Weeks, “Robust deconvolution of correlation
functions,” in IEEE Ultrasonics Symposium Proceedings, pp. 1031 — 1035, 1987.

[48] A. Yamada, “On-line deconvolution for the high resolution ultrasonic pulse—echo mea-
surement with narrow-band transducer,” in IEEE Ultrasonics Symposium Proceedings,
pp. 1027 — 1030, 1987.

[49] D. F. Elliott, ed., Handbook of digital signal processing, ch. 10. San Diego: Academic
Press, Inc, 1987.

[50] W. Mendenhall, R. L. Scheafl'er, and D. D. Wackerly, Mathematical Statistics with
Applications. Boston: Duxbury Press, 3rd ed., 1986.

[51] E. E. Hundt and E. A. Trautenberg, “Digital processing of ultrasonic data by deconvo-
lution,” IEEE Transactions on Sonics and Ultrasonics, vol. su-27, pp. 249 - 252, Sept.
1980.

[52] T. Yokota and Y. Sato, “Super-resolution ultrasonic imaging by using adaptive focus-
ing,” Journal of Acoustical Society of America, vol. 77(2), pp. 567 — 572, Feb. 1985.

[53] D. P. Shattuck, M. D. Weinshenker, S. W. Smith, and O. T. V. Ramm, “Explososcan:
a aprallel processing technique for high speed ultrasound imaging with linear phased

arrays,” Journal of Acoustical Society of America, vol. 75(4), pp. 1273 — 1282, Apr.
1984.

[54] S. W. Smith, J. H. G. Pavy, and O. T. V. Ramm, “High-speed ultrasound volumetric
imaging system - part I: Transducer design and beam steering,” IEEE Transactions on
Ultrasonics, Ferroelectrics, and Frequency Control, vol. 38, pp. 100 ~— 108, Mar. 1991.

[55] O. T. von Ramm, S. W. Smith, and J. H. G. Pavy, “High-speed ultrasound volumetric
imaging system - part 11: Parallel processing and imaging display,” IEEE Transactions
on Ultrasonics, Ferroelectrics, and Frequency Control, vol. 38, pp. 109 — 115, Mar.
1991.

[56] J. F. Greenleaf and R. C. Bahn, “Clinical imaging with transmissive ultrasonic comput-
erized tomography,” IEEE Transactions on Biomedical Engineering, vol. 28, pp. 177 -
185, 1981.

[57] J. F. Greenleaf, “Computerized tomography with ultrasound,” Proceedings of the
IEEE, vol. 71, pp. 330 — 337, 1983.

[58] H. Ermert and G. Rohrlein, “Ultrasound refelction-mode computerized tomography for
in-vivo imaging of small organs,” in IEEE Ultrasonics Symposium Proceedings, pp. 825
— 828, 1986.

 

149

[59] W. D. Richard, “A new time-gain correction method for standard B-mode ultrasound

imaging,” IEEE Transactions on Medical Imaging, vol. 8, pp. 283 — 285, Sept. 1989.

R. Momenan, M. H. Loew, R. F. Wagner, M. F. Insana, and B. S. Garra, “Applica-
tion of pattern recognition techniques in ultrasound tissue characterization,” in IEEE

Engineering in Medicine 6? Biology Society Conference, pp. 411—412, 1989.

D. Brzakovic, X. M. Luo, and P. Brzakovic, “An approach to automated detection of
tumors in mammograms,” IEEE Transactions on Medical Imaging, vol. 9, pp. 233 —
241, Sept. 1990.

S. M. Lai, X. Li, and W. F. Bischof, “On techniques for detecting circumscribed masses
in mammograms,” IEEE Transactions on Medical Imaging, vol. 8, pp. 377 - 386, Dec.
1989.

M. Hashimoto, P. V. Sankar, and J. Sklansky, “Detecting the edges of lung tumors by
classification techniques,” in IEEE Ultrasonics Symposium Proceedings, pp. 276 - 279,
1982.

M. Bomans, K. Hohne, U. Tiede, and M. Riemer, “3-D segmentation of MR images of
the head for 3-D display,” IEEE Transactions on Medical Imaging, vol. 9, pp. 177 —
183, June 1990.

S. P. Raya, “Low-level segmentation of 3-d magnetic resonance brain images - a rule-
based system,” IEEE Transactions on Medical Imaging, vol. 9, pp. 327 - 337, Sept.
1990.

S. S. Trivedi, G. T. Herman, and J. K. Udupa, “Segmentation into three classes using
gradients,” IEEE Transactions on Medical Imaging, vol. 5, pp. 116 — 119, June 1986.

J. W. K. Jr., C. L. Vaughan, T. D. F. Jr., and L. T. Andrews, “Segmentation of echocar-
diographic images using mathematical morphology,” IEEE Transactions on Biomedical
Engineering, vol. 35, pp. 925 — 934, Nov. 1988.

D. J. Michael and A. C. Nelson, “Handx: a model-based system for automatic seg-
mentation of bones form digital hand radiographs,” IEEE Transactions on Medical
Imaging, vol. 8, pp. 64 — 69, Mar. 1989.

[69] A. Rosenfeld, “Computer vision: a source of models for biological visual processes,”

IEEE Transactions on Biomedical Engineering, vol. BME-36, pp. 93 — 96, Jan. 1989.

B. A. Auld, Acoustic Fields and Waves in Solids, Vol. I. New York: John Wiley &
Sons, Inc., 1973.

[71] V. M. Ristic, Principles of Acoustic Devices. New York: John Wiley & Sons, Inc.,

1983.

150

[72] S. Rama, J. R. Whinnery, and T. V. Duzer, Fields and waves in communication elec-
tronics. New York: John Wiley & Sons, Inc., 2nd ed., 1984.

[73] N. R. Draper and H. Smith, Applied Regression Analysis. New York: John Wiley &
Sons, Inc., 2nd ed., 1981.

[74] L. W. Johnson and R. D. Riess, Numberical analysis. Reading, Massachusetts:
Addison-Wesley Publishing Company, 2nd ed., 1982.

[75] N. H. Wang, J. T. Sheu, and B. Ho, “Tissue characterization by hierarchical clustering
techniques,” in International Symposium on Ultrasonic Imaging and Tissue Charac-
terization, 1991.

[76] A. K. Jain and R. C. Dubes, Algorithms for Clustering Data. New Jersey: Prentice-
Hall, Inc., 1988.

[77] J. C. Birnholz, “An approach to specific tumor diagnosis by ultrasound imaging,” in
IEEE Ultrasonics Symposium Proceedings, 1972.

NSTR TE UNI V.

H]|][|[[|]|]||Hll] 0I] [[1]] [||[ll[[ll]][[l