) A Hand-Held Probe for Vibro-Elastography by Hassan Rivaz B.A.Sc, Sharif University of Technology, Tehran, Iran, 2003 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in THE FACULTY OF GRADUATE STUDIES (Mechanical Engineering) THE UNIVERSITY OF BRITISH COLUMBIA September 2005 © Hassan Rivaz, 2005 Abstract Vibro-elastography is a new medical imaging method that identifies the mechanical properties of tissue by measuring tissue motion in response to a multi-frequency external vibration source. Previous research on vibro-elastography used ultrasound to measure the tissue motion and system identification techniques to identify the tissue properties. This thesis describes a hand-held probe with a combined vibration source and ultrasound transducer which is operational in the 5-25 Hz range to cover a significant bandwidth of the tissue response. The design uses a vibration absorption system to counter-balance the reaction forces from contact with the tissue. Four different types of vibration absorption are briefly described and among them, active dy-namic vibration absorption is selected for the hand-held device. A proportional integrator control is designed for the active vibration absorption and is compared with a well-established active vi-bration absorption method, the delayed resonator. An electromagnetic actuator is selected for active vibration absorption with a single accelerometer providing the feedback. The dynamics of the electromagnetic actuator as well as the response of the different niters are considered in the control law. The stability of the both absorber system and the combined system are considered to find the operational frequency range. The concept of tuning speed for different vibration absorbers is elaborated and is related to the vibration absorption speed. The design of the hand-held device, which includes the active vibration absorber, is presented next. The design utilizes another electromagnetic actuator which is used to vibrate the tissue and is operated with displacement feedback. The design of the vibration absorber is elaborated next. In this design the effect of different absorber parameters on the actuator stroke, control force, stability of the combined system and the tuning speed of the absorber is studied. After ii Ill manufacturing the device, the absorber parameters are identified in order to optimize the vibration absorption performance. Simulation results are provided next that verify the theories presented on the stability and tuning speed. The results show that 100% vibration absorption can be achieved in steady state for both the proportional integrator and the delayed resonator controllers. They also show that nonlinearities in the absorber system decrease the amount of vibration absorption. Experimental results are presented next, showing approximately 85% vibration absorption in the operational range. The existing parameter identification methods that identify the mechanical properties of the tissue using the ultrasound are modified for the hand-held device. Simulation results are provided to validate the revised parameter identification methods. The first elastograms obtained with the hand-held device are finally presented. Contents Abstract ii Table of Contents iv List of Tables viii List of Figures ix Glossary xvi Acknowledgements xix 1 Introduction 1 1.1 Ultrasound Imaging 2 1.2 Static Elastography Methods 4 1.3 Dynamic Elastography Methods, using Doppler ultrasound 6 1.4 Vibro-Elastography 7 1.4.1 Equation of Motion Method with Mass-Spring-Damper Model 7 1.4.2 Transfer Function Method 9 1.4.3 Hand-Held Probe for Vibro-Elastography 11 1.5 Dissertation at a Glance 13 2 Theory of Vibration Absorption 14 2.1 Vibration Reduction Methods 14 2.1.1 Vibration Isolation 14 2.1.2 Passive Dynamic Vibration Absorption 15 iv CONTENTS v 2.1.2.1 Passive DVA with Zero Damping 17 2.1.3 Semi-Active Dynamic Vibration Absorption 18 2.1.4 Active Dynamic Vibration Absorption 18 2.2 Mechanical and Electrical Components 19 2.2.1 Mechanical Components 20 2.2.2 Actuators 20 2.2.3 Sensors 21 2.3 Active DVA Control 21 2.3.1 Delayed Resonator 22 2.3.1.1 DR-DVA with Electromagnetic Actuator 24 2.3.2 PI Vibration Absorber Controller 26 2.3.2.1 PI-DVA with Electromagnetic Actuator 28 2.4 Stability of the Combined System 29 2.4.1 DR Method 29 2.4.2 PI Method 30 2.5 Tuning Speed 31 2.6 Summary 32 3 Design of the Hand-Held Device 34 3.1 Hand-Held Probe Design 34 3.2 Design of the Vibration Absorber Parameters 36 3.2.1 Effect of Absorber Parameters on Control Force and Absorber Mass Displace-ment 37 3.2.2 Effects of Absorber Parameters on the Combined System Stability 42 3.2.2.1 DR Stability 42 3.2.2.2 PI Stability 46 3.2.3 Vibration Absorber Tuning Speed 50 3.2.3.1 PI Tuning Speed 51 3.2.3.2 DR Tuning Speed 52 3.2.4 Designed Values for the ma, ka and ba 54 3.3 Identification of the Absorber Parameters 55 CONTENTS vi 3.4 Summary 57 4 Simulation and Experimental Results 59 4.1 Simulation Results 59 4.1.1 Passive DVA 60 4.1.2 DR-DVA 61 4.1.2.1 Absorber System 61 4.1.2.2 Combined System 61 4.1.3 PI-DVA 61 4.1.3.1 Absorber System 63 4.1.3.2 Combined System 64 4.1.4 The Effect of Tuning Speed on the Vibration Absorption 65 4.1.5 The Effect of Coulomb Friction 67 4.1.5.1 Effect of Excitation Amplitude on the Amount of Vibration Ab-sorption 71 4.2 Experimental Results 72 4.2.1 Absorber System 72 4.2.2 Combined System 73 4.3 Summary 79 5 Vibro-Elastography with the Hand-Held Device 81 5.1 Mass-Spring-Damper Model 81 5.1.1 Equation of Motion from a Stationary Reference 81 5.1.2 Equation of Motion from a Moving Reference 83 5.1.3 Parameter Identification 84 5.1.3.1 Parameter Identification Iterations 89 5.1.4 Simulation Results 90 5.2 Transfer Function Method 93 5.2.1 Calculating Elasticity from Transfer Functions 95 5.2.1.1 Transfer Functions from a Stationary Reference 95 5.2.1.2 Transfer Functions from a Moving Reference 96 5.2.2 Simulation Results 97 CONTENTS vii 5.2.3 Experimental Results 98 5.2.4 Summary 99 6 Conclusions and Future Work 101 6.1 Conclusions 101 6.2 Future Work 102 Bibliography 104 Appendices 109 A Voice Coil Actuator Selection 109 A . l Theory of Voice Coil Actuators 109 A.2 Sizing Voice Coil Actuators 110 List of Tables 3.1 Properties of the electromagnetic actuator stated by the manufacturer. R and L are the DC resistance and the inductance of the coil respectively, Kj is the force sensitivity and Kb is the back emf constant 42 3.2 Absorber system parameters 44 3.3 Primary system parameters 44 3.4 Filters used 48 4.1 Absorber system parameters 66 4.2 Primary system parameters 67 4.3 Amplitude of the hand-held device (primary system) with different vibra-tion absorbers. Vibration of the device with passive DVA, DR-DVA and PI-DVA subjected to the forces depicted in Figure 4.10 68 viii List of Figures 1.1 A typical RF signal. Obtained from a plastic phantom 3 1.2 Ultrasound probe 5 1.3 Static elastography. The depth x inside tissue is related to the time of the RF signal g(t) by x = | where v is the speed of sound 5 1.4 The discrete one dimensional visco-elastic tissue model, (from [63]) 8 1.5 Discrete tissue model in two dimensions.(from [63]) 8 1.6 Transfer function described between two points of tissue, (from [63]) . . . . 9 1.7 A typical transfer function and coherence function. The solid line is transfer function and the dotted line is the coherence function, (from [63]) 10 1.8 A proposed probe for dynamic elastography: (a) ultrasound probe; (b) and (d) audio transducer; (c) water-proof acoustic membrane; (e) vibration damping pad; (f) tissue surface, (from [43]) 12 2.1 Vibration isolation, (a) Vibration isolation to reduce transmitted displacement, (b) Vibration isolation to reduce transmitted force, (c) Vibration isolation in the hand-held device 16 2.2 Passive vibration absorption, (a) Passive dynamic vibration absorber attached to a 1 DOF system (mp, bp and kp represent primary system mass, viscosity and density respectively), (b) Passive vibration absorber attached to a general primary system. (c) Passive vibration absorber in the hand-held device 17 2.3 Passive DVA with zero damping connected to a 1 DOF system 18 2.4 Semi-active vibration absorption 18 2.5 Active vibration absorption attached to the hand-held device 19 ix LIST OF FIGURES x 2.6 Tissue exciter, vibration absorber and signal flow diagram: (a) and (b) amplifier; (c) and (d) signal conditioners; (e) data acquisition board; (f) host PC. . 20 2.7 A typical DR root locus plot when r is fixed, (from [49]) 23 2.8 Gain and frequency vs. delay. 24 2.9 Block diagram representation of the DR-DVA with electromagnetic actu-ator. Kf and Kj, are the force sensitivity and the back emf constant of the actuator, and L and R are the inductance and DC resistance of the actuator. The primary system is not shown in this diagram 25 2.10 Block diagram representation of the PI-DVA with electromagnetic actu-ator. The primary system is not shown in this diagram 28 2.11 A semi-active DVA with zero damping 31 3.1 Exploded view of the hand-held probe: (a) moving head; (b) ultrasound probe; (c) slide rod; (d) machined shell; (e) spacer; (f) linear potentiometer; (g) and (h) primary actuator's coil and magnet; (i) accelerometer; (j) and (k) absorber actuator's coil and magnet; (1) absorber spring 35 3.2 Simplified imaginary model of active DVA. (a) Active dynamic vibration ab-sorber, (b) Equivalent model in steady state condition considering a control force Fc applied according to DR-DVA or PI-DVA control law 39 3.3 Spring-damper model of the tissue 40 3.4 Control force and absorber mass displacement vs. excitation frequency. (a) i 7 ^ as a function of ut. The solid and the dashed curves correspond to kt » btui and the dotted and dashed-dotted curves correspond to btto » kt. Two different absorber system damping coefficients of £ a i =0.1 (the solid and the dotted curves) and (a2 — 0.5 (the dashed and the dashed-dotted curves) are shown. In all curves una = 15 Hz. (b) Minimizing the maximum of Fc by moving ujna to 9 Hz. It is assumed that kt « btco and £a = 0.5. (c) Xa as a function of to. The solid lines and dotted lines represent kt » btcj and btcj » kt respectively. Fo all graphs, Fe = 1 N at w = 5 Hz 41 LIST OF FIGURES xi 3.5 DR stability chart. Gain and frequency versus delay. The plots show the 1st and 2nd branches of the absorber system and the combined system. The second branch plots of the absorber system and the combined system are very similar (the dashed curve and the dashed-dotted curve) 45 3.6 Variation of the maximum and minimum frequency limit of the DR-DVA with absorber system damping 45 3.7 Block diagram representation of the PI-DVA with the lowpass filter added. Duplicated from Figure 2.10 with the addition of the lowpass filter 46 3.8 PI stability chart. Kp and frequency vs. Kp 49 3.9 Maximum frequency limit variation with system properties, (a) copimax as a function of absorber damping ba. (b)wp/m a l as a function of primary system damping bp 49 3.10 Maximum frequency limit variation with lowpass filter cutoff frequency. . 50 3.11 PI tuning speed, (a) Active DVA. (b) Simulation of PI-DVA tuned to resonate at 250 Hz. (c) Simulation of PI-DVA tuned to resonate at 500 Hz. (d) A segment of part (b) from 3 s to 3.1 s. (e) A segment of part (c) from 3 s to 3.1 s 53 3.12 Dominant poles real part vs. tuning frequency. The dynamics of the actuator is not considered in this plot in order to simply show the effect of DR infinite poles and not the actuator pole 54 3.13 Free vibration experiment results 56 3.14 20 Hz forced vibration experiment results, (a) Absorber mass acceleration vs. time.(b) Absorber mass acceleration vs. time, smaller excitation force 57 4.1 Passive vibration absorption simulation results. Absorber mass displacement (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies: (a) 8.77 Hz, (b) 7 Hz, and (c) 11 Hz. The amplitude of excitation force is 1 N for all three simulations 60 4.2 Absorber mass-spring-damper trio free vibration simulation, (a) No control force, (b) DR-DVA tuned to 10 Hz by selecting ra = 0.0157 s and ga = 0.069 kg. (c) Increasing ga to ga = 0.069 x 1.05 kg. (d) Tuning DR-DVA to 6.5 Hz 62 LIST OF FIGURES xii 4.3 Combined system simulation results with the DR-DVA. Absorber mass dis-placement (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies (the DR-DVA tuned to the excitation frequency): (a) 8.77 Hz, (b) 11 Hz, and (c) 6.5 Hz. The amplitude of excitation force is 1 N for all three simulations 62 4.4 Simulation of the PI controlled mass-spring-damper trio vibration, (a) Tuned to 6.5 Hz. (b) Tuned to 2 Hz 63 4.5 PI and DR controlled mass-spring-damper trio vibration simulation, (a) PI. (b) DR. Both plots involve changing the resonance frequency from 15 Hz to 7.5 Hz at t=2 s 64 4.6 Combined system simulation results with the PI-DVA. Absorber mass dis-placement (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies (the PI-DVA tuned to the excitation frequency): (a) 8.77 Hz, (b) 11 Hz, and (c) 6.5 Hz. The amplitude of excitation force is 1 N for all three simulations 64 4.7 Control force of the DR-DVA and PI-DVA when the combined system is excited at 8.77 Hz. (a) DR-DVA. (b) PI-DVA 66 4.8 Stability chart of the DR-DVA. The parameters of the absorber system and the combined system are according to Tables 4.1 and 4.2, and the lowpass filter of Table 3.4 is used. Looking at the lower graph, the operational range of frequency is 3.98-10.28 Hz 67 4.9 Stability chart of the PI-DVA. The parameters of the absorber system and the combined system are according to Tables 4.1 and 4.2, and the lowpass and highpass filters of Table 3.4 are used. The asterics shows the point [Ki Kp]=[0 0] (passive DVA) and indicates the stable part of the chart. Looking at the lower graph, the operational range of frequency is 0-11.19 Hz 68 LIST OF FIGURES xiii 4.10 Comparison of the vibration absorption performances of the DR-DVA and the PI-DVA. (al) Excitation frequency is changed from 4.8 Hz to 4.1 Hz at t=4 s. (a2) and (a3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (al). (bl) Excitation amplitude is changed from 1 N to 1.5 N at t=4 s with a constant frequency of 4.3 Hz. (b2) and (b3) Vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (bl). (cl) Excitation frequency is a constant 6 Hz for t=0-4 s and is decreased linearly from 6 Hz at t=4 s to 4 Hz at t=10 s. (c2) and (c3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (cl). (dl) Excitation amplitude is a constant 1 N for t=0-4 s and is increased linearly from 2 N at t=4 s to 8 N at t=10 s. The frequency of excitation is kept constant at 4.2 Hz.(d2) and (d3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (dl) 69 4.11 Coulomb friction models (a) The discontinuous model, (b) The approximate continuous model 70 4.12 PI controlled vibration absorption with Coulomb friction (a) The value of the Coulomb friction is 0.1 N and the excitation frequency is 5 Hz. (b) The value of the Coulomb friction is 0.2 N and the excitation frequency is 5 Hz. (c) The friction force of the part (b) 70 4.13 The effect of the amplitude of the excitation on absorption achieved with the PI-DVA in the Coulomb friction damped system (a) Amplitude of probe motion with and without vibration absorption, (b) The percentage of the absorption level 71 4.14 "Free vibration" experiments on the absorber system (a) DR controlled mass-spring-damper system tuned to resonate at 10 Hz. (b) PI controlled mass-spring-damper system tuned to resonate at 10 Hz. (c) Unstable DR controlled ab-sorber system caused by choosing a value of ba that is too large 73 4.15 Absorber mass 3 Hz "free vibration". PI controlled mass-spring-damper system tuned to resonate at 3 Hz 73 4.16 The hand-held device, (a) Held by hand, (b) Held by the mechanical arm 74 LIST OF FIGURES xiv 4.17 Vibration absorption experiments at four different frequencies, (al) Vi-bration absorption with the DR-DVA at 8 Hz. (a2) Vibration absorption with the PI-DVA at 8 Hz. (bl) Vibration absorption with the DR-DVA at 10 Hz. (b2) Vi-bration absorption with the DR-DVA at 10 Hz. (cl) Vibration absorption with the DR-DVA at 12 Hz. (c2) Vibration absorption with the PI-DVA at 12 Hz. (dl) Vi-bration absorption with the DR-DVA at 15 Hz. (d2) Vibration absorption with the PI-DVA at 15 Hz 76 4.18 Vibration absorption experimental results with the PI-DVA at 5 Hz. . . . 77 4.19 Absorption performance of the DR-DVA and PI-DVA with variable am-plitude excitation force, (a) DR-DVA. (b) PI-DVA. The excitation amplitude is increased from 1.7 mm to 2 mm at t=1.3 s. Dashed and solid curves show the motion of the probe without and with the active vibration absorber respectively. . . 77 4.20 Absorption performance of the DR-DVA and the PI-DVA with variable frequency excitation force. The excitation frequency is 8 Hz for t=0-2 s and is increased linearly to 16 Hz in 5 s. Dashed and solid curves show the motion of the probe without and with the active vibration absorber respectively, (a) DR-DVA. (b) PI-DVA 78 4.21 Probe vibration for experiments with a human operator on real tissue. The PI controller is used in this experiment and the vibration absorber is turned on at t=0.7 s 78 4.22 Probe vibration for experiments with a human operator on real tissue. The PI controller is used in this experiment and the vibration absorber is turned on at t=0.7 s 79 5.1 The discrete one dimensional visco-elastic tissue model, (from [63]) 82 5.2 Imaging tissue with a moving ultrasound probe 83 5.3 Overall structure of the simulation environment, (from [63]) 91 5.4 Actual coefficients of the MSD model, (a): values, (b): bi values, (c): fcj values 93 LIST OF FIGURES xv 5.5 Tissue displacement and probe displacement, (a): Tissue displacement at the 1s t, 3 r d and 5th nodes of the 6 node MSD model, (b); Displacement assumed for the ultrasound probe 93 5.6 Estimated parameters with 1 mm amplitude xcs- (a), (b) and (c): rrii, bi and ki respectively. No motion for the ultrasound probe is considered, xcs — 0- (d), (e) and (f): rrii, bi and ki respectively at different iterations. The probe is supposed to vibrate as depicted in Figure 5.5 (b) 94 5.7 Convergence variable for 1 mm amplitude xcs- The variable 5m is defined by equation (5.12), and 5b and 5k are defined similarly for b and k 94 5.8 Estimated parameters for 0.2 mm amplitude xcs- (a), (b) and (c): rrii, bi and ki respectively at different iterations. The probe is supposed to vibrate as depicted in Figure 5.5 (b) but with 0.2 mm aplitude 95 5.9 Described procedure's convergence for 0.2 mm amplitude xcs- The variable 5m is defined by equation (5.12), and 5b and 5k are defined similarly for b and k. . . 95 5.10 Estimated elasticity using the transfer function method, (a) Estimated ki values. No motion for the ultrasound probe is considered, (b) Estimated ki values. The probe is assumed to vibrate with a 1 Hz sinusoidal at 1 mm amplitude 97 5.11 B-mode image and static elastograms of a three layered gelatin phantom, (a) B-mode image, (b) Best-case example of static elastogram. (c) Worst-case static elastogram. The geometric locations of the boundaries between layers are indicated with triangles 98 5.12 Vibro-elastograms of the three layered gelatin phantom obtained with the hand-held device mounted on the mechanical arm. (a) Vibro-elastogram at 5 Hz, no vibration absorption is used, (b) Vibro-elastogram at 5 Hz, vibration of the probe cancelled with PI-DVA. The geometric locations of the boundaries between layers are indicated with triangles 99 5.13 Vibro-elastogram of the three layered gelatin phantom obtained with an approximately fixed probe at 5 Hz excitation. The geometric locations of the boundaries between layers are indicated with triangles 100 A . l Voice coil move profile I l l Glossary Elastography Techniques for imaging mechanical properties of tissue. Elastogram An image obtained by elastography. ID One dimensional. 2D Two dimensional. 3D Three dimensional. v Poisson's ratio. E Young's modulus. RF Radio frequency. B-mode image Brightness mode ultrasound image. MR Magnetic resonance. CT Computed tomography. Vibration absorption Cancelling unwanted vibration. DVA Dynamic vibration absorption. DR Delayed resonator. DR-DVA Delayed resonator dynamic vibration absorption. PI Proportional integral. PI-DVA Proportional integral dynamic vibration absorption. PID Proportional, integrator and derivative L-MSD Lumped mass spring damper. xvi GLOSSARY xvii ma Absorber mass. ba Viscosity of the absorber system. ka Stiffness of the absorber system. mp Primary mass. bp Viscosity of the primary system. kp Stiffness of the primary system. Displacement of the absorber mass. Xp Displacement of the primary mass. xa Steady state amplitude of the absorber mass. fc Control force. Fc Steady state amplitude of the control force. Fe Excitation force. ADC Analog to digital converter. DAC Digital to analog converter. LVDT Linear variable differential transformer. 9 Gain term in the delayed resonator. T Delay term in the delayed resonator. L Inductance. R DC resistance. Kb Back emf constant. Kf Force sensitivity. I Electrical current. Resonance frequency of the absorber system. OJp Resonance frequency of the primary system. Kp Proportional coefficient in the proportional integral controller Ki Integral coefficient in the proportional integral controller. GLOSSARY xviii Xt Displacement of the tissue skin. Ca Damping coefficient of the absorber system. CP Damping coefficient of the primary system. Ff Coulomb friction. 6 Parameter vector. <f> Measurement matrix. p Density. ai(t) Acceleration measurement of node i at time t. Vi(t) Velocity measurement of node i at time t. Xi(t) Displacement measurement of node i at time t. CS Coordinate system CS. a.ir(t) Acceleration measurement of node i at time t. Vir(t) Velocity measurement of node i at time t. Xir(t) Displacement measurement of node i at time t. ki Stiffness between nodes i and i + 1. bi Damping between nodes i and i + rrii Mass at node i. Acknowledgements This thesis would not have been completed without the appreciable help of many great people who made my life in UBC an enjoyable and enlightening experience. First of all I would like to express my gratitude to my supervisor, Dr. Robert Rohling. I am greatly indebted to him for his continuous support and patience from the time I entered UBC to the current time. I would like to thank Dr. Tim Salcudean for the idea of using active vibration absorption in the hand-held device. I gratefully acknowledge Simon Bachman for providing the Pro Engineer drawings and helping me with the design of the hand-held device. I would like to thank Dr. Joseph Yan for the invaluable discussions regarding parameter identification. I gratefully acknowledge the support of my colleagues at the robotics and control lab in UBC. I'm specially thankful to Reza Zahiri for helping me with the ultrasound and ultrasound motion tracking. I would like to thank Ehsan Dehghan for the idea of Coulomb friction approximation and Hani Eskandari for helping me with electronics. I would like to thank my wonderful friends at UBC especially Pirooz Darabi, Hamed Mohsenian, Amin Karami, Reza Molavi, Sara Khosravi, Shahrzad Mazlouman, Sara Badiei, Mehran Motamed, Majid Khabbazian, Hosna Jabbari, Ramin Khorasani and Mohsen Abbaspour. I would like to thank my great friends at Green College, who have always cared for me and have made my life here an enjoyable one. Finally, this work is dedicated to my mother and my father's soul and my dear sisters and brothers for all their love and help. They have been my biggest motivation with their support throughout my studies. xix Chapter 1 Introduction Elastography Elastography, the imaging of the elastic modulus1 of tissue, is an emerging imaging method. The motivation of elastography is that soft tissue abnormalities are often correlated to a local change in their mechanical properties. For example, about half of breast cancers detected in US are first dis-covered by the patient feeling a hard lump in her breast [55]. Currently, clinical palpation is widely used to examine breast abnormalities. Unfortunately this technique is subjective and is limited to lesions near the skin surface [34]. To overcome this problem, tissue motion measurement, as the response to an excitation, is proposed as a method to determine the local mechanical properties of tissue. Numerous medical applications of elastography include tumor detection [24], evaluation of vascular health [28], study of skeletal muscle contraction [41], assessment of fetal lung maturity [1], renal transplant rejection and characterization of vascular plaques [12] Different parameters can be extracted by measuring the mechanical properties of tissue. The longitudinal stress is related to the longitudinal strain2 by Young's modulus 3 . Poisson's ratio (v) is the ratio of transverse strain to longitudinal strain, given a stress applied in the longitudinal direction. Density (mass per unit volume) is a mechanical property that is also of interest. Another 1 The terms "elastic modulus", "elasticity" and "stiffness" will be used interchangeably throughout the thesis when referring to tissue. 2Strain (shown by symbol e) is defined as the spatial differentiation of displacement. 3Young's modulus or elastic modulus (shown by symbol E) is defined as the ratio of stress over strain in a homogenous isotropic one dimensionally loaded medium. 1 1.1 Ultrasound Imaging 2 property of tissue, a visco-elastic material, is viscosity. Viscosity is the ratio of the friction shear stress (tangential friction force per unit area) to the velocity gradient perpendicular to the direction of movment. A Young's modulus image is well known to be useful for the diagnosis of cancer tumors in prostate and breast tissue [36]. It is believed that viscosity and density images of tissue can reveal more information about tissue, although the values for normal and abnormal tissue are rare in the literature. A value of 0.5 for the Poisson's ratio means the material is incompressible. Tissue has a poisson's ratio value of between 0.49 and 0.499 which means it is almost incompressible. This small range suggests that it is more difficult to use Poisson's ratio as measurement for diagnosis. The property of tissue elasticity has been the most widely researched of all properties. The basic steps of tissue elasticity imaging can be summarized as follows: 1. Take an image of the region of interest before excitation. 2. Take one image after excitation or some images during excitation. 3. Compare images and investigate tissue response to the excitation at multiple locations. 4. Using tissue response, infer mechanical properties. 5. Create an image where intensity or color is related to the mechanical properties. The excitation can be static or dynamic, leading to static and dynamic elastography. The imaging method is usually Magnetic Resonance (MR), Computed Tomography (CT) or ultrasound. We will focus mainly on ultrasound imaging because of its low cost, speed of image capture and non-ionizing nature. The basic concepts of ultrasound imaging are provided next, followed by a brief description of static and dynamic elastography. 1.1 Ultrasound Imaging Ultrasound is defined as acoustical waves with frequencies from about 20 kHz to several hundred MHz. Medical ultrasound machines typically use a frequency range of 1 MHz to 10 MHz. The waves are generated by small piezoelectric transducers that are placed on the surface of the skin and are triggered electrically. The waves propagate into the tissue where a portion of them is reflected from numerous interfaces between tissue types with different acoustic properties. The reflected 1.1 Ultrasound Imaging 3 waves (echoes) travel back to the piezoelectric crystals and are converted to electrical voltages by the same piezoelectric effect. Ultrasound machines are capable of several modes of operation. We briefly explain A-mode, B-mode and Doppler here. A-mode ultrasound, like most ultrasound modes, uses the pulse-echo technique. In this tech-nique, a short ultrasound pulse is transmitted by a piezoelectric transducer into tissue. The reflec-tion echos are recorded as electrical voltages and then digitalized. The raw unprocessed electrical signal is usually called the radio-frequency (RF) signal. Figure 1.1 shows a typical RF signal obtained from an ultrasound machine. 5001 1 1 1 1 1 1 0 2 4 6 8 10 MM Figure 1.1: A typical RF signal. Obtained from a plastic phantom. The depth of a tissue boundary is calculated by the total transit time from initial pulse trans-mission to reception of the echo, using the propagation speed of sound as the conversion factor. The amplitude of the echo shows the relative difference of acoustical properties at the boundary. Therefore, a single RF signal produces a one dimensional image of the tissue boundaries along the propagation path of the pulse. Usually the amplitude of the RF signal is displayed, hence the letter A in A-mode. B-mode imaging is the natural extension of one dimensional A-mode to two dimensions by using many piezoelectric crystals arranged in a line. These piezoelectric crystals are a part of a hand-held ultrasound probe (Figure 1.2). The various RF signals are combined together in columns. Again the amplitude is measured and converted to brightness of the pixels in the 2D image, hence the 1.2 Static Elastography Methods 4 letter B in B-mode. Doppler ultrasound is mainly used to measure blood flow velocity by measuring the Doppler shift of the ultrasound. The frequency of the wave reflected from a moving body is shifted by an amount proportional to the velocity of the body. As before, the wave is propagated into tissue by the piezoelectric transducers and the echo is received at the piezoelectric transducer, and converted to an RF signal. The frequency shift is calculated to measure the velocity of a local region of tissue, converted to velocity, and displayed as a color overlay. Colour hue and intensity are coded as velocity and direction. 1.2 Static Elastography Methods In static elastography using ultrasound, the tissue is compressed slowly (quasi-statically) with a hand-held probe and ultrasound images are taken before and after compression [51]. The contact force between the probe and skin is controlled by the operator. Suppose that the RF echo signals acquired before and after compression are g(t) and g'{t) respectively (Figure 1.3). Since time in the RF signals of g(t) and g'(t) is related to the depth of tissue through the speed of sound conversion factor, the tissue may be divided into blocks by taking small windows in time of the RF signals. The local displacement of tissue at a point corresponding to time t along the sound beam is estimated from the time lag r obtained by maximizing the cross-correlation of signal g with signal g' over r where w(t) is a window function. The peak location of the correlation coefficient predicts the displacement of that window. This process is repeated for all windows. Therefore, by exploiting correlation techniques between ultrasound RF signals acquired before and after compression (first introduced by [13]) displacement measurements of the tissue are obtained. A displacement, strain or elastic modulus image can be obtained from the displacement vector field obtained by static elastography. A displacement image clearly is not desirable since it is difficult to directly infer the tissue mechanical properties. Therefore a strain image is obtained by calculating the spatial differentiation of the displacement field. The ultimate goal, however, is to oo (1.1) — O O 1.2 Static Elastography Methods 5 Figure 1.2: Ultrasound probe. Figure 1.3: Static elastography. The depth x inside tissue is related to the time of the RF signal g(t) by x = j where v is the speed of sound. produce an image of elastic modulus. Finding elasticity modulus from the displcement or strain field is a challenging inversion problem given the boundary conditions are unknown. So instead, the strain distribution is normally displayed directly. For cases such as small hard nodules in much larger soft tissue, this image can be interpreted well, but more complicated anatomy makes interpretation of strain images more difficult. From an image quality perspective, modulus images reduce the boundary effects, yet they seem to have lower spatial resolution than strain images [25]. Also elasticity calculations from strain data require sophisticated algorithms and there are issues about the uniqueness of the solution [5]. 1.3 Dynamic Elastography Methods, using Doppler ultrasound 6 We have shown that a static elastogram can be obtained by two RF signals, acquired before and after compression. The assumption with cross-correlation based displacement measurement is that the shape of the waveform is relatively unaltered and only the locations are changed. To increase the contrast between hard and soft tissue, the difference in compression level should be increased. On the other hand, as the compression increases, it becomes more difficult to correlate the two RF signals. When strains go beyond 0.01, the uncertainty in determining the cross correlation peak with noisy RF signals increases significantly because the two RF signals appear less similar [18]. To overcome this problem and allow for large compressions, several images are acquired during the compression and motion tracking is performed between every two consecutive RF signals [44]. Local displacements between sequential signals are then summed to calculate the total strains. 1.3 Dynamic Elastography Methods, using Doppler ultrasound In dynamic elastography, tissue is excited with an external vibrator and the excitation penetrates the tissue. In one approach, called sonoelasticity, tissue motion is detected using the Doppler signals [40]. The concept is that the Doppler shift of an ultrasonic wave scattered from a point with sinusoidal motion is given by a Fourier-Bessel series of equally spaced harmonics above and below the Doppler carrier frequency. The amplitude of vibration of the point is also equal to the ratio of the first and second harmonics [39]. If the excitation frequency is swept over a frequency range, amplitudes of vibration at different frequencies are obtained [42]. One problem with this method is that a very stiff region inside tissue has very a small vibration amplitude and therefore is represented by a hole in the sonoelasticity image. On the other hand, the simple absence of the signal also produces a hole, making it difficult to distinguish stiff regions from signal drop-out. Also, this method does not consider any particular dynamic tissue model and therefore only a qualitative estimate of tissue elasticity can be provided. Yamakoshi et al. [66] developed a technique in which both the amplitude and the phase distrib-ution of a low frequency vibration were obtained, aiming at estimating the visco-elastic properties of tissue. In their method, low frequency sinusoidal vibrations under several hundred Hz were applied externally and the resulting movement in tissue was measured from the Doppler frequency shift of the simultaneously transmitted probing ultrasonic waves. In this method, stiff regions are 1.4 Vibro-Elastography 7 also represented by holes because they have very small vibration amplitude also. Therefore it is difficult again to distinguish stiff regions from signal drop-out. 1.4 Vibro-Elastography Vibro-elastography (VE) is a new imaging technique that vibrates tissue externally over a range of low frequency vibrations [64]. What makes VE different from other elastography methods is that tissue response at different frequencies and different time instances are obtained and stiffness is calculated from the combined set of these measurements. The basic principle is to apply techniques from system identification to the problem of tisue property estimation. The VE method also has the potential to extract viscosity and density from the dynamic response of tissue. In VE, ultrasound images are captured while tissue is vibrating and the local displacements are extracted from RF echoes using standard correlation techniques [68]. These displacement fields are then used to determine the visco-elastic properties of tissue. To extract the visco-elastic properties of tissue, two different approaches have been previously developed [63], and described briefly in the following two sections. 1.4.1 Equation of Motion Method with Mass-Spring-Damper Model In the first VE approach to find soft tissue mechanical properties, the tissue is modeled as a set of mass elements inter-connected with dampers 4 and springs [63], as shown in Figure 1.4. This is a one dimensional model that includes both viscous and elastic properties of tissue. It also takes the density variation of tissue into account via the mass values. The two dimensional version of this model is depicted in Figure 1.5. The following assumptions are made in order for this two dimensional model to be valid. 1. Each mass-spring-damper chain is totally decoupled from its neighbors, i.e. tissue properties can only vary in the axial direction. 2. Tissue does not bulge under axial stress, i.e. tissue has zero Poisson's ratio. 4 T h e terms damping and viscosity will be used interchangeably in this thesis when referring to tissue. 1.4 Vibro-Elastography 8 kn . . kn-1 k2 ki . I j — W W W H m n L—VWWVH m n , I , , _|— W A W - ^ | — W W W ^ ^ bn bn-1 b 2 - = E -bi •F(t) Figure 1.4: The discrete one dimensional visco-elastic tissue model, (from [63]) vibration UUlUUUlii Figure 1.5: Discrete tissue model in two dimensions.(from [63]) The equations of motion to this tissue model can be extracted easily, using Newton's laws, Hooke's law and the definition of a damper [63]. The displacement measurements are acquired in real time by a set of successive ultrasound scans, with a typical rate of 100 Hz. These displacement measurements are then differentiated to obtain velocity and acceleration. Then, using parameter identification techniques, elasticity, viscosity and density can be determined. In other words, m,, bi and ki are found such that they fit the equations of motion. It is straightforward to show that the excitation source should contain at least two frequencies in order to extract all rrii, bi and ki [63]. Simulation results show good estimation of stiffness, viscosity and density for homogenous material even when the measurements are noisy. Yet initial results for real inhomogeneous material, like a three layered gelatin phantom with middle layer stiffer that the two surrounding layers, gives less reliable results for viscosity and density. One problem is the poor velocity and acceleration estimations, obtained respectively by differentiation and double differentiation of the displacement values. However, this method gives good results for elasticity even for inhomogeneous materials. The method has an advantage over previous methods because it uses displacement data at many different time instances and solves an over-determined set of equations by the least squares method. 1.4 Vibro-Elastography 9 The least squares fit of the redundant measurements makes the method more robust to measurement noise [63]. 1.4.2 Transfer Function Method A n alternative to the method using the equations of motion is the transfer function method. The transfer function method does not employ the mass-spring-damper model. Here the tissue is again excited externally with multiple frequencies, such as a swept-sine or a band limited white noise input. The transfer function between the measured displacements of two arbitrary points in the tissue represents the dynamics of the local tissue between the two selected points (Figure 1.6). The low frequency values of this transfer function can be considered quasi-static, so it depicts the elasticity aspects of the tissue. Therefore by carefully selecting a small range over low frequencies and calculating the transfer function magnitude over this range, the desired elasticity properties of tissue can be extracted. The higher frequency values provide the viscosity and density properties [63]. vibration X] (input) Xj (output) Figure 1.6: Transfer function described between two points of tissue, (from [63]) One advantage of this method to calculate elasticity is that it performs averaging on the low frequency values of the transfer function, i.e. low frequency values of displacements. This averaging step makes the method more robust to the measurement noise. Another advantage is that the only assumption needed to calculate the transfer function is linearity of the tissue response. The validity of this assumption can be checked by computing the coherence function for each calculated transfer function. Figure 1.7 shows a typical transfer function and coherence function between two points of tissue. A coherence function value of 1 between points Xi and Xj at frequency fo indicates that if Xi has only one frequency component (i.e. pure sinusoidal motion), the Xj would also have only 1.4 Vibro-Elastography 10 one frequency component. Lower values of the coherence function means that applying a pure sinusoidal excitation at point rcj at the frequency of /n, point Xj will move with a spectrum around /n. Perfectly linear systems have a coherence function value of 1. Previous experiments were performed on a phantom of tissue mimicking gel [63]. In those experiments, the phantom was vibrated by a source using a band-limited white noise and the tissue response was imaged with ultrasound at 499 frames per second. It was observed that the coherence function value drops after 30 Hz, indicating that only the transfer function values between 0 and 30 Hz are reliable. However, it should be mentioned that a low value coherence function does not necessarily mean that the tissue is not linear at that frequency. A lack of signal, and therefore a low signal to noise ratio, can also result in a low coherence function. Given a fast excitation, and displacement tracking algorithms that can track motions at higher frequencies, the viscosity and density values should also be identifiable from the higher frequency components of the transfer function. This makes the VE method one of the few elastography methods that has the potential to create viscosity and density images. 1.4, 1.2 0 . 6 0 . 8 0 . 4 0 . 2 r Frequency Figure 1.7: A typical transfer function and coherence function. The solid line is transfer function and the dotted line is the coherence function, (from [63]) 1.4 Vibro-Elastography 11 1.4.3 Hand-Held Probe for Vibro-Elastography Currently, the VE method provides good resolution elastograms from gel phantoms. The transfer functions contain useful information on tissue mechanical properties up to 30 Hz. In [63], the tissue mimicking gels were excited on one side and the ultrasound probe was placed on the other side of tissue. This means having access to both sides of the gel. This is not feasible for most clinical applications on human subjects, where we often have access to only one part of tissue. The work on gels also clamped the ultrasound probe, but clinical ultrasound uses a hand-held probe. This is necessary to keep ultrasound easy to use, convenient and fast. Clinical acceptance of vibro-elastography would be aided by the development of a hand-held probe that includes both the imaging device (the ultrasound probe) and the excitation source. Figure 1.8 shows a probe designed previously for dynamic elastography [43]. Parts labeled (b) and (d) are the heads of two audio transducers connected with vibration damping pads (labelled as (e)) to the ultrasound probe (labelled (a)). Part (c) is a water-proof acoustic membrane which seals the audio transducers from contaminates such as ultrasound coupling fluid (coupling fluid is normally used between the ultrasound probe and tissue skin to translate ultrasound waves more efficiently from the probe to the tissue). Label (f) refers to the tissue surface. The right image shows ultrasound waves propagating into tissue and returning to the probe after reflection from moving tissue. It also shows the mechanical vibration induced by the two audio transducers. When the audio transducers vibrate tissue, reaction forces are exerted from the tissue to the excitation sources and eventually to the probe. On the other hand, the whole device is designed to be light, to reduce operator fatigue. Therefore the light weight hand-held device will vibrate with large amplitudes, sometimes more than the amplitude of tissue excitation. Vibration of the device is undesired because it can produce errors in the VE calculations (described in chapter 5) and increase operator fatigue. Therefore, vibration absorption is needed in the design of such a hand-held device for clinical use [65]. The other problem with the previous design emanates from having two actuators. There are differences in the phase and amplitude of the two audio transducers forces resulting from the differences of signals sent to them and differences in the two transducers' responses. Therefore, the device will experience unbalanced forces, producing a torque to rotate the device. Even small rotations of the probe will result in a change in the imaging plane. This change will produce 1.4 Vibro-Elastography 12 Figure 1.8: A proposed probe for dynamic elastography: (a) ultrasound probe; (b) and (d) audio transducer; (c) water-proof acoustic membrane; (e) vibration damping pad; (f) tissue surface, (from [43]) decorrelation of the RF signals and reduce the accuracy of the displacement measurements, and eventually deteriorate the quality of the elastograms. The main contribution of this thesis is the design of a hand-held device as the next step of VE development. The device should contain a vibration absorber (VA). The goals of the VA are to provide control of the vibration amplitude of the tissue, reduce operator fatigue, ease of scanning (avoid loss of contact between the ultrasound probe and the coupling gel due to probe vibrations) and make easier the parameter identification of tissue properties (described in chapter 5). The following are some desired properties and the aspects of the hand-held probe. 1. The excitation source should contain multiple frequency components. The required frequency range is 5-25 Hz. 2. The amplitude of the excitation force varies, depending on the requested vibration amplitude, the force from the operator on the probe, and the tissue's mechanical properties. 3. The probe can be held by various operators with various inertia, damping and stiffness prop-erties of their hand and arm. This means that the properties of the overall system, whose vibration is to be cancelled, can change substantially from operator to operator. 4. The designed vibration absorber, with all of its components, should be small and light enough to be hand-held. 1.5 Dissertation at a Glance 13 1.5 Dissertation at a Glance In this chapter, a brief review of ultrasound and the different approaches to elastography were provided. The motivations for the design and the requirements of a hand-held device for VE were then presented. The key is to include a method for vibration absorption. In chapter 2, Theory of Vibration Absorption, vibration absorption methods that are potentially applicable to the hand-held device are briefly reviewed. A vibration absorption method is selected and the required components for the selected method are listed. The theory of two feedback control approaches for the chosen vibration absorber is presented next. At the end of the chapter, issues associated with the stability of the combined system are addressed and the tuning speed of a vibration absorber is defined. In chapter 3, Design of the Hand-Held Device, the components of the designed hand-held device are described. The design of the vibration absorber is presented in detail. The stability of the combined system with the designed vibration absorber is examined next. The performance of the vibration absorber is then studied by deriving the settling time of the absorber. Having values for all parameters of the system, the device is constructed. This requires identifying the parameters of the device to optimize the performance of the vibration absorber, which comes at the end of this chapter. Chapter 4, Simulation and Experimental Results, starts with the simulation results that validate the claims made in chapter 2 and chapter 3. Specifically, the claims regarding the set-tling time of the two control methods are verified by the results obtained from the simulation of the absorber system and the combined system. The effect of Coulomb friction on the vibration absorption performance is then investigated by simulation. Finally the experimental results of the absorber system and the combined system are presented. In chapter 5, Vibro-Elastography with the Hand-Held Device, the effect of device vi-bration on the elastograms obtained from real tissue is studied. It is shown that solving for tissue properties, using the equation of motion method with mass-spring-damper model, with the tissue response obtained from vibrating probe requires solving a set of nonlinear equations. The transfer function method using the displacement data obtained from vibrating probe is also investigated. Finally, the first real elastograms obtained with the device from a tissue mimicking phantom are presented, as well as static elastograms obtained from the same phantom. Chapter 2 Theory of Vibration Absorption Many mechanical systems are disturbed by periodic loading and undergo unwanted vibration. Reducing unwanted vibration arises in mechanical and aerospace engineering and civil engineering. This chapter starts with a brief overview of different methods of vibration reduction. According to the requirements of the hand-held device, a vibration reduction method is then chosen. Finally, the theory of the chosen method is presented. 2.1 Vibrat ion Reduction Methods Four different methods for reducing the vibration are well known: vibration isolation, passive dynamic vibration absorption, semi-active dynamic vibration absorption and active dynamic vi-bration absorption. There are numerous variations of these vibration absorption methods. Many of those variations cannot be implemented in a miniature design that must be hand-held. Among all variations, those that are potentially feasible for the VE device are elaborated here. 2.1.1 Vibration Isolation In vibration isolation, a spring-damper combination is placed between the source of the excitation force and the system to be isolated. Figure 2.1 shows a depiction of vibration isolation and how it can be exploited in the hand-held device. In this figure, ki is the stiffness of the isolating material, bi is the viscosity of the isolating material, Fe is the excitation force, Ft is the transmitted force 14 2.1 Vibration Reduction Methods 15 and the primary system is the hand-held probe. To minimize the force transmitted to the primary system, the spring and damper values should be selected according to the excitation frequency [62]. This vibration reduction method is widely used in industry to reduce the transmitted force. As an example, force transmission reduction to the base of mechanical devices with periodic motion, like CNC (Computer Numerically Controlled) milling machines, can be achieved by using visco-elastic materials in their foundation. Another well established application of vibration isolation is to reduce transmitted displacement. The suspension system of automobiles, which reduces the transmitted displacement from the wheels to the chassis, is another example. These vibration isolators can be designed for only a narrow band of excitation frequencies. Also, they only reduce the transmitted vibration and, even theoretically, they never can be designed to isolate the primary device from vibration perfectly. Tunable vibration isolators were introduced recently [16]. These systems are able to work over a wider frequency range. Tuning can be achieved by using either variable dampers or variable springs. Electro-rheological fluid materials, which undergo significant instantaneous changes in mechanical properties by applying voltage to them, are employed in variable dampers [4], [53]. Changing the number of effective coils in the absorber spring can be used to tune the system to absorb a wider frequency range of excitations [22]. Many recent applications of tunable vibration isolators are found in the automobile industry. Examples are variable stiffness suspension [67], electro-rheological fluid damping suspension [15], electro-hydro-pneumatic suspension [57], tunable gas suspension [29] and tunable suspension for passenger trains [60]. In order to control variable dampers and springs, different classical control methods, like optimal control [19], have been used. 2.1.2 Passive Dynamic Vibration Absorption For passive dynamic vibration absorption (DVA), an additional mass-spring-damper trio is attached to the primary system (Figure 2.2). In this figure, ma, ka and ba are the absorber mass, stiffness and viscosity respectively and xa and xp are the position of the absorber mass and primary mass. There are numerous types of passive dynamic vibration absorbers. Reference [35] has a fairly complete list of these methods. Unlike the previous method, passive vibration absorption adds one degree of 2.1 Vibration Reduction Methods 16 (a) isolated device y(t)=Ytsin(wt) moving base ; ^ j | isolating material 21 Ty(t)=Ysin(wt) (b) vibration source fixed base f(t)=Fsin(wt) (c) Primary sys. F t = Frsin(e)t) isolating material Actuator I • • i | 1 S 0 l a t i n g m a t e n a l ^ F o s i n O o t ) Jf(t)=F tsin(wt) j~ Tissue "| Figure 2.1: Vibration isolation, (a) Vibration isolation to reduce transmitted displacement, (b) Vibration isolation to reduce transmitted force. (c) Vibration isolation in the hand-held device. freedom (DOF) to the system. It can be shown that steady state amplitude values for xp and xa are [30] x _ [(fca - mau)2) + bau>j] Fe p ~ det{K - u 2 M + ujjB) Xa (ka + baujj)Fe (2-1) (2.2) det(K - OJ2M + UJJB) where det(...) refers to the determinant function and M, B and K are the system coefficient matrices: M mp 0 0 m„ bp + ba -ba -ba ba K = kp -f- ka ka kn. kn_ (2.3) Xp and Xa are complex values meaning that they have two terms, one in phase with the applied force and one out of phase. To minimize the vibration of the primary system (xp), the absorber's mass, spring and damper values should be selected according to the frequency range of vibration of the primary system [30]. Yet these absorbers can be designed for only a narrow band of excitation frequencies . There is a trade-off between the frequency bandwidth and the absorption level: the wider the bandwidth (or simply bigger the ba), the lower the absorption. Another tradeoff exists between the absorber mass and its displacement [61]. 2.1 Vibration Reduction Methods 17 (a) t—. Fesin(wt) I ] l x , t ) (0 Absorber sys. Primary sys. Absorber svs S i Primaiy sys. F e= Fosin(o)t) L J Figure 2.2: Passive vibration absorption, (a) Passive dynamic vibration absorber attached to a 1 DOF system (rap, bp and kp represent primary system mass, viscosity and density respectively), (b) Passive vibration absorber attached to a general primary system. (c) Passive vibration absorber in the hand-held device 2.1.2.1 Passive DVA with Zero Damping Passive DVA with zero damping has several interesting properties that gives intuition to many vibration absorption methods. Figure 2.3 shows a passive DVA with zero damping attached to a 1 DOF system. The equation of motion of the mass ma is ka(xp — xa) = maxa or equivalently ka%p = maxa + kaxa. Taking the Laplace transform, the transfer function between xa and xp can be obtained: Xp(s) = mas2 + ka Xa(s) ka • [ • > where Xp(s) and Xa(s) are the Laplace domain representation of xp and xa. Suppose that the excitation force Fe has the frequency LU. The absorber system is said to be tuned to the excitation frequency if ka = maw2. Equation 2.4 indicates that the steady state value of xp is zero provided that the absorber system is tuned to the excitation frequency. This means that a resonator (undamped mass-spring system) can cancel perfectly the vibration of the primary system if it is tuned to the frequency of excitation. A passive absorber introduces other peak response frequencies in the system, typically one below and one above the vibration frequency it suppresses. Therefore the frequency of the excitation should be within a very small frequency range that the passsive DVA is designed for. In the hand-held device, it is required to excite the tissue in the 5-25 Hz frequency range. Therefore a passive DVA is not suitable. 2.1 Vibration Reduction Methods 18 ]lx(t, Absorber sys. Fesin(wt) 1 Primary sys. Figure 2.3: Passive DVA with zero damping connected to a 1 DOF system. The second problem is that in order for these vibration absorbers to be able to absorb the vibrations perfectly, they should exhibit resonance behavior which can only be achieved with zero damping. Therefore, even for single frequency excitations these absorbers cannot suppress the vibration perfectly because every physical system has some degree of damping. 2.1.3 Semi-Active Dynamic Vibration Absorption Semi-active vibration absorption follows the same idea as passive absorption, but uses a tunable spring and damper (Figure 2.4). Varying the spring and damper values allows the absorber to be tuned to different excitation frequencies. Jalili [32] provides a survey of different semi-active vibration absorption systems . Semi-active vibration absorbers still have the problem of not being able to suppress the vibrations perfectly because of inherent damping. T Tl.->IIC: Figure 2.4: Semi-active vibration absorption. 2.1.4 Active Dynamic Vibration Absorption In active DVA, an actuator is added to the mass-spring-damper trio of a passive dynamic absorber (Figure 2.5). Active DVA differs from semi-active DVA in that it needs more control effort to sup-2.2 Mechanical and Electrical Components 19 press the vibration, requiring an external energy source. Its advantage over semi-active methods is that it can cancel unwanted vibrations perfectly provided that the excitation force has some char-acteristics, described later in Section 2.3. Also, for swept-sine excitations, it can cancel vibrations much faster, described later in Section 2.5. We choose this method for the hand-held probe. Figure 2.5: Active vibration absorption attached to the hand-held device. 2 .2 M e c h a n i c a l a n d E l e c t r i c a l C o m p o n e n t s Given active DVA, the hand-held probe and the vibration absorber can be designed as depicted in Figure 2.6. The moving head is the part that vibrates the tissue and is driven by force Fe, applied by the first actuator. To provide a smooth vibration, this actuator works in a closed loop control system, with the feedback fa. This feedback should provide the relative position of the moving head and the hand-held device to the controller since position control is required. The feedback fa is amplified by the signal conditioner (part (c)) and then converted to a digital signal by the ADC (Analog to Digital Converter). The digital signal is then fed to the host PC (labelled (f)). A conventional PID (Proportional Integrator Differentiator) controller is used here. The digital control signal of the host PC is then convereted to an analog signal by the DAC (Digital to Analog Converter, part (e)) and drives the first actuator after being amplified at the amplifier (a). The active DVA is controlled by force Fc, applied by the second actuator, fa is the feedback that is used to drive the second actuator. A similar signal flow as the first actuator can be used to drive the second actuator. Primaiy sys. T Ti—l ie 2.2 Mechanical and Electrical Components 20 Figure 2.6: Tissue exciter, vibration absorber and signal flow diagram: (a) and (b) ampli-fier; (c) and (d) signal conditioners; (e) data acquisition board; (f) host PC. 2.2.1 Mechanical Components Regular coil springs are chosen for the absorber system. The required damping value depends on the control method. It will be shown later that an approximate linear damping value of 3 kg/s is beneficial to the vibration absorption performance, but not necessary. A suitable miniature damper (dashpot) with such a small damping value is impractical. Other types of dissipation (like Coulomb friction or Eddy currents) provide damping but their nonlinearities affect performance. 2.2.2 Actuators As depicted in Figure 2.6, two actuators are required: one to excite the tissue and the other for vibration absorption. Either piezoelectric actuators or electromagnetic actuators are suitable here since they offer fast response. The multiplication of the absorber mass vibration amplitude and the absorber mass value is a design parameter (shown later in equation (3.2)). Therefore a piezoelectric actuator's small stroke requires the absorber mass to be very heavy (for a 2 mm stroke piezoelectric actuator the absorber mass should be approximately 4 kg). Therefore an electromagnetic actuator is selected for the absorber actuator. The disadvantages of an electromagnetic actuator compared to a piezoelectric actuator are weight and response speed. The magnet comprises 90% of the weight of the electromagnetic actuator 2.3 Active DVA Control 21 selected for the device. By mounting the coil to the probe, the magnet can be used as the absorber mass of the active DVA. It is shown in Section 3.2.3.1 that the electromagnetic actuator used in the device is fast enough for the 5-25 Hz range of frequency. A similar voice coil is also chosen for the tissue exciter for convenience. 2.2.3 Sensors It will be shown later that the vibration absorption control methods used here only require one feedback from the absorber mass. Therefore only two sensors are required in the device: one displacement sensor to control the tissue exciter actuator and one acceleration (or displacement) sensor to control the vibration absorber actuator. Possible choices for sensors include an LVDT (linear variable differential transformer), capaci-tive displacement sensor, Eddy current position sensor, magneto resistive, optical sensor (like laser triangulation sensor), linear encoder and potentiometer. None of the sensors that rely on a mag-netic field can be used because of the very strong magnetic fields produced by the two actuators. Considering size and price, linear encoder and optical sensors are eliminated. A potentiometer is chosen for the tissue excitation where large forces axe present. It cannot be used to measure the displacement of the absorber mass because the friction force is a significant problem. Therefore an accelerometer is used to provide the feedback for vibration absorption. 2.3 Active D V A Control Active DVA control methods span a wide range of feedback control laws including linear quadratic regulator [6], H ^ [46], neural networks [38] and Delayed Resonator [49]. Using these modern control technologies, active vibration absorbers are finding their way in new applications like flexible space structures [8], structural control [58], [59] and helicopter vibrations [2]. DVA design usually comprises an analysis of the full system vibration, including the DVA. Examples include bandpass vibration absorber [21], delayed feedback vibration absorber [33], LQG active vibration control [10], optimal dynamic absorber [3], and dynamic absorber for general multi degree of freedom systems [69]. Such design approaches cannot be used for the hand-held 2.3 Active DVA Control 22 device because the probe can be held by various operators with various inertia, damping and stiffness properties. This means that the properties of the system, whose vibration is to be cancelled, can change substantially. A stand-alone active DVA should be used in the hand-held device. That is, the vibration absorption control law is decoupled from the primary system and does not require any information from it. Two such implementations of active DVA are the Delayed Resonator (DR), first introduced by Olgac and Holm-Hansen in 1994 [49], and the Proportional Integrator (PI) controller, developed in this thesis. In order to assist the reader, we present a brief summary of the DR active DVA (DR-DVA). We then also introduce the PI controller as an alternative stand alone DVA controller (PI-DVA). Both of these methods can cancel vibrations perfectly provided that the excitation frequency has a single frequency component1. Therefore, we can cover the 5-25 Hz required frequency range for VE by exciting the tissue by a swept-sine signal whose frequency increases linearly (or nonlinearly) from 5 Hz to 25 Hz. An example of acquisition time in VE is 10 s, resulting in a 2 Hz/s sweep rate. 2.3.1 Delayed Resonator The control force, Fc, (Figure 2.5) in the DR active DVA (DR-DVA) with acceleration of the absorber mass as feedback is Fc = gxa(t - T) (2.5) where g is the feedback gain and r is a time delay. Having this force as the control command, the equation of motion of the absorber mass is [49] maxa + baxa + kaxa + gxa(t - r) = 0 (2.6) The Laplace domain representation leads to the following transcendental characteristic equation mas2 + bas + ka + gs2e~TS = 0 (2.7) This equation has infinitely many complex roots for s when g 0 and r ^ 0. Figure 2.7 depicts 1DR-DVA can be used to cancel vibrations with two frequency components [50] also. However the two frequencies are fixed and are determined by the parameters of the absorber system. 2.3 Active DVA Control 23 the root locus plot for this equation. To achieve resonator behavior at the frequency u>a, two roots of the equation (2.7) should be placed at ±jcoa while the other roots remain in the left-half of the complex plane. Subscript a refers to the absorber's "natural" frequency. Solving for g and r one has [49] 9a = -2\/(Va) 2 + {maujl - ka)2 (2.8) r a = — {atan2(bau>a, mauia - ka) + 2(1 - l)tr} ,1 = 1,2,- (2.9) The variable parameter I refers to the branch of root loci that the roots s = ±juia axe located on. The importance of this parameter is visualized in the ga(^a) vs. Ta(uja) plot. Figure 2.8 shows gain (ga, (kg)) and frequency (wa, (Hz)) vs. delay ( T 0 ) (S)) when the absorber parameters are selected to be ma = 0.2 kg, ba = 3.6 kg/s and ka = 637 N/m. These values are the actual values of the absorber system paxameters and they will be obtained in chapter 3. The first and the second branch cross approximately at gc = 0.15 kg and r c = 0.058 s (top plot). Subscript c refers to the critical gain and delay value. This point corresponds approximately to wc = 7 Hz (bottom plot). From the bottom plot, it can be inferred that to tune DR to a frequency less than uc using its first branch, the delay should be bigger than r c, say r+. From the top plot, suppose that for the delay value T+ the gain value required to put the roots of the first branch (the solid line) at ±jtuc is g+. On the other hand, at the same delay value of r+, the gain value that is required to put the roots of Figure 2.7: A typical DR root locus plot when r is fixed, (from [49]) 2.3 Active DVA Control 24 the second branch (the dotted line) on the imaginary axis is some value less than g+, say g2nd branch (because for T > r c, the dotted line is below the solid line). This means that for the delay value T + and gain value g^nd branch the two roots of the system are on the imaginary axis (the two roots on the second branch) and the two roots on the first branch have not yet reached imaginary axis. By increasing g, trying to put the two roots on the first branch on the imaginary axis, the two roots on the second branch move to the RHP (right half plane) of the imaginary plane, making the system unstable. Therefore, for I = 1 (i.e. the first branch of root loci crossing the imaginary axis at the ±ju>a), the delayed resonator can be tuned to any frequency between 7 Hz and oo. For I > 2, following similar reasoning, the frequency range of operation will be limited by both upper and lower bounds due to previous and next branch crossings respectively. It can be shown that the minimum frequency limit for the first branch is less than the lower bound of any other branch [17]. Therefore, the lower operational frequency limit of DR-DVA can be obtained by the crossing of the first and second branch (on a ga vs. r a plot) while there is no upper frequency limit. 2.3.1.1 DR-DVA with Electromagnetic Actuator The block diagram representation of a DR controlled mass-sprimg damper with an electromagnetic actuator is shown in Figure 2.9. With the electromagnetic actuator, the control command is the applied voltage. Therefore, the dynamics of the electromagnetic actuator must be considered. Modifying the formulation of [17] to fit DR-DVA with an electromagnetic actuator and acceleration feedback, the system is modeled as: Figure 2.8: Gain and frequency vs. delay. 2.3 Active DVA Control 25 The DR controller <H1 gain Delay accelerometer data 1 Ls*R Figure 2.9: Block diagram representation of the DR-DVA with electromagnetic actua-tor. Kf and Kb are the force sensitivity and the back emf constant of the actuator, and L and R are the inductance and DC resistance of the actuator. The primary system is not shown in this diagram. maxa(t) + baxa(t) + ka(t) — fc(t) LI(t) + RI(t) + Kbxa(t) = ec{t) fc(t) = KfI(t) ec(i) = -gx(t - T) (2.10) (2.11) (2.12) (2.13) where L and R are the inductance and the DC resistance and Kf and Kb are the force sensitivity and back emf constants of the actuator. Taking the Laplace transform and forcing s = ±ju>a, the gain and delay values can be calculated. Reference [17] provides the gain and delay values for the DR-DVA with electromagnetic actuator and displacement feedback. Modifying this to DR with the acceleration feedback, we obtain: 9a (Ls + R)(mas2 + bas + ka) + KfKbs Kfs2 S=]bJa ra = ±{(2l- IK ~ Z { l S + R ) { m a S 2 +tS + K ) + KfKbS } ,1 = 1,2,-w a I Kf I (2.14) (2.15) In contrast with the previous section where there was no upper bound for the tuning frequency of DR-DVA, equation (2.15) shows that there is an upper bound for tuning frequency. Equation 2.3 Active DVA Control 26 2.15 indicates that ra is a monotonously decreasing function of tuning frequency, wa. Therefore zero delay corresponds to the highest frequency that DR-DVA can be tuned to. Solving equation (2.15) for tja letting r a = 0, one has To obtain the minimum operation frequency, ga vs. ra can be plotted similar to the procedure that was used to obtain the minimum operation frequency without considering the dynamics of the electromagnetic actuator (Section 2.3.1). Therefore, the DR-DVA using an electromagnetic actuator has both lower and upper limits. Considering equation (2.16), the upper frequency bound clearly increases by increasing damping and decreasing absorber mass (although decreasing absorber mass will eventually decrease the operation frequency by impacting the combined system stability). The main advantages of this method are the minimum number of feedbacks required (one sensor) and being decoupled from the primary system (i.e. no information of the primary system is required to determine the gain and delay control parameters). 2.3.2 PI Vibration Absorber Controller A PI controlled active DVA (PI-DVA) with acceleration feedback is now developed for an elec-tromagnetic actuator. For an intuitive view of the final controller (which considers the dynamics of the electromagnetic actuator and the filters used) we start with a simple model neglecting ac-tuator dynamics and assuming smooth acceleration data with zero drift. Like DR-DVA, we then investigate the electromagnetic actuator in the analysis. As mentioned in Section 2.1.2.1, a passive DVA has two limitations: first it can only cancel the vibration at a very small range of frequency, and second the damping exists in every physical system and prevents achieving perfect vibration absorption. The control force in an active DVA should obviate the two limitations of the passive DVA: first by acting such that the natural frequency of the absorber system changes with the excitation frequency, allowing for vibration absorption in a much wider range of frequency, and second by cancelling the damping force, allowing for perfect (2.16) 2.3 Active DVA Control 27 vibration absorption. The control force, Fc (Figure 2.5), of such a PI controller has the form Fc = Kpxa + Ki±a (2.17) where Kp and Ki are the proportional and integral gain coefficients respectively and xa is the position of the absorber mass. The first control parameter, Kp, will be used to allow for changes in the natural frequency of the absorber system, and the second control parameter, Ki, will be used to cancel the damping force. Having this force as the control command, the equation of motion of the absorber mass is maXa + baXa + kaXa = KpXa + Ki±a (2.18) Taking the Laplace transform one has mas2 + bas + k = Kps2 + KiS (2.19) By selecting the Kp and Ki such that the roots of the above equation are placed on the imaginary axis at ±ju>, the absorber system will mimic a resonator at the frequency of u>a. This system can cancel the vibrations of the primary system at the frequency uia following the discusion in Section 2.1.2.1. Solving equation (2.19) with s = ±juJa, one has Kp = ma - ka/Jl (2.20) Ki = ba (2.21) The advantage of this controller over DR is that it doesn't have any lower bound on the operation frequency. The upper frequency bound, like DR, is oo. Another advantage of this method over DR is the absence of other roots. In DR, there are infinitely many roots and the closest root to the imaginary axis in the left-half of the complex plane determines the speed at which DR cancels the vibration. These roots need to be considered in the design process. 2.3 Active DVA Control 28 2.3.2.1 PI-DVA with Electromagnetic Actuator The block diagram representation of the absorber system with the electromagnetic actuator is provided in Figure 2.10. In this figure, L and R are the inductance and the DC resistance and Kf and Kb are the force sensitivity and back emf constants of the actuator. To obtain the velocity data from the acceleration data, a first order Butterworth highpass filter with a 0.5 Hz cutoff frequency is used before integration to eliminate the low frequency drift of the accelerometer data. Again, to have a vibration absorber at the frequency uia, Kp and Ki should be selected such that the roots of the closed loop system are at ±jcoa- By finding the closed loop transfer function of the whole system and forcing it to have poles at ±ju>a, the integrator and proportional gains can be found L(ka - macol) + Rba + KfKb K f Kp = _ R{K/ul - ma) - Lba ujhp Kf Ki (2.22) (2.23) where ujhp represents the highpass filter cutoff frequency (0.5 x 2 x )• The PI Controller Ki lb utter Integrate ' Highpass Filter 1* accelerometer data 1 L .s*R Figure 2.10: Block diagram representation of the PI-DVA with electromagnetic actua-tor. The primary system is not shown in this diagram. 2.4 Stability of the Combined System 29 2.4 Stability of the Combined System Although it is guaranteed that both DR and PI controllers make the absorber system marginally stable, the stability of the combined system should be investigated. In this section, formulas required to examine the stability of the combined system will be presented. We will show in the following chapters that this range agrees with the simulation and experimental results. 2.4.1 D R Method The dynamic equations of the combined system (Figure 2.6) can be written as mpxp(t) + (bp + ba)xp(t) + (kp + ka)xp(t) - baxa(t) - kaxa(t) = fe - KfI(t) (2.24) maxa(t) + baxa(t) + kaxa(i) - baxa(t) - kaxp(t) = KfI(t) (2.25) LI(t) + RI{t) + kbxa(t) - Kbx(t) = gxa(t - T ) (2.26) Taking the Laplace transform, one obtains: G(s) mpS 2 + (bp + ba)s + (kp + ka) -(bas + ka) Kf X Fe -(bas + ka) mas 2 + bas + ka -Kf Xa = 0 —Kbs Kbs — gs2e~TS Ls + R I 0 (2.27) The characteristic equation of the system is det(G(s)) = 0. Doing so, one has: N(s) + D(s)ge-TS = 0 where N(s) = (mps2 + bps + kp) {(mas2 + bas + ka)(Ls + R) + KfKbs} +mas2 {(Ls + R)(bas + ka) + KfKbs}) (2.28) (2.29) 2.4 Stability of the Combined System 30 D(s) = Kcs2(mps2 + bps + kp) (2.30) The combined system characteristic equation is transcendental, like the absorber system charac-teristic equation. Increasing gain for such a system ultimately causes instability as the roots move to the right half plane. The combined system becomes marginally stable when the first two roots from the "root cloud" cross the imaginary axis, say at s = ±jwp (subscript p refers to the primary system). Introducing this condition into equation equation (2.28), one obtains the gain and delay values that make the combined system marginally stable: 9P N(s) T p up { (21 - 1 ) T T - Z D(s) N(s) s=±juip I = 1,2, (2.31) (2.32) 'D(s)ja=±ju equation (2.31) and equation (2.32) will be used later to determined the stability of the combined system. 2.4.2 PI Method The dynamic equations of the combined system are very similar to equation (2.24), equation (2.25) and equation (2.26) from the DR-DVA. The only difference is that the control force of the DR-DVA (gxa(t — T)) should be replaced with the control force of the PI-DVA (Kpxa + KiXa). Repeating the same procedure as the previous section, the characteristic equation is: N(s)+D(s)(Kps2 + KiS) =0 (2.33) where N(s) is given by equation (2.29) and D(s) is: D(s) = Kf(mps2 + bps + kp) (2.34) Again, the marginally stable point of the combined system is an important point in studying the system's stability. Finding Kp and Ki from equation (2.33) gives: *>=-4MfL„ <2-35) 2.5 Tuning Speed 31 Ki = I mag jN(s)\ (2.36) \D(s)j We refer to these equations later when we want to investigate the stability of the combined system. Tuning speed of an active DVA or a semi-active DVA can be defined as the speed that the resonance frequency of the mass-spring-damper system of the DVA can be changed. Consider the passive DVA of Figure 2.3. Assume that the absorber spring ka can be changed simultaneously (with no delay) with the excitation frequency (Figure 2.11). Regardless of the speed of the variation of the excitation frequency, the absorber system tunes itself to the excitation frequency. Therefore, the tuning speed of this system is theoritically infinite provided that the variable spring ka can be changed with an infinitely large speed. The low speed that the variable spring can be changed in the semi-active DVA (because of the mechanical difficulties of the mechanisms and the actuators that change the stiffness value) makes these absorbers not suitable for excitations with variable frequency since they cannot be tuned fast enough to the excitation frequency. However this kind of DVA is suitable for excitations with variable amplitude since the properties of the absorber system should not be changed when the amplitude of excitation changes. Figure 2.7 (or equation (2.7)) indicates that there are infinitely many roots in the characteristic equation of a DR-DVA other than the two roots on the imaginary axis. These many roots can be neglected in the steady state since they all have negative real values, but they determine the transient response of the DR-DVA. A DR-DVA system tuned to resonate at the frequency w is similar to the mass-spring resonator of Figure 2.3 only after a settling time. Therefore, usuaslly 2.5 Tuning Speed Primary sys. Absorber sys. Figure 2.11: A semi-active DVA with zero damping. 2.6 Summary 32 a longer time is required for a DR-DVA (depending on the place of its extra poles) to cancel the vibrations compared to a PI-DVA. It should be clear that the theoretical steady state vibration absorption amount for all resonant DVA systems (including the passive DVA with zero damping, the DR-DVA and the PI-DVA) is 100% according to equation (2.4). However, the time required to reach the steady state is different for different methods. We explore the tuning speed of the DR-DVA and PI-DVA in Section 3.2.3. It is worth noting here that, in general, the tuning speed is not only important when dealing with the excitations with variable frequency, but also for excitations with variable amplitude. For example, imagine a DR-DVA attached to a primary system that is subjected to an excitation force with continuously variable amplitude (or frequency). This system never reaches the steady state and therefore the extra roots of the DR-DVA deteriorate the vibration absorption performance. We also mention in practice, excitations with constant frequency and constant amplitude can contain other frequency components with small power compared to the power of the excitation at its main frequency. No active DVA that makes the absorber system a resonator (including DR-DVA and PI-DVA) is capable of cancelling vibration caused by these other frequency components. How-ever, these other frequency components deteriorate the performance of the DR-DVA by preventing the system from reaching stady state. Therefore the vibration absorption of the DR-DVA at the main frequency will be affected by the power of these other frequency components. 2.6 Summary-Variations of four different vibration absorption methods were studied in this chapter. Among those vibration absorption methods, the active DVA was selected and two control laws for this vibration absorber, DR-DVA and PI-DVA, were elaborated. The required electrical and mechanical components for an active DVA were described in this chapter. The required formulations for stability analysis and some intuition about the tuning frequency of active and semi-active DVA were provided at the end. Using PI-DVA with a single acceleration feedback is a novel approach to vibration absorption. All PI formulations, including the combined system stability analysis, are new contributions. The DR formulation is also modified, as mentioned before, to fit the specific type of sensor and actuator 2.6 Summary 33 that we use: acceleration feedback and electromagnetic actuator. The analysis of tuning speed is also novel. Chapter 3 Design o f the Hand-Held Device The design of the hand-held device is described in this chapter. The absorber system parameters are first selected according to design issues. A stability analysis is then performed to investigate the operational range of the device. The performance of the absorber is studied next. Having the final values of the absorber parameters, the device is assembled and parameter identification of the absorber properties is performed to optimize the performance. 3.1 Hand-Held Probe Design Figure 3.1 shows the hand-held device. Part labeled as (a), (c) and (d) are the moving head which vibrates tissue, the slide rods, and the CNC machined aluminum shell respectively. Parts (g) and (h) are the primary actuator's (LA-14-17-000 A, BEI Kimco, San Marcos, USA) coil and magnet. Appendix A gives details of the actuator selection process. This actuator vibrates the tissue, while the ultrasound probe (b) images the tissue. The absorber actuator is mounted in an opposite way: its coil (j) is fixed and the magnet (k) moves. Therefore the absorber mass is simply the magnet of the actuator. This choice obviates the need for an additional absorber mass, making the hand-held device lighter. Part (1) is one of the two absorber springs. There is no damper in the absorber; the inevitable Coulomb friction between the magnet and the housing and also the rods and the bushings play the role of a damper. Parts (f) and (i) are a linear potentiometer (LCP8S-10, ETI Systems, Carlsbad, USA) and a piezoelectric accelerometer (4508-002, Bruel and Kjaer, Naerum, 34 3.1 Hand-Held Probe Design 35 Denmark) respectively. The accelerometer is mounted on the absorber mass (the magnet) and provides the only feedback for the vibration absorption algorithm. The potentiometer signal is only used to operate the first actuator under closed loop control for displacement control. Here a conventional P ID controller with the coefficients Kp = 1, Ki = 0.2 and Kd = 0.002 is used. Since the mathematical model of the vibrator is unknown, the coefficients of the P I D controller are obtained by the Ziegler-Nicholas method [47]. To prevent the magnetic fields of the two linear actuators interfering, especially since the absorber mass (k) is the magnet, a hollow spacer (e) is placed between them. Figure 3.1: Exploded view of the hand-held probe: (a) moving head; (b) ultrasound probe; (c) slide rod; (d) machined shell; (e) spacer; (f) linear potentiometer; (g) and (h) primary actuator's coil and magnet; (i) accelerometer; (j) and (k) absorber actuator's coil and magnet; (1) absorber spring. 3.2 Design of the Vibration Absorber Parameters 36 3.2 Design of the Vibration Absorber Parameters The required operational frequency range, as mentioned in the device requirements in the first chapter, of the hand-held device is 5-25 Hz. The design process involves selecting values for ma, ka and ba such that the probe is as light and small as possible (i.e. small absorber mass (ma) and small absorber mass displacement). The design procedure of a passive absorber will be briefly described next since it can be used to obtain an estimate of the active absorber parameters. These values are different from the optimal values obtained for an active DVA, but they can be used as an initial guess for the design of an active DVA. To design a passive DVA, equation (2.1) and equation (2.2) should be used. Mass and spring elements determine the natural frequency of the absorber. This frequency usually should be selected close to the frequency of excitation force. Damping is added to a passive DVA to prevent resonance or to improve the effective bandwidth of operation of a vibration absorber. The mass, stiffness and viscosity values can be determined by minimizing the maximum of the combined system's transfer function over the operation frequency bandwidth (a possible regularization). This optimization is performed over the physical bounds of the absorber mass, absorber damping and absorber stiffness. Optimized values of the stiffness and damping usually are between their physical bounds, whereas the optimum value of the absorber mass is always found at the upper bound ( [31], [52] and [45]). In active DVA, the design procedure starts with finding how the required control force, Fc, and the absorber mass displacement, Xa, change with the selection of ma, ka and ba over the 5-25 Hz frequency interval. Fc and Xa are important factors in actuator selection. The method that we introduce here can be used for all active DVA methods that make the absorber mass resonate, including the DR-DVA and PI-DVA methods. The selected values for ma, ka and ba should ensure the stability of the combined system throughout the operation frequency range. Having the final choices for ma, ka, ba for the control methods, the choice of DR-DVA or PI-DVA can be made by comparing the absorption level achieved, and the tuning speed. 3.2 Design of the Vibration Absorber Parameters 37 3.2.1 Effect of Absorber Parameters on Control Force and Absorber Mass Dis-placement The required control force and the absorber mass displacement are determined by the absorber system parameters ma, ka and ba. As mentioned before, the best results always happen at the upper physical bound of the absorber mass. Therefore, we select a maximum possible mass, ma = 0.2 kg 1 , and keep it constant. The ka and ba values should be designed according to the operational frequency bandwidth requirements. The objective is, of course, to minimize the control force, Fc (Figure 3.2 (a)) and the absorber mass displacement, Xa- Fc varies with time but, for actuator selection, the steady state value is most important. Similarly, the steady state value of Xa is desired since a stable system soon reaches steady state. To find Fc and Xa, DR-DVA and PI-DVA should be treated separately. For DR-DVA, starting from the gain value (equation (2.14)), Fc and Xa can be found using equation (2.24), equation (2.25) and equation (2.26). For PI-DVA, starting with the proportional and integrator coefficients (equa-tion (2.22) and equation (2.23)), Fc can be found using equation (2.24), equation (2.25) and an equation similar to equation (2.26) with the right hand side replaced with the PI-DVA control force (equation (2.17)). A computer simulation for solving the dynamic equations of the system (equation (2.24), equation (2.25) and equation (2.26)) can be used for this purpose. We use a much easier way that actually gives both Fc and Xa for both DR-DVA and PI-DVA. Figure 3.2 (a) shows an active DVA attached to a 1 DOF primary system. Using either DR or PI control to drive the actuator, the system will be equivalent to Figure 3.2 (b) at steady state. In this figure, the spring ka, the damper ba and the actuator that connect the absorber mass to the primary mass are replaced by a single spring with the value of k = macj2 where CJ is the excitation frequency and the frequency that the active DVA is tuned to. This is true because both DR-DVA and PI-DVA methods make the absorber part a resonator at the frequency of excitation2 and therefore the absorber system will be identical to the absorber system depicted in Figure 2.3. 1lt will be shown that this value is the mass of the magnet of the electromagnetic actuator, which is used as the absorber mass of the active DVA. 2Although DR-DVA has infinitely many roots, this statement is still true. That is because all of the roots in the stable DR-DVA, of course except the two on the imaginary axis, have negative real part. Therefore they can all be neglected in steady state. The same argument applies to the extra roots that are introduced to the DVA system by the dynamics of the electromagnetic actuator and required filters that are used to modify feedback signals. They all have negative real values provided that the system is designed to be stable. 3.2 Design of the Vibration Absorber Parameters 38 In Figure 3.2 (b), the primary mass is fixed, which happens at the steady state condition provided the DVA is tuned to the excitation frequency. Writing Newton's law for the mass mp in Figure 3.2 (b), one obtains Fs = Fesin(wt) (3.1) where F3 represents the force of the imaginary spring k = maw2 and Fe is the amplitude of the excitation force. The displacement of the absorber mass Xa can be found as x « = T- = (3-2) k mau>z Coming back to Figure 3.2 (a), the force of the damper ba is ba(xa — xp) and the force of the spring ka is ka(xa — xp), which simplify to baxa and kaxa at the steady state. At the steady state, the summation of the forces of the spring ka, damper ba and the actuator in Figure 3.2 (a) should be equal to the force of the imaginary spring k = maui2 in Figure 3.2 (b): kaxa + baxa + fc = mau>2xa or equivalently fc = (maw2 — ka)xa — baxa. Therefore, the amplitude of fc (called Fc) can be calculated as Fc = V ( m a o ; 2 - ka)2 + b2L02 Xa. (3.3) Or, equivalently, knowing that xa(t) = Xasin(u)t) in the steady state condition and using equation (3.2) to replace Xa with Fe in equation (3.3) FK = i2 UI where una = \fka/ma and £a = ba/2y/kama are respectively the absorber system natural frequency and damping coefficient. In equation (3.4), u>na and Ca are the parameters of the absorber system and OJ is the excitation frequency. The question is now how to determine the value of Fe and whether it changes. Fe depends on many factors, the most important ones are: the required tissue vibration amplitude, the mechanical properties of tissue to be excited, and the frequency of excitation. Fe has both transient and steady state components. Similar to the arguments regarding Fc and Xa, our interest is only in the steady state part. In VE, a predetermined amount of strain is required (amplitude of excitation). Too little strain 3.2 Design of the Vibration Absorber Parameters 39 '7M m„ 1 m Fesin(wt) k = maw2 ]///// (a) (b) Figure 3.2: Simplified imaginary model of active DVA. (a) Active dynamic vibration absorber, (b) Equivalent model in steady state condition considering a control force Fc applied according to DR-DVA or PI-DVA control law. results in poor displacement estimation and too much strain leads to nonlinear tissue behaviour. Therefore for a specific tissue we can assume the vibration amplitude is fixed through the operational frequency range. Regarding the mechanical properties of the tissue to be excited, the values can change over several orders of magnitude depending on the subject and tissue type. Again to have an estimate of the excitation force as a function of excitation frequency, we consider tissue as a simple spring-damper model [23] (Figure 3.3). The differential equation of motion for such a system is btxt + xtkt = Fesin(u>t) where xt is the displacement of the tissue, and bt and kt are tissue equivalent damping and stiffness values. Solving the inhomogeneous differential equation, the steady state answer for the xt can be found as Xtsin(u)t) where Xt = , n 1 „F e (3.5) The elasticity modulus of muscle and liver tissue are reported to be approximately 1.5 kPa and 0.35 kPa respectively [9]. Unfortunately the value for tissue viscosity is difficult to measure and is rarely found in the literature. [37] reports the shear viscosity of a gelatin phantom (which has the elasticity value of 9 kPa) as 15.3 Pa.s. Although they are very rough estimates, they are sufficient for the design step of investigating excitation force over frequency. Plugging in these values results in btu> « kt for the 5-25 Hz excitation frequency. Nevertheless, 3.2 Design of the Vibration Absorber Parameters 40 Figure 3.3: Spring-damper model of the tissue. to show that ratio of btu) and kt does not have significant result on the answers, both extremes are shown in the plots of Fc vs. u> and Xa vs. ui. In other words we look at the cases where the stiffness force is much bigger than damping force (kt » btuj) and where the damping force is much bigger than the stiffness force (btu) » kt). Figure 3.4 (a) shows the required control force (equation (3.4)) for both extremes for two damping coefficients of the absorber system. As this figure indicates the results of the two extremes are not very different. In this figure, the natural frequency of the absorber is considered to be (5+25)/2=15 Hz (ujna = 15 Hz) as the first guess. As we can see, this guess results in bigger Fc at the beginning of interval (5Hz) than at its end (25 Hz). Minimizing the maximum of Fc results in ujna = 9Hz with the maximum Fc value of approximately 2.8 N (Figure 3.4 (b)). Accordingly, the absorber spring value can be determined: ka = mau>la = 640AT/m The spring is implemented as two Spaenaur B-646 springs in series each having a rated stiffness of 265 N/m, resulting in ka< nominai — 2x 265N/m = 530 N/m (It will be shown later that the actual stiffness value of the real spring is 637 N/m, slightly larger than the stated value). Figure 3.4 (c) shows the required stroke of the actuator (equation (3.2)) as a function of the excitation frequency. This value does not depend on ka or ba and is only a function of ma (equa-tion (3.2)). To plot these diagrams, the value of Fe in equation (3.2) and equation (3.4) at ui = 5 Hz is obtained by assuming a strain (e) of 0.01 for a E=10 kPa tissue with a 10 cm x 10 cm cross section area (A) and length I: EA Fe = ktAl = — x el = EAe = 104 x 0.12 x 0.01 = IN at UJ = 5Hz 3.2 Design of the Vibration Absorber Parameters 41 The viscosity of tissue is neglected in this equation and the values for E and A are overestimated since an approximate value of Fe suffices. It should be mentioned again that the Fc and Xa plots are independent of the control method (DR or PI). Also, the calculated value for Fe acts as a scale in Figure 3.4 and has no effect on determining the best u>na value. Figure 3.4 (a) indicates that a smaller control force corresponds to smaller valus of ba and, therefore, a small value for ba is desired (theoretically zero). In fact, a significant value of ba is inevitable because of Coulomb friction between the moving mass and the housing, as well as the rods and the bushings. Given these initial design choices, the next section discusses stability, both for the absorber system and the combined system. To obtain the stability results, the properties of the actuator are needed. A BEI LA-14-17-000A linear voice coil is selected as the absorber actuator. Appendix A provides a detailed guide on how the actuator was selected, and includes issues like peak force requirement, RMS force requirement and total stroke. The electrical parameters of the actuator are listed in Table 3.1. — k » o w C,- c w » k Cj, k » c w C^, c w » k w(Hz) v» (Hz) w(Hz) Figure 3.4: Control force and absorber mass displacement vs. excitation frequency. (a) Fc as a function of to. The solid and the dashed curves correspond to kt » hto and the dotted and dashed-dotted curves correspond to btto >> kt. Two different absorber system damping coefficients of Cal = 0.1 (the solid and the dotted curves) and Ca2 — 0.5 (the dashed and the dashed-dotted curves) are shown. In all curves tona = 15 Hz. (b) Minimizing the maximum of Fc by moving tona to 9 Hz. It is assumed that kt w btto and Ca = 0.5. (c) Xa as a function of LO. The solid lines and dotted lines represent kt » btuo and bTLO » kt respectively. Fo all graphs, Fe = 1 N at LO = 5 Hz. 3.2 Design of the Vibration Absorber Parameters 42 R L Kf Kb 6.5 n 1.07 mH 5.56 N/A 5.58 vs/m Table 3.1: Properties of the electromagnetic actuator stated by the manufacturer. R and L are the DC resistance and the inductance of the coil respectively, Kf is the force sensitivity and Kb is the back emf constant. 3.2.2 Effects of Absorber Parameters on the Combined System Stability The ma = 0.2 kg and ka = 637 N/m values designed in the previous chapter should provide absorption in the 5-25 Hz frequency interval. Therefore the stability of both the absorber system and combined system should be ensured in this interval. Also, the effect of the remaining absorber parameter, ba, on the stability of the combined system will be explored. 3.2.2.1 D R Stability In Section 2.3.1, it was shown how the stability of the absorber system dictates a minimum tuning frequency for the DR-DVA. Section 2.3.1.1 indicates there is a maximum tuning frequency for the DR absorber. In Section 2.4.1 the stability of the combined system was studied. This stability condition gives a new maximum and minimum tuning frequency for DR that, in general, is not the same as the maximum and minimum tuning frequencies obtained by studying only the absorber system. The minimum tuning frequency of the DR absorber system is 7.7 Hz with the gain (ga) set to 0.2134 vs2/m and the delay (ra) set to 0.0393 s and Ca = 0.16 (ba = 3.6 kg/s). This value changes a small amount to 7 Hz for Co = 0 and to 9.5 Hz for Ca = 0.3 (ba = 6.8 kg/s). These values are obtained from the crossing point of the first and second branches of equation (2.15) , similar to the procedure that was performed on equation (2.9) in Section 2.3.1. For the properties of the absorber system and the actuator given above, equation (2.16) in Section 2.3.1.1 gives the maximum tuning frequency 81 Hz with C = 0.16. This equation gives the minimum of the DR absorber maximum tuning frequency equal to 61 Hz, happening at Ca = 0 (equation (2.16)). The stability of the absorber system does not guarantee the combined system stability. The "D-subdivision" stability method is exploited here to analyze the combined system stability [20]. It 3.2 Design of the Vibration Absorber Parameters 43 originates from the continuity of the roots of quasi-polynomials with variable parameters g and r. That is, if there exists at least one path between the two sets of the roots of the quasi-polynomials p(s,gi,Ti) and p(s,52,^2) that does not cross the imaginary axis, then either both polynomials are stable or both are unstable. The stability is then defined in the boundary crossing theorem [7]: if the parameters of a stable (unstable) quasi polynomial with variable parameters change in an interval and its roots do not cross the imaginary axis, this polynomial is stable (unstable) within the parameter interval. In DR-DVA, it is guaranteed that the combined system is stable for g — 0 since it corresponds to zero control force (passive DVA). Therefore any point [g T] = [0 r°] in the g — T plane can be a starting point to look for stable regions. If a point [g* T*] can be connected to point [0 r°] without crossing the imaginary axis, then this point makes the polynomial stable. Before plotting the stability chart, one last practical modification to the DR formulation should be made. It was observed experimentally that the voltage applied to the absorber actuator has a high frequency noise component (several hundred Hz). Therefore the accelerometer data should be passed through a lowpass filter. An 80 Hz first-order lowpass Butterworth filter is chosen here. Accordingly, the equation (2.30) will be changed to D(s) = Kcs2(mps2 + bps + kp)F(s) (3.6) where F(s) = 1 + 3/( 1 2 7 r 8 0) represents the lowpass filter. Equations (2.14), (2.15), (2.31) and (2.32) are also modified similarly. Figure 3.5 shows the plot of gain and frequency versus delay for both the absorber system alone and the combined system. All plots include the dynamics of the actuator and the filter. The absorber and primary system properties are from Tables 3.2 and 3.3 and the parameters of the electromagnetic actuator are from Table 3.1. The properties of the primary system depend on how the probe is held and also on the operator. The mass of the whole hand-held device (excluding the absorber mass) is about 500 g and it is assumed the effective mass of the operator's hand is 500 g 3 . The 140 N/m value of kp (the stiffness of the hand) is obtained experimentally by seeing 3Since the hand is a distributed mass system and it vibrates around the arthroses at wrist, elbow and shoulder, the effective mass is only a percentage of the whole mass of the hand. As an example, for a rod with uniform mass distribution which rotates around a hinge at one of its ends, the effective mass is 1/3 of the total mass. 3.2 Design of the Vibration Absorber Parameters 44 that when the probe is held in normal condition, adding 1 N force (100 g) causes a displacement of approximately 1 cm. The bp = 2 kg/s value is a low value of damping (C,p = 0.08), since low damping values usually tend to make the system less stable4. The minimum operational frequency is 7.7 Hz which corresponds to g = 0.2216 vs2/m and T = 0.0388 s. The minimum operational frequency is a result of the second branch curve in the root-locus diagram of the combined system crossing the first branch curve in the root-locus diagram of the absorber system (Figure 3.5). This result is similar to the 7.7 Hz absorber's minimum tuning frequency; this is because of the specific values we chose for the primary system and is not true in general. The maximum tuning frequency is 24 Hz which corresponds to the zero delay of the absorber system (similar to the limit for equation (2.16) but now considering the lowpass filter). The corresponding gain and delay values for this operation point are g = 0.219 vs2/m and, of course, r = 0. K (u>a) ba (Ca) 0.2 kg 637 N/m (9 Hz) 3.6 kg/s (0.16) Table 3.2: Absorber system parameters. mp kp (wp) bp (CP) 1 kg 140 N/m (1.9 Hz) 2 kg/s (0.08) Table 3.3: Primary system parameters. Figure 3.6 shows other tuning frequency intervals that are found for different damping coeffi-cients of the absorber system. For each value of the ba in this figure, the stability chart (Figure 3.5) is replotted and the crossing points are found in order to find the minimum and maximum tuning frequency. Since the properties of the primary system can change, depending on how the device is held, a variations of kp in a range of 100-800 N/m and bp in a range of 2-20 kg/s (0.08 < Cp < 0.85) were investigated. Similar plots as Figure 3.5 were obtained by various kp and bp in these ranges. The results show that varying kp and bp has no significant effect on the tuning frequency interval (less than 0.5 Hz). 4These values of the properties of the primary system will be allowed to vary in some intervals since they vary for different persons holding the device. 3.2 Design of the Vibration Absorber Parameters 45 Branch 1, absorber Branch 2, absorber Branch 1, combined Branch 2, combined 0.08 0.08 Figure 3.5: DR stability chart. Gain and frequency versus delay. The plots show the 1st and 2nd branches of the absorber system and the combined system. The second branch plots of the absorber system and the combined system are very similar (the dashed curve and the dashed-dotted curve). - minimum frequency I • maximum frequency | Figure 3.6: Variation of the maximum and minimum frequency limit of the DR-DVA with absorber system damping. As a summary for the results obtained in this section, the range of tuning frequency of the absorber system of the DR-DVA with the properties from Table 3.2 is 7.7 Hz and 81 Hz respectively. The range of tuning frequency of the combined system with the DR-DVA, with the properties from 3.2 Design of the Vibration Absorber Parameters 46 Table 3.2 and Table 3.3 is 7.7 Hz to 24 Hz. These values change with the absorber system damping, as depicted in Figure 3.6. However, the variations of the properties of the primary system in the 2-20 kg/s for bp and 100-800 N/m for the kp do not have a significant effect on the minimum and maximum tuning frequency. 3.2.2.2 PI Stability Unlike DR-DVA, there is no lower limit for the PI-DVA tuning frequency in terms of stability. Of course, there are other limitations, like the amplitude of the absorber mass, that limit the operational range of any active DVA including the PI-DVA method (Figure 3.4 (c)). Theoretically, no upper frequency limit exists for PI-DVA either. Here, the frequency tuning speed becomes the key factor at high frequencies, as described in the Section 3.2.3. Therefore, the tuning frequency range is theoretically from 0 Hz to oo for the PI-DVA absorber system. However, the combined system may not be stable in this range, and like DR-DVA, the stability of the combined system must be ensured. Before plotting stability charts for the combined system, one last practical modification will be made to the PI-DVA formulation. It was experimentally observed that when the PI-DVA was tuned to absorb low frequency excitations, say 5 Hz, the first term of the control force (Kpxa term) has very high frequency components (around 900 Hz). Although the sound of the voice coil vibration due to these very high frequency voltages could be heard, these high frequency components did not deteriorate the vibration absorption performance. Nevertheless, to avoid damage to the voice coil and electronics, the accelerometer data is lowpass filtered before feeding it to the proportional gain term (Figure 3.7). The PI Controller P H Highpass accelerometer data -<3 Figure 3.7: Block diagram representation of the PI-DVA with the lowpass filter added. Duplicated from Figure 2.10 with the addition of the lowpass filter. 3.2 Design of the Vibration Absorber Parameters 47 With a first-order lowpass Butterworth filter, the characteristic equation of the combined system becomes N(s)+D(s) [Ulp + S e S+UJhp = 0 (3.7) where N(s) and D(s) are defined in equation (2.29) and equation (2.34) and Jf'+S and are the lowpass and highpass filters respectively. To solve this equation for Kp and Ki, we make these approximations: - ^ _ « ( 1 - S / W l p ) (3.8) —±—^(1-CJ^/S) (3.9) o I ^hp assuming that s/u>ip and uJhp/s « 1. With these approximations, equation (3.7) can be solved using both of its real and imaginary parts: i „ ) c > -t- u i „ i m n n < -=4-4 ^ = (1 _ ^ ) W 3 ( 3 - 1 0 ) u)comreal { ^ 14 } + wipimag { j ^- _ I. ' J s=ju>com l \ > J s=juip ^ ^ (u>ip - Uhp)uCom where ju)com is the crossing of the characteristic equation of the combined system with the imaginary axis. Also, the Kp and Ki values that make the absorber system resonate at uia (equation (2.22) and equation (2.23)) should be changed if the lowpass filter is used in the system: K _ "«real {A(s)}s=ju)a +UhPimag{A(s)}s=jula Kfd-^Pl (3"12) R = "» r e a l {^ (g)}s=jo,a + "Ipirnag {A(s)}s=jula Kf(oJip - UhP)ua where juja is the crossing of the characteristic equation of the absorber system with the imaginary axis and A(s) is defined as A(s) = (mas2 + bas + ka)(Ls + R) + KfKbs (3.14) 3.2 Design of the Vibration Absorber Parameters 48 There are many possible ways of determining the stability of the combined system for a PI controller. Here, the D-subdivision method is employed again. Figure 3.8 shows the plot of Kv versus Ki for both the absorber system and the combined system. In this figure, the properties of the absorber system and the primary system axe from Tables 3.2 and 3.3. The properties of the actuator are fom Table 3.1 and the highpass and lowpass filters axe as Table 3.4. The point [Kp Kj]=[0 0] (indicated with an asterisk) can be the starting point for the stability analysis: it represents zero control force or, equivalently, a passive absorber whose stability (both the absorber system and the combined system) is guaranteed. Moving from the point [0 0] in the Ki-Kp plane, the stability of the combined system is guaranteed as long as the imaginary axis is not crossed. This crossing happens at any point on the dotted line (the dotted line is obtained by plotting Kp in equation (3.10) versus Ki in equation (3.10) by varying u>Com- These equations give the Kp and Ki values that make the combined system marginally stable.). Therefore, the section of the absorber marginal stability (the solid line) that can be connected to the point [0 0] without crossing the dotted line corresponds to the stable region. This is anywhere on the solid curve befoxe its intexsection with the dotted curve, which happens at a frequency of 24.6 Hz frequency, with Ki = 0.58 vs/m and Kp = —0.2023 vs2/m corresponding values. No other crossing exists between the two curves, meaning that the combined system is stable in the frequency range of 0-24.6 Hz. DR-DVA lowpass (uip) PI-DVA highpass (uhP) PI-DVA lowpass (wip) 80 Hz first-order Butterworth 0.5 Hz first-order Butterworth 80 Hz first-order Butterworth Table 3.4: Filters used. Since the absorber system parameters (especially ba) and the primaxy system paxametexs (es-pecially c and k) can change ovex time and from experiment to experiment, the stability chart is also replotted, allowing for these vaxiations (the intexvals fox ba, bp and k are investigated for PI-DVA, are the same as those investigated for DR-DVA). The intersection point determines the high frequency limit of the PI-DVA, upimax- Figure 3.9 shows how wp/ m a i vaxies with ba and bp. The tJpimax changes by only 1 Hz when the kp vaxies fxom 100 N/m to 800 N/m, so is not plotted. Figuxe 3.9 (a) indeed reveals the upper frequency limit of the PI-DVA. Provided perfect lubrica-tion on the absorber part (ba = 0, the woxst case fox the combined system stability), the maximum fxequency of opexation is 19.3 Hz. 3.2 Design of the Vibration Absorber Parameters 49 0.5 absorber combined 120 100 120 Figure 3.8: PI stability chart. Kp and frequency vs. Kp. bp(ks(,) Figure 3.9: Maximum frequency limit variation with system properties, (a) wp/ m o l as a function of absorber damping ba. (b)wprmai as a function of primary system damping bp. The operational frequency of both DR-DVA and PI-DVA is highly influenced by the cutoff frequency of the lowpass filter, uip. To obtain the range of operational frequency for DR-DVA and PI-DVA versus wjp, stability charts (Figure 3.5 and Figure 3.8) are replotted for each value of the u>ip and the crossing points should be determined. Figure 3.10 shows how the maximum operational frequency increases by increasing the cutoff frequency. In this figure, toip is varied over 3.2 Design of the Vibration Absorber Parameters 50 a large frequency interval in order to see the asymptotical behaviour of the DR-DVA and PI-DVA; in practice, a cutoff frequency of more than 150 Hz is not appropriate for the hand-held device. The minimum operational frequency of the DR-DVA changes slightly (about 0.5 Hz) by changing the wip over the range of 80 Hz to 150 Hz. 11 Or Figure 3.10: Maximum frequency limit variation with lowpass filter cutoff frequency. As a summary for this section, the range of tuning frequency of the absorber system of the PI-DVA with the properties from Table 3.2 is 0 to oo (as opposed to 7.7 Hz to 81 Hz in DR-DVA). The range of tuning frequency of the combined system with the PI-DVA, with the properties from Table 3.2 and Table 3.3 is 0 to 24.7 Hz (as opposed to 7.7 Hz to 24 Hz in DR-DVA). These values change with the damping of the absorber system and the damping of the primary system, as depicted in Figure 3.9. However, the variations of the properties of the primary system in the 100-800 N/m for the kp does not have significant effect on the minimum and maximum tuning frequency, similar to the DR-DVA. 3.2.3 Vibration Absorber Tuning Speed So far, only the steady state performance of the active DVA has been considered, but the tran-sient response is also important when dealing with excitations with variable frequency or variable amplitude. The transient response is also important when dealing with excitations with constant frequency and constant amplitude excitations, as described in Section 2.5. Transient responses are usually very difficult to obtain from the equations of dynamics. The transient response of an active DVA is studied by the tuning speed of the active DVA here. 3.2 Design of the Vibration Absorber Parameters 51 The PI-DVA tuning speed is studied first, because it is easier to show that, with a PI controller, the chosen electromagnetic actuator is sufficiently fast for the frequency range of the device. Then, the tuning speed of the DR-DVA is investigated. 3.2.3.1 PI Tuning Speed Neglecting the dynamics of the electromagnetic actuator, the characteristic equation of the PI-DVA (equation (2.19)) has no roots other than the two roots on the imaginary axis (the two roots at zLjuJa that make the absorber mass resonate). Therefore, the tuning speed is theoretically infinity in this case. The dynamics of the actuator are more important when the absorber is tuned to resonate at high frequencies. To study the effect of the electromagnetic actuator dynamics on the frequency tuning speed, we study the characteristic equation of the absorber system. It is easy to show that the characteristic equation of the absorber system is: CEpi{s) = maLss + (Lba + Rma - KfKp)s2 + (LKa + Rba + KfKb - KfKJs + Rka = 0 (3.15) where CEpi(s) is the characteristic equation of the PI-DVA, neglecting the 0.5 Hz highpass filter5 (Figure 3.11 (a)). To tune the absorber to frequency ua, Ki and Kp should be calculated according to equation (2.22) and equation (2.23), again neglecting the highpass filter (term ^%-Ki). Inserting the obtained Ki and Kp into equation (3.15), the CEpi(s) will have three roots: S i ^ = ±jw a a n d S3 = —a. The third root determines the tuning speed of the PI-DVA at the frequency wa. After inserting the numerical values of the system parameters from Tables 3.2 and 3.1 we find where toa is in Hz. With this first order pole, the absorber has a time constant of 1/a and a settling time of 4/a. The decrease in the absorber tuning speed by increasing the tuning frequency can be best observed in the simulation of the PI-DVA with the electromagnetic actuator (Figure 3.11). In these simulations, the absorber mass has an initial displacement (x) of 0.1 m and no 5The highpass filter can be neglected since we only want to investigate the effect of the dynamics of the electro-magnetic actuator a = 4.9 x 105 x w, a - 2 (3.16) 3.2 Design of the Vibration Absorber Parameters 52 external excitation force is applied. The absorber parameter values are as of equation (3.16) . In the parts (b) and (c), Ki and Kp are set to 250 Hz and 500 Hz respectively, using equation (2.22) and equation (2.23). The settling time increases from 0.49 s to 1.99 s when the tuning frequency is increased from 250 Hz to 500 Hz, as equation (3.16) indicates. Equation 3.16 indicates that the third root is far enough from the imaginary axis (in the other words the actuator is fast enough) in the frequency range that the absorber operates. For example, at the frequency of 25 Hz, the settling time would be A/a — 4/784 ss 0.005 s which is negligible for a excitation with a 2 Hz/s sweep rate. 3.2.3.2 DR Tuning Speed For PI-DVA, neglecting the dynamics of the actuator (L = 0 in equation (3.15)) makes CEpi a 2nd polynomial and, therefore, the tuning speed is oo. This is not the case however in DR-DVA: even neglecting the dynamics of the actuator, the system has infinitely many roots which slow the absorber tuning speed. Equation 2.7 gives the characteristic equation of the DR absorber system neglecting the actuator's dynamics (repeated here for convenience) rnas2 + bas + ka + gs2e~TS = 0. (3.17) Among all of these roots, the two roots with the smallest real part (closest to the imaginary axis) are the dominant poles, except, of course, the two on the imaginary axis. In order to analyze the transient response of the absorber, these dominant poles should be found. There are two reasons that make finding these roots difficult with existing numerical methods: first, the roots are complex roots of a transcendental equation. Second, there are infinitely many roots and it's not easy to make the numerical method converge to the specific root that we are interested in. By looking at the DR root-locus plot, knowing that the dominant poles are on the second branch (Figure 2.3.1), we select the starting point as [<r° w°] = [a* 27r/ra] where [a w\ is the root of the equation (3.17) and ra is the delay value that makes the absorber resonate at the frequency of ua. The real part of the root, a, determines the tuning speed of the DR-DVA. Figure 3.12 shows how the real part of the dominant poles changes with the tuning frequency. The absorber parameters here are the same as Table 3.2. The minimum tuning frequency to keep 3.2 Design of the Vibration Absorber Parameters 53 - 2 1 ' ' ' ' 1 O1 1 ' ' ' 1 3 3.02 3.04 3.06 3 08 3.1 3 3.02 3.04 3 06 3.08 3.1 time (s) timo (s) Figure 3.11: PI tuning speed, (a) Active DVA. (b) Simulation of PI-DVA tuned to resonate at 250 Hz. (c) Simulation of PI-DVA tuned to resonate at 500 Hz. (d) A segment of part (b) from 3 s to 3.1 s. (e) A segment of part (c) from 3 s to 3.1 s. the absorber system with these parameters stable is 7 Hz, as mentioned in Section 2.3.1. In Figure 3.12 the point 7 Hz is associated with a = 0, which is in agreement with the results of Section 2.3.1. Figure 3.12 suggests that it is better to avoid the minimum operational frequency range of the DR-DVA since the dominant poles are very close to the imaginary axis. For example, at the frequency of 8 Hz, the dominant poles have the real value of cr = 3.46 l/s which translates into a time constant of 4/<7 = 4/3.46 = 1.16 s. This settling time is not only important in dealing with excitations with variable frequency, but also for the excitations with variable amplitude (which theoretically means that the excitation is not a single frequency excitation). Also, this large settling time is important even for excitations with constant frequency and constant amplitude, as was explained in the Section 2.5. We therefore see that there is a tradeoff between the settling time and the range of operation: the smaller the desired settling time, the farther the operation frequency is from the minimum operational frequency. 3.2 Design of the Vibration Absorber Parameters 54 Figure 3.12 is obtained when the first branch of the root locus plot (Figure 2.3.1) is used as the loci that carries the resonator poles (±jwa). A similar analysis show that the settling time increases dramatically when other branches are used to carry the resonator poles. Therefore only the first branch crossing with the imaginary axis (I = 0 in equation (2.15)) is shown here. 3.2.4 Designed Values for the ma, ka and ba The following is a summary of the previous sections. Section 3.2.1: The absorber mass, ma, was assumed to be 0.2 kg, a maximum tolerable value. The absorber spring, ka, was selected to be approximately 640 N/m. This value was obtained by minimizing the maximum control force, Fc. There is no control over the absorber mass amplitude of vibration: it is determined by the absorber mass and is independent from absorber spring and damper. Section 3.2.2: The stability charts revealed the relationship between the absorber system damp-ing coefficient, £a, and the stable frequency range. For DR-DVA, the stable frequency range is 8-25 Hz for Co = 0.2 or equivalently ba = 4.5 kg/s. For PI-DVA, the stable frequency range is 0-25 Hz for the same damping value (using an 80 Hz lowpass filter). Since no damper is incorporated in the device and the uncontrollable friction and Eddy currents play the role of damping, there is no direct control of the £a value. By increasing the lowpass filter cutoff frequency the maximum operable frequency can be increased for both DR-DVA and PI-DVA methods (Figure 3.10). Fortunately, it was observed from the experiments that the lowpass filter is mainly required for low frequency excitations. Therefore if the damping of the device is insufficient to give the maximum 25 Hz tuning 0 20 40 GO 80 Wa (Hz) Figure 3.12: Dominant poles real part vs. tuning frequency. The dynamics of the actuator is not considered in this plot in order to simply show the effect of DR infinite poles and not the actuator pole. 3.3 Identification of the Absorber Parameters 55 frequency and a 25 Hz tuning frequency is desired, the lowpass filter cutoff frequency (80 Hz in Figure 3.8) can be increased such that 25 Hz falls into the system's stable frequency interval. Section 3.2.3: The performance of an active DVA is determined by its tuning speed (for single frequency excitation as well as sweeps). The DR has a large settling time close to its lower frequency bound, that makes it a slow absorber at these frequencies. That is because of the fact that the characteristic equation of the DR-DVA has infinitely many roots , and, close to the lower stability bound, two of these roots are very close to the imaginary axis. On the other hand, the PI-DVA absorber is much faster because the only extra roots are from the electromagnetic actuator dynamics (which creates the same extra root in DR-DVA) and the highpass and lowpass filters used (which can be shown that create poles close to the poles to the lowpass filter used in DR-DVA creates). 3 .3 I d e n t i f i c a t i o n o f t h e A b s o r b e r P a r a m e t e r s After manufacturing the device, the exact values of raa, ka and ba need to be measured. The value of ma was measured by a scale and is 0.199 kg. To measure ka and ba, a free vibration test on the absorber mass was done and the acceleration was measured by the accelerometer (Fig-ure 3.13). Since the accelerometer data is noisy, its noise should be removed before parameter identification. The data is integrated twice, and the drifts are subtracted in each integration, re-sulting in absorber mass displacement. This method is preferred over filtering since it does not distort the phase or change the amplitudes. Assume that the free vibration is of the standard form x(t) = Ae~^aUnatsin(<jJdat-r4>) where ujna = yj(ka/ma) is the absorber system natural frequency and Ca is the absorber system damping coefficient and is the frequency of absorber system damped oscillation. The two unknowns ujna and Ca can be found by matching the measured motion to the free vibration equation. This gives una = 8.9 Hz and Ca = 0.16. These translate into ka = 637 N/m, bali = 3.6 kg/s, and ui^a = 8.9 Hz. ka can also be found using the equation k = 8 ^ ja ( [27]) where n e, G, d and D are the number of effective coils, shear modulus, wire diameter and coil diameter of the spring respectively. For the Spaenaur B-646 spring n = 8, d = 0.018 in, D = 1/4 — 0.018 in and G = 11.6 x 106lb/in2. Using two springs in series and using ne = 8 — 1.5 = 6.5 because of the two clamps at the two ends of the springs, gives ka = 653 N/m which agrees closely with the ka = 637 N/m value obtained from the free vibration test. 3.3 Identification of the Absorber Parameters 56 10 5 0 t ~ 5 -10 -15 -20 1.8 1.9 2 2.1 2.2 1(s) Figure 3.13: Free vibration experiment results. The Coulomb friction and the Eddy currents induced by the moving magnet in the housing are the sources of damping, and the given value for ba is the result of linearization. Coulomb friction force induces several undesirable nonlinear effects to the system, such as nonlinear instability phenomena [26] and chaos for small excitation frequency and amplitude force vibration [11]. To investigate the contribution of Coulomb friction and Eddy currents to the damping observed in the absorber system, a forced vibration test at 20 Hz is performed on the absorber mass. Figure 3.14 shows the accelerometer data obtained from the test. Strange behaviors can be found look-ing carefully at the maximum (or minimum) points of the acceleration graph. Since these points correspond to the ends of the stroke, they can be associated with actuator saturation or other nonlinearity sources in the system like Eddy currents or Coulomb friction. To eliminate the possi-bility of actuator saturation (at its stroke limits), a similar test with smaller excitation amplitude is carried out (Figure 3.14 (b)). Since the nonlinear behavior at the stroke limits is still visible, the possibility of actuator saturation is eliminated. Also, since the speed is zero at the stroke limits, Eddy forces are zero, so we conclude the nonlinearities axe caused mainly by Coulomb friction. An approximately 0.85 m/s2 discontinuity can be seen in both graphs at the stroke limits. With an absorber mass of 0.2 kg, Coulomb friction force is approximately: 2 x F/ ~ 0.85 x 0.200 = 0.17 N or Ff w 0.09 N. Ff is multiplied by 2 because the discontinuity is the result of a direction change of the friction force. To find the equivalent viscous damping, two methods exist. The first method equates the energy dissipated by the nonlinear and equivalent linear damping in one cycle [62], and the second method is the describing function linearization [54]. Both methods lead to a ba = 0.4 kg/s equivalent viscous damping (at Xa = 5 mm and toa = 8 Hz which are the approximate amplitude 3.4 Summary 57 <•) (b) Figure 3.14: 20 Hz forced vibration experiment results, (a) Absorber mass acceleration vs. time.(b) Absorber mass acceleration vs. time, smaller excitation force. and frequency of the free vibration experiment). This shows that about 10% of the damping emanates from Coulomb friction and therefore Eddy currents is responsible for the reminder. 3 .4 S u m m a r y The components of the designed hand-held device were described in this chapter. A conventional PID controller was used to provide smooth tissue vibration. The design of the DVA mass, spring and damper parameters was presented in detail. It was shown that the control force and the absorber mass displacement are the same in DR-DVA and PI-DVA. The stability of the combined system with both DR-DVA and PI-DVA (considering the filters used in each method and the dynamics of the electromagnetic actuator) was studied next. DR-DVA and PI-DVA have similar upper bounds on the operational frequency, but different lower bounds. The tuning speed of the DR-DVA and PI-DVA was then examined. It was shown that the electromagnetic actuator selected in this work is sufficiently fast for the 5-25 Hz frequency range, and that its dynamics do not decrease the tuning speed of the DR-DVA or PI-DVA. However, it was shown that the extra roots of the DR-DVA characteristic equation deteriorate the vibration absorption performance significantly in some frequencies. Parameter identification of the parameters of the absorber system was performed last to optimize the performance of the vibration absorber. The easy derivation of control force and absorber mass displacement using the imaginary equiv-alent model of Figure 3.2 (b) is novel in this work. A lowpass filter was added to the PI-DVA in 3.4 Summary 58 order to make it a suitable control approach for noisy feedback. The comparison of the tuning speed between PI-DVA and DR-DVA is another novel aspect. The tuning speed of PI-DVA (Figure 3.8 (b) and (c)) and the relation of DR-DVA dominant poles with respect to tuning frequency (Section 3.2.3.2) are other novelties of this work. C h a p t e r 4 Simulation and Experimental Results The performance of the DR-DVA and PI-DVA systems are investigated in this chapter by simulation of both the absorber system and the combined system with different kinds of excitation forces. The simulation environment is implemented by Simulink and MATLAB. The regions of stability obtained in the previous chapters will also be explored by simulation. The claims on performance and the obtained poles of the characteristic equations of the DR-DVA and PI-DVA will be verified here by two different approaches: first by looking at the resonance behavior of DR-DVA and PI-DVA, and second by investigating their vibration absorption perfor-mance. Experiments with the overall system are also performed on different types of gelatin phantoms with different mechanical properties mimicking different tissue types. 4.1 Simulation Results In this section, the primary system is assumed to be a 1 DOF mass spring damper system and the vibration absorber is a mass-spring-damper trio. The absorber system and primary system properties for all simulations are from Tables 3.2 and 3.3 unless otherwise specified. 59 4.1 Simulation Results 60 4.1.1 Passive D V A Figure 4.1 shows the simulation results of the probe motion and absorber mass motion assuming a passive DVA is used to cancel probe vibration at three excitation frequencies of 8.8 Hz, 7 Hz and 11 Hz. As expected, because of the damping ba, even at the absorber system's frequency of resonance (uja « 8.8 Hz) the vibration suppression is not perfect and residual vibration can be observed in the primary mass (Figure 4.1 (a)). Changing the frequency of excitation from the absorber resonance frequency decreases the level of vibration suppression (the ratio of the amplitude of the primary mass vibration with and without vibration cancellation). The corresponding vibration absorption is 32% at coe = 8.77 (the maximum absorption value, obtained at the resonance frequency of the absorber system), 26% at uje = 7 Hz and -14% at cue = 11 Hz. Note that at the excitation frequency of 11 Hz the passive DVA works as a vibration amplifier. First, the resonance behaviour of the DR-DVA is studied. The vibration absorption performance of the PI-DVA is studied next. Figure 4.1: Passive vibration absorption simulation results. Absorber mass displacement (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies: (a) 8.77 Hz, (b) 7 Hz, and (c) 11 Hz. The amplitude of excitation force is 1 N for all three simulations. 4.1 Simulation Results 61 4.1.2 D R - D V A 4.1.2.1 Absorber System The absorber mass-spring-damper trio is considered first. Figure 4.2 (a) shows the free vibration of the trio with no control force, which is clearly dissipative and stable. To make the trio marginally stable with two roots crossing the imaginary axis at ±ju>a = ±j27rl0, DR control parameters ga and ra are calculated according to equation (2.14) and equation (2.15) and the DR control force is applied to the absorber mass, resulting in a 10 Hz resonator (Figure 4.2 (b)). From the DR-DVA root locus plot it is clear that increasing the gain moves the poles on the imaginary axis to the RHP and makes the system unstable. Figure 4.2 (c) shows the instability of the system by increasing g from the ga value by 5%. Also, there is the minimum frequency of 7.7 Hz for the DR-DVA operational range as mentioned before. By tuning the DR-DVA below this level, the roots on the second branch move to the RHP of the complex plane making the absorber system unstable. Figure 4.2 (d) shows the simulation results when the DR-DVA is tuned to resonate at 6.5 Hz. The roots on the second branch which have a imaginary part of approximately ±j27rl5 have moved to the RHP of the complex plane. 4.1.2.2 Combined System A DR-DVA can completely cancel the primary system motion provided it is properly tuned to the excitation frequency. Figure 4.3 (a) and (b) show the simulation results of 8.77 Hz and 15 Hz excitations. In both simulations complete vibration suppression is achieved. The combined system stability must be ensured for a practical system. Figure 4.3 (c) shows the system becomes unstable when the DR-DVA is tuned to the frequency of 6.5 Hz, outside the 7.7-24 Hz range of stability. Theoretically 100% vibration absorption is achieved through the operational range of frequency at the steady state. 4.1.3 PI-DVA Again, the resonance behaviour of the PI controlled mass-spring-damper is studied first, followed by the vibration absorption performance of the PI-DVA. 4.1 Simulation Results 62 A A A 0.1 0.2 0,3 0.4 0.5 t(S) 0.2 0.3 0.4 0.5 U S ) Figure 4.2: Absorber mass-spring-damper trio free vibration simulation, (a) No control force, (b) DR-DVA tuned to 10 Hz by selecting ra = 0.0157 s and ga = 0.069 kg. (c) Increasing ga to ga = 0.069 x 1.05 kg. (d) Tuning DR-DVA to 6.5 Hz. | Xa | | xp | HI T. (c) Us) Us) Figure 4.3: Combined system simulation results with the DR-DVA. Absorber mass displace-ment (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies (the DR-DVA tuned to the excitation frequency): (a) 8.77 Hz, (b) 11 Hz, and (c) 6.5 Hz. The amplitude of excitation force is 1 N for all three simulations. 4.1 Simulation Results 63 4.1.3.1 Absorber System The performance of the DR-DVA and PI-DVA are very similar since they both place the roots of the absorber system characteristic equation on the imaginary axis. There are two major differences: unlike DR-DVA, there is no lower limit for PI-DVA operational range of frequency. Figure 4.4 (a) shows the mass-spring trio tuned to resonate at 6.5 Hz, which was impossible to do with the DR control law (Figure 4.2 (d)). It is possible to tune the absorber system in the PI-DVA to even lower frequencies. Figure 4.4 (b) shows PI controlled trio tuned to vibrate at the frequency of 2 Hz. The second difference arises from the fact that DR-DVA has infinitely many roots and the two dominant roots among these roots slow the tuning speed. Figure 4.5 shows the difference between the DR and PI controlled mass-spring-damper trio when their tuning frequency is suddenly changed by changing their control parameters (g and r for DR and Kp and Ki for PI). The tuning frequency is changed from 15 Hz to 7.5 Hz by changing the control parameters (the delay and gain in DR, and the proportional and integrator coefficients in PI) and the response of the two systems is observed. The DR-DVA becomes tuned to 7.5 Hz from 15 Hz after a settling time, a transition that requires much less time for the PI-DVA (too small to be observable in the graph). The very small settling time required for the PI-DVA is in fact associated with the actuator dynamics and niters used, unlike DR-DVA that has an intrinsic settling time. (t>> 1 .5 1 1.5 2 t ( S ) Figure 4.4: Simulation of the PI controlled mass-spring-damper trio vibration, (a) Tuned to 6.5 Hz. (b) Tuned to 2 Hz. E " 1 o a X-0,2 -0,4 -0.6 -0.8 -1 —. 0.2 E E o x -0.2 -0.4 -0.6 -0.8 t(S> 4.1 Simulation Results 64 t(s) l(s) Figure 4.5: PI and DR controlled mass-spring-damper trio vibration simulation, (a) PI. (b) DR. Both plots involve changing the resonance frequency from 15 Hz to 7.5 Hz at t=2 s. 4.1.3.2 Combined System Like DR-DVA, 100% vibration absorption of the primary system can be achieved by PI-DVA in the steady state. There are two differences in their vibration absorption performance however. The first difference is that, as mentioned before in the stability analysis, there is no lower operational limit for PI. Figure 4.6 shows the simulation results of the combined system with PI-DVA. The primary mass is subjected to a 8.77 Hz, 11 Hz and 6.5 Hz excitation force in parts (a), (b) and (c) respectively, and its vibration is perfectly cancelled with the vibration absorber. Figure 4.7 shows the DR and PI control force required to suppress the primary system vibration Figure 4.6: Combined system simulation results with the PI-DVA. Absorber mass displace-ment (xa, dotted line) and primary mass displacement (xp, solid line) at three different excitation frequencies (the PI-DVA tuned to the excitation frequency): (a) 8.77 Hz, (b) 11 Hz, and (c) 6.5 Hz. The amplitude of excitation force is 1 N for all three simulations. 4.1 Simulation Results 65 when it is excited with a 8.77 Hz excitation force. As anticipated before (Section 3.2.1), the steady state control force is the same for DR-DVA and PI-DVA. The only difference is in the transient control force which is because of extra roots in the DR-DVA. The same behavior is observed in the absorber mass displacement (xa) in the DR-DVA and PI-DVA (Figure 4.3 and Figure 4.6, a different transient but same steady state response) again as anticipated in Section 3.2.1. 4.1.4 The Effect of Tuning Speed on the Vibration Absorption The second difference between PI-DVA and DR-DVA is the difference between their tuning speed. We showed in Figure 4.5 that PI-DVA can be tuned to different frequencies almost instantaneously, while DR-DVA requires some time for tuning. This difference comes from the fact that DR-DVA has infinitely extra roots (other than the two roots that make the absorber system a resonator). Comparing Figure 4.3 and Figure 4.6, the vibration absorption performance of the two methods are very similar. In order to best reveal the difference between the two methods, we choose another system for simulation in this section (Table 4.1 and Table 4.2). The value of the mp (the mass of the primary system excluding the absorber mass) is decreased to 0.75 kg assuming that the device is remanufactured lighter, the value of ba is decreased to 1.26 kg/s assuming that the Eddy current and Coulomb friction dissipations are minimized by making the housing of the absorber part from plastic and replacing the bushings with linear bearings, and the value of the ka is decreased to 200 N/m 1 . Figures 4.8 and 4.9 show the stability charts of the DR-DVA and PI-DVA of this system using the filters of Table 3.4. In order to make the simulations faster, the stability charts are plotted without considering the dynamics of the electromagnetic actuator2. The operational range of frequency for DR-DVA is 3.98-10.28 Hz and for PI-DVA, this range is 0-11.19 Hz, using the stability charts. As mentioned in the Introduction, the VE device is required to vibrate tissue with a range of frequencies in order to extract the mechanical properties of the tissue. This can be done by sweeping the frequency over the required range. Also, since the probe is hand-held, the amplitude of the excitation force can vary. Figure 4.10 shows the performance of DR-DVA and PI-DVA for various 1lt will be noted in chapter 6 that future work will focus on developing a light device with less friction. 2The simulation results can be generalized to the DR-DVA and PI-DVA with the electromagnetic actuator because it can be shown that the dynamics of the electromagnetic actuator affects equally the performance of both DR-DVA and PI-DVA. 4.1 Simulation Results 66 (b) - 0 . 3 -- 0 . 4 • - 0 . 5 -t(s) Ms) Figure 4.7: Control force of the DR-DVA and PI-DVA when the combined system is excited at 8.77 Hz. (a) DR-DVA. (b) PI-DVA. kinds of excitation that the hand-held device can experience. All simulations reveal the effect of the extra DR poles through system transient responses. The PI controller performs better in all cases. Table 4.3 compares the amount of vibration absorption achieved with DR-DVA and PI-DVA. 7 s is chosen for the comparison of passive DVA, DR-DVA and PI-DVA in parts (a), (b) and (d) because the system is still passing its transient response at this time and has not reached the steady state yet. An analysis of the DR-DVA dominant poles, similar to the Section 3.2.3.2, shows that the DR-DVA has large settling time if it is tuned a frequency close to the 3.98 Hz minimum operational frequency, say between 3.98 Hz and 5.5 Hz. The excitation frequencies are selected in this range in order to to best reveal the differences between DR-DVA and PI-DVA. At higher frequencies (say between 5.5 Hz and 10.28 Hz), when DR dominant poles become less influential (Figure 3.12), DR and PI perform quite similarly. Also, as mentioned before, 100% vibration absorption throughout the operational ranges (3.98-10.28 Hz for the DR-DVA and 0-11.19 Hz for the PI-DVA) is achieved with either method at the steady state. ma K (W0) ba (Ca) 0.2 kg 200 N/m (5 Hz) 1.26 kg/s (0.1) Table 4.1: Absorber system parameters. 4.1 Simulation Results 67 mp kp (ujp) Cp (Cp) 0.75 kg 600 N/m (1.9 Hz) 5 kg/s (0.1) Table 4.2: Primary system parameters. — Branch 1 .absorber Branch 2, absorber - Branch 1, combined Branch 2, combined 0.2 0.25 0.3 delay (s) Figure 4.8: Stability chart of the DR-DVA. The parameters of the absorber system and the combined system are according to Tables 4.1 and 4.2, and the lowpass filter of Table 3.4 is used. Looking at the lower graph, the operational range of frequency is 3.98-10.28 Hz. 4.1.5 The Effect of Coulomb Friction As mentioned before, Coulomb friction exists between the moving magnet and the housing and also between the bushings and the rods. It was impossible numerically to simulate the system with Coulomb friction because of its discontinuities at the times when xa — xp changes its direction. Therefore, the Coulomb friction model of Figure 4.11 (a) is approximated with the continuous friction force model of Figure 4.11 (b) where A = Ff/M and M is a suitably large number. Figures 4.12 shows the simulation results of PI-DVA considering the Coulomb friction between the absorber mass and the housing. A value of 100 kg/s is assumed for the M in these simulations. This value was obtained by increasing the M value until the results remain unchanged, ie. increasing 4.1 Simulation Results 68 absorber combined rji 1 1 1 1 1 1 -10 -8 -6 - 4 - 2 0 2 Ki(kg/s) Figure 4.9: Stability chart of the PI-DVA. The parameters of the absorber system and the combined system are according to Tables 4.1 and 4.2, and the lowpass and highpass filters of Table 3.4 are used. The asterics shows the point [Ki Kp]=[0 0] (passive DVA) and indicates the stable part of the chart. Looking at the lower graph, the operational range of frequency is 0-11.19 Hz. Excitation force Figure 4.10 (al) Figure 4.10 (bl) Figure 4.10 (cl) Figure 4.10 (dl) t 7 s 7 s 10 s 7 s Passive DVA 3 mm 3.05 mm 6 mm 7.5 mm DR-DVA 0.83 mm 0.7 mm 0.1 mm 3 mm PI-DVA 0.07 mm 0.01 mm 0.08 mm 0.7 mm Table 4.3: Amplitude of the hand-held device (primary system) with different vibration absorbers. Vibration of the device with passive DVA, DR-DVA and PI-DVA subjected to the forces depicted in Figure 4.10. M to 1000 kg/s has no effect on the results. The parameters of the absorber system and the combined system are as Table 3.2 and Table 3.3. In Figures 4.12 (a) the amount of Coulomb friction is 0.1 N. the equivalent linear damping required in the PI control law is calculated as 4Ff 1T(jJeX (4.1) 4.1 Simulation Results 69 0 2 4 6 8 10 0 2 4 6 8 10 0 2 4 6 8 10 I (a) I (8) 1 (s) Figure 4.10: Comparison of the vibration absorption performances of the DR-DVA and the PI-DVA. (al) Excitation frequency is changed from 4.8 Hz to 4.1 Hz at t=4 s. (a2) and (a3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (al). (bl) Excitation amplitude is changed from 1 N to 1.5 N at t=4 s with a constant frequency of 4.3 Hz. (b2) and (b3) Vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (bl). (cl) Excitation frequency is a constant 6 Hz for t=0-4 s and is decreased linearly from 6 Hz at t=4 s to 4 Hz at t=10 s. (c2) and (c3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (cl). (dl) Excitation amplitude is a constant 1 N for t=0-4 s and is increased linearly from 2 N at t=4 s to 8 N at t=10 s. The frequency of excitation is kept constant at 4.2 Hz.(d2) and (d3) Residual vibration of the primary system with DR-DVA and PI-DVA respectively, subjected to the excitation force of (dl). 4.1 Simulation Results 70 (a) (b) V F 1 f - A n 4 V Figure 4.11: Coulomb friction models (a) The discontinuous model, (b) The approximate continuous model. where cafiq is the equivalent viscous damping, Fj is the Coulomb friction, uie = 27T.5 rad/s is the excitation frequency and X is the absorber mass relative motion with respect to the probe (xa — xp). This equivalent viscous damping dissipates the same energy in each cycle as the Coulomb friction [62]. As we can see, vibration absorption is not 100% in the steady state. This residual motion of the primary mass in the steady state is a result of model error (approximating the nonlinear Coulomb friction with a linear viscous friction). Figures 4.12 (b) shows the simulation results of the displacement of the hand-held device when the Coulomb friction is increased to 0.2 N. This graph shows that the performance of the vibration absorber declines by increasing the amount of the Coulomb friction force. Figures 4.12 (c) is the Coulomb friction force obtained from the simulation, while the approximation of the Figure 4.11 (b) is assumed. The friction value is very close to the discontinuous model of the Figure 4.11 (a), verifying the approximation made. Figure 4.12: PI controlled vibration absorption with Coulomb friction (a) The value of the Coulomb friction is 0.1 N and the excitation frequency is 5 Hz. (b) The value of the Coulomb friction is 0.2 N and the excitation frequency is 5 Hz. (c) The friction force of the part (b). 4.1 Simulation Results 71 The performance of the DR-DVA with Coulomb friction can be examined using the same ap-proximation of Figure 4.11. The results show that the two methods have similar performance (about 2% difference) in the presence of the Coulomb friction. 4.1.5.1 Effect of Excitation Amplitude on the Amount of Vibration Absorption The transfer function for a linear system is independent of the input, a fact that does not hold true for nonlinear systems. Figures 4.13 shows how the amplitude of the excitation force can affect the absorption value achieved with the PI-DVA. The excitation force is a 10 Hz sinusoidal force and its amplitude changes between 0.1 N and 1 N. Figure 4.13 (a) shows the amplitude of probe vibration with and without active vibration absorption (amplitudes measured when the system reaches its steady state condition). As we can see, the amplitude of probe motion with the active vibration absorption is almost independent of the excitation amplitude, while amplitude of probe vibration expectedly increases linearly when the excitation force amplitude increases. Figure 4.13 (b) shows the absorption level achieved, with the PI-DVA at each amplitude. There is no dissipation except Coulomb friction, the amount of Coulomb friction is 0.05 N for all excitation amplitudes, and the nonlinear Coulomb friction damping is approximated using the linear viscous damping value from equation (4.1). Similar results (about 2% variations) were obtained using the DR-DVA. (•> 0.35; 100 0.3 W i t h o u t o b r o i p t k x . W i t h a b a o i p l l o n 0.25 0.1 0.05 0.2 0.6 0.8 Fe (N) Fe (N) Figure 4.13: The effect of the amplitude of the excitation on absorption achieved with the PI-DVA in the Coulomb friction damped system (a) Amplitude of probe motion with and without vibration absorption, (b) The percentage of the absorption level. 4.2 Experimental Results 72 4.2 Experimental Results Experimental results of the absorber system vibration are presented first. Then, the combined system behavior and vibration absorption results are presented. 4.2.1 Absorber System Both PI-DVA and DR-DVA methods make the stable absorber system marginally stable by plac-ing two poles of the characteristic equation of the absorber system on the imaginary axis. This marginal stability was shown in the simulation results (Figure 4.2 and Figure 4.4). This resonance behaviour of the absorber system is the key point to the vibration absorption (Section 2.1.2.1) and is investigated experimentally in this section. Results of a free vibration test on the actual device with an approximately 3 mm absorber mass initial displacement were shown before in Figure 3.13 (b). The system is clearly damped and vibrates approximately 1.5 cycles before dissipates all energy. Figures 4.14 (a) and (b) show the lowpass filtered acceleration results of the absorber system under a 10 Hz "free vibration" controlled with DR and PI respectively. We call this "free vibration" although it is a forced vibration since the actuator applies force on the absorber mass. This is because of the similarity between the DR-DVA and PI-DVA systems and a mass-spring system (Figure 3.2). Clearly the highly damped vibration of Figure 3.13 (b) is changed to slightly damped vibrations. Theoretically, the damping of the absorber system vibration should be zero when DR-DVA or PI-DVA are properly tuned. However because of uncertainties in the system and model errors, it is not possible to make the absorber system marginally stable. To obtain the best results, the absorber system properties and especially ba which has the most uncertainty (these properties are required for both DR-DVA and PI-DVA control laws) were changed manually in some interval around the values obtained by the parameter identification to accommodate some error in parameter identification. Clearly, increasing ba in the control law eventually makes the system unstable, as we can see in Figure 4.14 (c). Figure 4.15 (a) shows the "free vibration" experimental results of the absorber system tuned to 3 Hz controlled with PI controller. This frequency is far below the DR-DVA minimum operational frequency range and trying to tune the DR controller to this frequency makes the absorber system 4.2 Experimental Results 73 I (s) 1 (s) 1 (s) Figure 4.14: "Free vibration" experiments on the absorber system (a) DR controlled mass-spring-damper system tuned to resonate at 10 Hz. (b) PI controlled mass-spring-damper system tuned to resonate at 10 Hz. (c) Unstable DR controlled absorber system caused by choosing a value of ba that is too large. highly unstable. Again, the "free vibration" is slightly damped because of the uncertainties in the absorber parameters. 4.2.2 Combined System Several experiments were carried out on the probe to show the effectiveness of the active vibration absorber. A mechanical arm was built to hold the probe in order to perform repeatable and comparable experiments (Figure 4.16). The arm has a stiffness of 140 N/m, mimicking the stiffness of the human arm holding the probe. This value is selected according to the argument in Section 3.2.2.1. An OPTOTRAK 3020 (NDI, Waterloo, Canada) tracker with a measurement rate of 500 Hz is used to track the motion of the probe. The OPTOTRAK is only used to investigate the performance of the vibration absorber and is not used as a sensor in the vibration absorption control algorithm. Figure 4.15: Absorber mass 3 Hz "free vibration". PI controlled mass-spring-damper system tuned to resonate at 3 Hz. 4.2 Experimental Results 74 (a) (b) Figure 4.16: The hand-held device, (a) Held by hand, (b) Held by the mechanical arm. It was observed experimentally that more noise exists at lower frequencies, so a lowpass filter with a lower cutoff frequency was required to cancel the noise. The PI-DVA stability chart (Figure 3.8) shows that decreasing the lowpass filter cutoff frequency brings the two absorber system and the combined system stability curves closer to each other. It means that small changes in the Kp and Ki values that make the absorber system marginally stable moves the combined system to instability. These variations are inevitable especially with the uncertainties in the absorber model. Therefore the system is not robust at low frequencies. Experimental results also showed that although the PI-DVA can be tuned to frequencies as low as 3 Hz, the combined system can become unstable below 5 Hz. A possible explanation is that the range of stability is sensitive to changes in the absorber properties. This can be seen in Figure 3.8 where the crossover point of the two curves can vary greatly for small changes in the shape of the curves. This frequency is still below the minimum operational frequency range of the DR which was experimentally 8 Hz. 4.2 Experimental Results 75 Vibrating the tissue at very high frequencies, say 20 Hz, vibrates the probe at very low am-plitudes3. Experimentally, it was measured around 0.2 mm. Low amplitude vibrations are very difficult to cancel mainly because of the Coulomb friction: systems with Coulomb friction tend to become chaos when they vibrate with very small amplitude [11]. In the other words, during low excitations, the absorber mass vibrates with very small amplitudes and is essentially stuck to the housing. Therefore the vibration absorption experiments are performed up to 15 Hz, although the vibration absorber works in theory up to around 25 Hz. Vibration absorption experiments are performed in the 8 - 15 Hz range for the DR-DVA and in the 5 - 15 Hz range for the PI-DVA. Fo each experiment, ba is allowed to be varied from the value obtained by parameter identification to achieve the best performance. This is done because the friction varies slightly over time as the device is operated. Figure 4.17 shows the probe vibration before and after turning the active vibration absorber on. The probe is held by the mechanical arm in all experiments. The excitation frequencies are 8 Hz, 10 Hz, 12 Hz and 15 Hz for parts (a), (b), (c) and (d) respectively. These experiments show around 85% vibration absorption with either DR-DVA or PI-DVA, which is promising especially considering the deteriorating effects of Coulomb friction and Eddy currents. The performance of both methods is influenced by Coulomb friction and Eddy current nonlineaxities. As mentioned, the amount of Coulomb friction was observed to be variable from experiment to experiment. Therefore direct comparisons should be done with caution. Figure 4.18 shows a vibration absorption experiment with the PI-DVA at the frequency of 5 Hz, a frequency far below the stability range of the DR-DVA . To cancel the noise, it was required to decrease the cutoff frequency of the lowpass filter to 20 Hz. To see the performance of the probe to absorb forces with variable amplitude, two experiments with and without active vibration absorption are done. In these experiments, the tissue excitation amplitude is increased from 1.7 mm to 2 mm at t=2 s (Figure 4.19). The amplitude of vibration is increased from 0.73 mm to 1 mm without vibration absorption and from 0.11 mm to 0.14 mm with the vibration absorption. The experiment shows an increase in the attenuation level when 3The device held by hand can be approximated by a 1 DOF mass-spring-damper model (Figure 2.6). The natural frequency (and the resonance frequency) of this system is about 2 Hz using kp = 200 N/m and mp = 1 kg and 5 Hz using kp = 1000 N/m and mp = 1 kg. In the forced vibration of the 1 DOF mass-spring-damper systems, the amplitude of vibration is maximum at the resonance frequency and it decreases by increasing the frequency of the excitation force [62]. Therefore the amplitude of the vibration of the device is very low at 20 Hz. 4.2 Experimental Results 76 DR-DVA PI-DVA i d ) t ( s ) Figure 4.17: Vibration absorption experiments at four different frequencies, (al) Vibra-tion absorption with the DR-DVA at 8 Hz. (a2) Vibration absorption with the PI-DVA at 8 Hz. (bl) Vibration absorption with the DR-DVA at 10 Hz. (b2) Vibration absorption with the DR-DVA at 10 Hz. (cl) Vibration absorption with the DR-DVA at 12 Hz. (c2) Vibration absorption with the PI-DVA at 12 Hz. (dl) Vibration absorption with the DR-DVA at 15 Hz. (d2) Vibration absorption with the PI-DVA at 15 Hz. 4.2 Experimental Results 77 8 10 Figure 4.18: Vibration absorption experimental results with the PI-DVA at 5 Hz. the excitation force is increased, as predicted before because of the absorber system model errors (Figure 4.13). - Without absorption - With absorption \\ i H n \i ii H V; Hi Hi ii ii i; ii ii ii ii P P ^ liU' 0 5 \ U UII ii It. h \.\ HI \. i ii ii \ ii ii ii i; =i Without absorption - With absorption 2 t ( s ) Figure 4.19: Absorption performance of the DR-DVA and PI-DVA with variable am-plitude excitation force, (a) DR-DVA. (b) PI-DVA. The excitation amplitude is increased from 1.7 mm to 2 mm at t=1.3 s. Dashed and solid curves show the motion of the probe without and with the active vibration absorber respectively. Another requirement is the ability to cancel vibration during a swept-sine excitation. Figure 4.20 shows the results of two experiments performed with and without active vibration absorption (solid curve and dotted curve respectively). The first 2 seconds represent a constant amplitude, constant frequency excitation at 8 Hz. At t = 2 s the excitation frequency start to increase linearly from 8 Hz to 16 Hz in 5 s. The PI controller is used in the active vibration absorption. The absorption performance deteriorates by increasing excitation frequency mainly because the amplitude of probe vibration decreases with increasing frequency (see footnote 2, page 74 for an explanation of why the amplitude of vibration of the probe decreases by increasing excitation frequency.). This will result in performance degradation, as predicted in Section 4.1.5. > 4.2 Experimental Results 78 Figure 4.20: Absorption performance of the DR-DVA and the PI-DVA with variable frequency excitation force. The excitation frequency is 8 Hz for t=0-2 s and is increased linearly to 16 Hz in 5 s. Dashed and solid curves show the motion of the probe without and with the active vibration absorber respectively, (a) DR-DVA. (b) PI-DVA. The mechanical arm is only made to perform comparable and repeatable results. The ultimate goal is to use the device hand held and perform vibration absorption. In Figure 4.21, the probe is held by a human operator on human tissue (the forearm of another human subject). (•) (b) 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8 Figure 4.21: Probe vibration for experiments with a human operator on real tissue. The PI controller is used in this experiment and the vibration absorber is turned on at t=0.7 s. To show the non-constant behavior of the friction in the device, Figure 4.22 shows a second experiment with the same setup as Figure 4.21. An even better level of vibration absorption is achieved but it is not always repeatable. In this experiment, the vibration absorber decreases the peak to peak motion from 1.28 mm to 0.13 mm. Because the probe is held by hand and applies a slightly uneven pressure to the tissue, a net motion at a low 2 Hz frequency (compared to the 10 Hz vibration frequency) is seen. This is a rigid body motion where both the vibrating forearm and hand-held probe move together and no vibration absorption is required for this vibration. 4.3 Summary 79 0.81 -1 0 0.5 1.5 2 2.5 3 t ( S ) Figure 4.22: Probe vibration for experiments with a human operator on real tissue. The PI controller is used in this experiment and the vibration absorber is turned on at t=0.7 s. 4 . 3 S u m m a r y The simulation results verified that the operational frequency range of the PI-DVA is wider than that of the DR-DVA, considering either the absorber system or the combined system. They also showed that 100% vibration absorption in the steady states can be achieved with either the DR-DVA or the PI-DVA. We showed through simulation results that the DR-DVA and the PI-DVA require the same control effort (same control force and same absorber mass displacement) to cancel the vibrations. We also related the resonance behaviour of the two active DVA methods (Figure 4.5) and their vibration absorption performance (Figure 4.10) in this chapter. We showed through the simulation results that the DR-DVA is very slow at some frequencies, so it should not be operated at those frequencies. A study of the effect of Coulomb friction on vibration absorption performance is made with a simple approximation, depicted in Figure 4.11. We showed that with Coulomb friction and linear control it is not possible to achieve 100% vibration absorption, so some residual vibrations always exists. The amount of residual vibration increases with increasing friction force. We also showed that with Coulomb friction, it is more difficult to cancel small vibrations. This was done using simulations, and was verified by the experimental results later in the chapter. The experimental results of the resonance of the absorber mass and the vibration absorption performance of the active vibration absorber was provided in this chapter. About 85% vibration absorption was achieved with the active DVA, but it required manual slight re-tuning of the parameters of the absorber system to accomodate changes in friction. By replacing the bushings with linear bearings and increasing the air gap between the magnet and the housing to decrease the eddy currents (or by making the housing from plastic), the performance and repeatability of the 4.3 Summary 80 vibration absorber should be improved. The device was designed for a 5-25 Hz frequency range. However the experiments were performed in the 5-15 Hz because the amplitude of vibration of the device is very small at frequencies more than 15 Hz (see footnote 2, page 74). The 5-15 Hz (a ratio of 3 between the maximum and minimum frequency) range of the experiments performed in this work is encouraging. Compared to other works (1040-1270 Hz in [48], 650-750 Hz in [56] and 14-32.5 Hz in [17]), it is especially encouraging considering the fact that the other devices are not hand-held. Chapter 5 Vibro-Elastography with the Hand-Held Device As mentioned in the Introduction, parameter identification of the tissue mechanical properties in VE can be performed by two methods. The first method uses the equations of motion of a mass-spring-damper model and the second method uses transfer functions. Both were reported first by Turgay et al. [63]. In this chapter, we reconsider these methods when the displacement information is measured with respect to a moving reference since the vibration absorption is not expected to be perfect. 5.1 Mass-Spring-Damper Model The equation of motion of a mass-spring-damper model of the tissue is presented first. The equa-tion of motion is modified from [63] for a moving reference. The parameter identification with the displacement data from a moving reference is explained next. Simulation results verify the parameter identification method at the end. 5.1.1 Equation of Motion from a Stationary Reference A discrete one dimensional mass-spring-damper model of tissue is considered here. Figure 1.4 shows such a model, shown again in Figure 5.1 for convenience. 81 5.1 Mass-Spring-Damper Model 82 kn bn kn-l - 3 -bn-1 |—WWAM M N i—wvw- m n , ,,, _|—WWWH M 2 i—v w w -k2 •F(t) Figure 5.1: The discrete one dimensional visco-elastic tissue model, (from [63]) Writing the equation of motion for each mass, from a stationary reference, using Hooke's law and the definition of linear damper [63] gives B M m\ 0 0 • 0 0 m2 0 • • 0 0 0 m 3 • • 0 • 0 0 0 0 • • mn + bi -bi 0 0 -h ( 6 1 + 6 2 ) - 6 2 0 0 -b2 ( 6 2 + 6 3 ) - 6 3 0 0 0 0 + 6„_i f l X2 £3 it) •En v(t) ± 1 ± 2 ± 3 (t) + K x(t) kl -h 0 0 0 Xl / -h (h + k2) -k2 0 0 x2 0 0 -k2 (k2 + fc3) -k3 0 (*) = 0 0 0 0 0 kn—\ kn + kn—i •£•71 0 (t) = u(t) Ma(t) + Bv(t) + Kx(t) = u(t) (5.2) 5.1 Mass-Spring-Damper Model 83 where m; is the mass parameter of element i, bi is the damping parameter between rrii and mj_i, ki is the elasticity parameter between m; and mj_i, Xi(t) is the displacement measurement of the ith element in the model at time t, Vi(t) = ±i(t) is the velocity measurement (ddj) at time t, ai(t) = Xi(t) is the acceleration measurement ( ^ i ) at time t, and f(t) is the force measurement at time t. The definitions of the matrices M, B, K and x(t), x(t) and x(t) are clear from equation (5.1). 5.1.2 Equation of Motion from a Moving Reference As mentioned, equation (5.1) is valid provided that X j , i j and Xi are measured from a reference with zero acceleration. Assume that the position of each mass is measured from a moving reference with position XQS while the first mass is vibrated at a amplitude of xt (Figure 5.2, subscript CS and t referring to coordinate system and tissue respectively). Equation 5.2 should be modified to become M x l r ( t ) Xlr(t) x l r ( t ) mi x2r(t) X2r(t) X2r(t) m 2 XSr(t) + B X3r(t) + K X3r{t) + m 3 Xnrif) x n r ( t ) Xnrif) m n acs{t) = u(t) (5.3) vibrator plate Fixed boundary (bone) Figure 5.2: Imaging tissue with a moving ultrasound probe. 5.1 Mass-Spring-Damper Model 84 where Xir(t), Xir(t) and Xir(t) are the acceleration, the velocity and the displacement of mass i relative to the coordinate system CS, acs is the acceleration of the coordinate system CS, and M, B, K and u(t) are as before. 5.1.3 Parameter Identification Suppose the system of Figure 5.2 with unknown values of m ,^ Cj and fcj is excited by force f(t) and the position of each mass is measured with respect to the coordinate system CS whose acceleration is the unknown function acs(t) = xcs(t). The goal is to find unknowns m,, C j and ki for i = 1,..., n given Xir(i), iir(i) and Xir(t). The coordinate system CS is the ultrasound probe and the unknowns mi, ki and Ci are the elements of the mass-spring-damper model of tissue. Xir(t) values are known by applying correlation techniques on the ultrasound measurements (Section 1.4) and Xir(t) and ±ir(t) are obtained from Xir(t) through differentiation. Equation 5.3 can be easily converted to the form of equation (5.2) if the acceleration of the coordinate system CS, acs, is known, acs can be found if the position of a fixed point (with respect to a stationary reference) is available in coordinate system CS. This is possible if a rigid boundary (which contains the fixed point) of tissue is visible in the ultrasound image, which is not always possible. Moreover, measurements in the far field of an ultrasound image are less reliable because the resolution decreases in the far field. Therefore acs(t) is assumed to be unknown. Writing equation (5.3) at a single time instant i i by inserting the values for Xir(t{), Xir(t\) and xir(H) , there will be 3n +1 unknowns (m ,^ 6j and ki for i = 1,..., n and acs(t\)) and n equations. This is an under-constrained problem. By inserting m values of Xir(t), ±ir(t) and Xir(t) at m different time instances tj, j — 1,... ,m into equation (5.3), the number of equations increase to mn and number of unknowns increase to 3n + m (mi, bi and ki for i = 1,..., n and acs(tj) for j = 1,..., m). Rearranging equation (5.3), 5.1 Mass-Spring-Damper Model 85 one obtains ' <pT(tl)' «(*i) e = u(t2) . VT{tm) _ _ u(*m) _ * u (5.4) where u contains the excitation forces, $ contains rrii, i — 1,... ,n and the measurements of motion of each element (tissue region) at a certain time instant as a response to the applied excitation, and 9 is a vector of entries of M, B, K and acs(tj), j = 1... m. Specifically, ipT(tj), 9 and u(ti) are defined as (using Vi = ii and aj = Xi) 0 0 0 0 0 0 0 0 0 0 a3r 0 0 0 0 0 a3, . (n-l r 0 0 0 0 0 <r 3 3 v2 — v\ 0 0 0 0 v2 - v3 3 3 0 0 0 0 0 • vLi - v i - vi-i 0 0 0 0 KT - XCS — X2 3 3 x2 ~~ x \ 0 0 0 0 x\ — x\ x\ — x\ 0 0 0 0 0 T3 - x3n - Xi-1 0 0 0 0 xnr—XCS (5.5) ML 5.1 Mass-Spring-Damper Model 86 where superscript j refers to time tj, j = 1 • • • m and ML = 0 0 ••• mi 0 ••• 0 0 0 ••• m 2 0 ••• 0 , j = 1 • - - m (5.6) 0 0 ••- mn 0 ••• 0 m-j and mi m 2 6i hi ki k2 kn and U(tj) = f(tj) 0 0 j = 1 • • • m (5.7) 5.1 Mass-Spring-Damper Model 87 Note that in equation (5.5), irj — vj+1 is used instead of vfr — v^i+1y since they are equal (and the same for xj — xj+1). Two modifications are required for equation (5.5) and equation (5.7) to solve for 6. First, the components (2n, n) and (3n, n) of the matrix <p(tj) in equation (5.5) contain ±cs and xcs which are unknown. This can be solved by removing the equation of motion of the mass mn from the equation set of equations 5.3. A second modification is required since the excitation force, in the vector u, is unknown. This is done by removing the equation of motion of the mass mi. This equation will be used to find the excitation force, / . This modification means we are not able to identify the parameters in absolute terms, they can only be identified relative to either ra^, bi or fcj. Relative values can still be displayed as a useful image. Here, all values will be identified relative to ki by assuming ki = 1. With these two modifications, equation (5.4) is changed to " v , r ( * l ) " e' = u'(ti) u'(t2) . u'(tm) _ u' (5.8) where 5.1 Mass-Spring-Damper Model 88 o?2r 0 0 0 0 0 0 0 0 0 0 a\r 0 0 0 0 0 a(n-2)r 0 0 0 0 0 °(n-l)r j 3 3 0 0 0 0 vz - v4 7 7 vi ~ vi • 0 0 0 0 0 3 3 • <_ 2 - < _ i 3 3 0 0 0 0 "(n-l) - 4 — x33 x3 ~ 4 0 0 0 0 T3 - T3 J _ J 0 0 0 0 0 4 - 2 — 4-1 r3 - T3 x n - l xn-2 0 0 0 4-i - 4 4 - 4-1 j -v{ 0 0 0 0 (5.9) where 0 0 ••• m 2 0 ••• 0 0 0 ••• m 3 0 ••• 0 0 0 ••• m n 0 ••• 0 V J m—j j = 1 • • • m (5.10) and 5.1 Mass-Spring-Damper Model 89 m 2 m 3 b2 h k2 h k-n—1 bi a\s a2 and 0 0 , j = 1 • • • m (5.11) Note that <ff does not contain xcs o r xcs (position and speed of the coordinate system C75) and that u' does not contain / (the excitation force). The equation (5.8), §'8' = u', is a nonlinear matrix equation since both and 6' contain unknowns. The following section describes how this equation can be solved. 5.1.3.1 Parameter Identification Iterations Depending on the amount of vibration of the ultrasound probe, initial estimates for M , B and K can be obtained. For example, if xcs, XCS and xcs are much smaller than xt, xt and xt respectively, M, B and K can be estimated by neglecting the term containing acs m equation (5.2), as is done 5.1 Mass-Spring-Damper Model 90 by Turgay et al. [63]. In the following procedure, another method is used that does not require this assumption. 1. Assume a value for m,, i = 2 • • • n. This value can be calculated by assuming the mass distribution is uniform in the tissue (Figure 5.2) and the density of tissue is equal to the density of water, 1000 kg/m3. 2. Insert the nii, i — 2---n values into equation (5.10). Insert the obtained Mel value into equation (5.9) at time instance j, j = 1 • • -m. Stack all ^ (tj) as in equation (5.8) to make 3. solve equation (5.8) for 6' with either least square method or instrumental variables method [63]. 4. rrii values are known from the vector 6'. Go back to step 2 if the ratio of the norm of the difference between rrii obtained in this iteration and rrij obtained in the previous section divided by the norm of m, obtained in the previous section is larger than a tolerance value. I.e. go to step 2 if 6m > tol. where 0~m = E K current iteration ^i, previous iteration) i=2 W^i, previous iteration i=2 (5.12) In equation (5.12) 5m is related to rrij. Similarly, 5b or 5k can be defined by replacing rrii with bi or ki respectively. These parameters can also be used to terminate the iterations. 5.1.4 Simulation Results In this section, an L-MSD (linear mass-spring damper) chain is simulated to see the result of probe vibration on the calculated visco-elastic properties of tissue. The overall structure of the simulation is divided into two blocks as shown in Figure 5.3. The first block defines the system and its output for different excitation types. The system properties, such as the number of elements to be used 5.1 Mass-Spring-Damper Model 91 in the MSD model chain and the coefficients of springs, dampers and masses, are defined in this block. The outputs are the displacements from each element, shown as X(t). The second block solves the parameter identification equation = u' (equation (5.8)). The displacements X(t) are filtered by the derivative and double derivative filters to obtain approx-imations for X(t) and X(t). A lowpass filter (( a + a )2 ) is used in the derivative filters to cancel high frequency noise in practical implementation. The cutoff frequency of the lowpass filters (a) is chosen to be same as the frequency that the white noise excitation is band-limited to. The L-MSD parameters 6 are then solved and displayed. The simulation environment is implemented by Simulink and MATLAB. The simulation parameters axe chosen to mimic the behaviour of real tissue as closely as possible. The MSD model coefficients are set to the estimated properties a 4 cm long section of human breast fat. The elasticity of human breast fat is given as 20 kPa in [63]. The density of fat is close to 1000 kg/m3. The stiffness and density of a cylinder (or cube) with cross section A, length /, elasticity E and density p is ^ and pAl respectively. Using these two relations and reasoning in [63], the following relation is derived: fc, F, 1 (5.13) ki mi E 1 p A /2 where Al is the distance between two neighboring elements of the L-MSD chain at rest. If a 5-element L-MSD chain is used to simulate a tissue sample having an axial length of 4 cm, then A/ is 8 mm, and fci/m; is 3 • 105. The ki/bi ratio is set to around 500 as is reported in [37]. This may not be an exact value for human tissue but will likely have the correct order of magnitude. Block : F(t) Band-limited white noise "•rT r^T~i^ = | X(t) Block 2 (s+ay ML (s+ay ML ML Solve 4>e = u Figure 5.3: Overall structure of the simulation environment, (from [63]) 5.1 Mass-Spring-Damper Model 92 According to these ratios for rrii, bi, and fcj presented above, the values shown in Figure 5.4 are assigned to the elements of the L-MSD chain. The excitation for the simulations is chosen to be white noise band limited to 15 Hz. This excitation is used instead of a sweep-sine to simplify the simulation. The resultant displacement of three of the six mass elements are shown in Figure 5.5 (a). The displacement of the ultrasound probe is assumed to be a 15 Hz sinusoidal vibration with the amplitude of 1 mm (Figure 5.5 (b)). Assuming that no vibration absorption exists, the amplitude of vibration of the probe, xcs, is chosen to be close to the amplitude of vibration of the skin surface, xt-The simulation results for m, b and fc, given that the ultrasound probe is stationary, are shown in Figures 5.6 (a), (b) and (c) respectively. Figures 5.6 (d), (e) and (f) show the simulation results for m, b and fc assuming that the probe moves as displayed in the Figure 5.5 (b). The initial guess for rrii is assumed to be 8 • 10 - 5 for all elements assuming a uniform mass distribution. The final values of fc (Figure 5.6 (f)) converge to the values obtained with a stationary coordinate system (Figure 5.6 (c)). However, the final values of 6 (Figure 5.6 (e)) and m ((Figure 5.6 (d)) do not perfectly converge to the values obtained with a stationary coordinate system (Figures 5.6 (b) and (a)). This is caused by the lowpass filter used to obtain velocity and acceleration from the displacement data. Figure 5.7 shows the how the variables Sm, 8b and Sk approach zero by each iteration. Figures 5.8 and 5.9 show estimated parameters and the convergence rate of the same system subject to the same 15 Hz frequency probe motion but with 0.2 mm amplitude. The amplitude of vibration of the coordinate system is decreased to 0.2 mm assuming that vibration cancellation on the motion of the coordinate system is performed. Comparing Figures 5.6 and 5.7 to 5.8 and 5.9, we can say that by decreasing the amplitude of the vibration of the probe from 1 mm to 0.2 mm, the number of iterations to converge to the solution decreases. For example, for a tolerance value of tol. = 4.5 • 10~3 for the 5m, the number of iterations to solve equation (5.8) decreases from 21 to 4. 5.2 Transfer Function Method 93 element index element index element index Figure 5.4: Actual coefficients of the MSD model, (a): rrii values, (b): bi values, (c): ki values (a) (b) 1 <s) I (s) Figure 5.5: Tissue displacement and probe displacement, (a): Tissue displacement at the 1 s t, 3rd and 5th nodes of the 6 node MSD model, (b); Displacement assumed for the ultrasound probe. 5.2 Transfer Function Method As mentioned in the Introduction, an alternative method for obtaining elasticity parameters is the transfer function method [63]. We discuss how this method is modified when the displacement measurements are obtained from a moving coordinate system. The first elastograms obtained from a gelatin phantom are provided next. The displacement data are captured from a moving ultra-sound probe, i.e. the coordinate system, with and without vibration absorption. The mechanical properties calculated from the displacement data are provided next. 5.2 Transfer Function Method 94 element index element index element index element index element index element index Figure 5.6: Estimated parameters with 1 mm amplitude xcs- (a), (b) and (c): rrii, bi and ki respectively. No motion for the ultrasound probe is considered, xcs = 0. (d), (e) and (f): rrii, bi and ki respectively at different iterations. The probe is supposed to vibrate as depicted in Figure 5.5 (b). Figure 5.7: Convergence variable for 1 mm amplitude xcs- The variable 5m is defined by equation (5.12), and 5b and 5k are defined similarly for b and k. 5.2 Transfer Function Method 95 1 1 1 8 t iter. - - - 2 n d iter. 5 • 3" 1 iter. — - a ' i ter . 2 V 5 rjl . . element index element index element index Figure 5.8: Estimated parameters for 0.2 mm amplitude xcs- (a), (b) and (c): rrii, h and ki respectively at different iterations. The probe is supposed to vibrate as depicted in Figure 5.5 (b) but with 0.2 mm aplitude. 0.2 0.18 0.16 0.14 0.12 «= 0.1 0.08 0.06 0.04 0.02 0 1 4 5 iteration no. Figure 5.9: Described procedure's convergence for 0.2 mm amplitude xcs- The variable 6m is defined by equation (5.12), and 5b and 5k are defined similarly for b and fc. 5.2.1 Calculating Elasticity from Transfer Functions 5.2.1.1 Transfer Functions from a Stationary Reference To find the elasticity from the displacement data using the transfer function method, low frequency components of the displacements are used [63]. At low frequencies, the magnitude of the transfer function (Figure 1.6) between two points is only a function of elasticity between two points and not density or viscosity. Therefore the local elasticity can be simply obtained from the low frequency values of the transfer function [63]. 5.2 Transfer Function Method 96 5.2.1.2 Transfer Functions from a Moving Reference By adding the first and fourth term on the left hand side of equation (5.3), one has xi(t) ±ir(t) xlr(t) x2(t) X2r(t) x2r(t) M x3(t) + B X3r(t) + K X3r(t) xn(t) xnr(t) By subtracting the displacement, velocity and acceleration of the mass m„ from the displacement, velocity and acceleration of the each mass rrii, i = 1 • • • n — 1, we wrtite equation (5.14) with respect to the coordinate system of the mass mn: xi(t) - Xn(t) x2(t) - xn(t) M n _ i x3{t) - xn{t) Xyi—l - xn(t) + Bn. x\{t) - xn(t) X2(t) - Xn(t) x3(t) ~ ±n(t) i „ _ l - Xn(t)(t) + K, n-l Xl(t) - Xn(t) X2(t) - Xn(t) x3r(t) - xn(t) xn-i(t) -xn(t) = u(t) (5.15) where Mn-i, Bn_i and Kn-i are (n— 1) x (n — 1) square matrices obtained by eliminating the last row and column of matrices M, B and K defined in equation (5.1). Note that in equation (5.15), ii(t)—xn(t) is used instead of xir(t) — xrn(t) since they are equal (and the same for Xi(t)—xn(t)). At low frequencies, rrii(xi(t)—xn(t)) << bi(xi(t)—xn(t)) « ki(xi(t)—xn(t)) (or equivalently m^w2 << biU « h where LO is a low frequency value). Therefore the first two terms in equation (5.15), containing terms M n _ i and Bn-i can be neglected. The simplified equation can be solved for the elasticity values using transfer functions as explained in [63]. Here, we wrote the equation of motion from the coordinate system of the last mass in the L-MSD chain (m„). However any other mass in the L-MSD chain can be used as the coordinate system. In fact, mass elements closer to the ultrasound probe are preferred over the mass elements 5.2 Transfer Function Method 97 farther from the ultrasound probe as the coordinate system since far field data in an ultrasound image has more error. 5.2.2 Simulation Results A L-MSD chain is considered here, with the mass, viscosity and elasticity values depicted in Figure 5.4. The excitation for the simulation is the same as before, a white noise band limited to 15 Hz. The displacement data X(t) are then calculated (Figure 5.3) and the elasticity values are calculated from the low frequency components of the transfer functions. It is now clear that the motion of the probe has no effect on the identification of the elasticity values. Figure 5.10 (a) shows the elasticity values obtained assuming the probe is stationary. Figure 5.10 (b) shows the elasticity values obtained assuming the probe has a 1 Hz sinusoidal vibration with the amplitude of 1 mm1. As the results indicate,the same answers are obtained with or without probe vibration. ( » (b) element index element index Figure 5.10: Estimated elasticity using the transfer function method, (a) Estimated ki values. No motion for the ultrasound probe is considered, (b) Estimated ki values. The probe is assumed to vibrate with a 1 Hz sinusoidal at 1 mm amplitude. Since only the low frequency components of the transfer functions are considered for calculating elasticity, the frequency of the vibration of the probe is decreased from the 15 Hz value of Figure 5.5 to 1 Hz. In the other words, to see the effect of probe vibration its frequency is decreased to 1 Hz because high frequency components of the displacement data are not considered in the elasticity identification 5.2 Transfer Function Method 98 5.2.3 Experimental Results A B-mode ultrasound image and two real static elastograms obtained from a three layered gelatin phantom are shown in Figure 5.11 2 . The middle layer of the phantom is stiffer than the other layers but exhibits the same acoustic properties. The stiff layer is not visible in the B-mode image but it is detectable in the traditional static elastograms. The static elastograms ((b) and (c)) are obtained by applying a 0.01 compression to the phantom slowly and rapidly acquiring RF data during the compression, to ensure high correlation coefficients between consecutive echo frames, as explained in the Introduction. The time domain cross-correlation with prior estimates [67] is used for calculating local displacements of the phantom. Best-case refers to the image with the highest contrast to noise ratio obtained over a set of repeated experiments, and worst-case refers to the image with lowest contrast to noise ratio. Figure 5.11: B-mode image and static elastograms of a three layered gelatin phantom. (a) B-mode image, (b) Best-case example of static elastogram. (c) Worst-case static elastogram. The geometric locations of the boundaries between layers are indicated with triangles. Figures 5.12 (a) and (b) show two vibro-elastograms obtained with the hand-held device, mounted on the mechanical arm, using a 5 Hz excitation without and with PI-DVA respectively. The transfer function method is used to calculate the elasticity. The accuracy of the obtained local displacements is reduced by probe vibration, resulting in a poor vibro-elastogram (Figure 5.12 (a)). A clearer image is produced with vibration absorption (Figure 5.12 (b)). The improvement in performance is mainly due to an improvement in the correlation based motion tracking algorithm that provides the displacement measurements. If one window is not correctly matched, subsequent 2The motion tracking software was obtained from Reza Zahiri, MASc student under the supervision of Prof. Salcudean. The motion tracking algorithm is currently under development and validation. 5.2 Transfer Function Method 99 windows are also mismatched, resulting in an incorrect set of displacement measurements for a whole column in an image. A 5 Hz excitation is chosen because it allows the vibration isolation to be demonstrated while maintaining the ability to obtain displacement measurements. The frame rate of the ultrasound currently prohibits higher frequency tests. This may be overcome in future tests with modifications to the ultrasound system. Currently, the motion tracking works well bellow 5 Hz. As explained in Section 4.2.2, the PI-DVA is theoretically capable of being tuned below 5 Hz. But the variations in the friction level produce instability. Even at 5 Hz, the PI-DVA performance is not as good as higher frequencies. In the future, the bushings guiding the absorber mass will be replaced by linear bearings. This will reduce friction and should make the system more stable and reduce the minimum operational frequency, and improve the level of vibration absorption at 5 Hz. To get an idea of the potential result of such changes, a vibro-elastogram is obtained of the same phantom at 5 Hz after clamping the probe fixed (Figure 5.13). This is an ad hoc test meant only to give insight into the possibility of further improving image quality. Figure 5.12: Vibro-elastograms of the three layered gelatin phantom obtained with the hand-held device mounted on the mechanical arm. (a) Vibro-elastogram at 5 Hz, no vibration absorption is used, (b) Vibro-elastogram at 5 Hz, vibration of the probe cancelled with PI-DVA. The geometric locations of the boundaries between layers are indicated with triangles. 5.2.4 Summary The equation of motion of the mass-spring-damper model of tissue from a moving reference was presented. We showed that solution of a nonlinear set of equations is required in order to solve Figure 5.13: Vibro-elastogram of the three layered gelatin phantom obtained with an approximately fixed probe at 5 Hz excitation. The geometric locations of the boundaries between layers are indicated with triangles. for the parameters of the mass-spring-damper model. An iterative approach to solve the inverse problem was presented. The simulation results show that the number of iterations decreases by decreasing the amplitude of the vibration of the moving reference. It was shown that no major modification from the method presented in [63] is required to calculate elasticity from the transfer functions when the ultrasound probe is not stationary. The first vibro-elastograms are obtained from a three layered phantom with the hand-held device. Elastograms obtained with and without probe motion show that calculating local displacements from the RF data becomes easier with small probe motion. Therefore vibration absorption can be used to enhance vibro-elastogram quality. Of course, vibration absorption also reduces operator fatigue. Chapter 6 Conclusions and Future Work 6.1 Conclusions A hand-held device that includes an active vibration absorber was designed and built. Active vibration absorption was selected among different methods that can be used to reduce the unwanted vibration. A new control method, PI-DVA, was exploited in the design along with an established method, DR-DVA. In most previous papers on DR-DVA, the emphasis was on the level of vibration absorption, and frequency range. Here for VE, we also focus on tuning speed, and obtaining the largest frequency range. Two interesting features of the PI-DVA are its wide frequency operational range and fast vibra-tion absorption, compared with the DR-DVA. The operational range of DR-DVA and PI-DVA was obtained using stability charts, and was verified later by simulation and experimental results. The speed of the vibration absorption of the two methods was compared by looking at the poles of the characteristic equations, and was verified by simulation results later. These two features suggest that PI-DVA is a suitable control method for the vibration absorption of the hand-held device. Guidelines for designing an active vibration absorption device were explained. Such a design should consider physical limits such as the allowed absorber mass and the possible stroke and force limits of the actuator. It should also satisfy the required operational range and absorption speed. Some tolerance should be considered for the deteriorating effects of nonlinearities in the system. In this device, Coulomb friction, Eddy currents, actuator saturation, and friction losses in the spring 101 6.2 Future Work 102 clamps are amongst these effects. Nonlinearities and errors in the models reduce the amount of vibration absorption, as well as produce a smaller operational frequency range. Identifying the tissue properties by VE was also shown to benefit from vibration absorption. Obtaining displacement measurements from RF ultrasound data is more difficult with large residual motions of the probe. Practical limitations of the current ultrasound measurement device limit operation to 5 Hz and lower. A re-design of the device is needed to reduce the deleterious effects of the nonlinearities of the low frequencies. 6.2 Future Work A new design for the device is envisioned. It will use smaller actuator, to reduce weight, a new shell design to minimize Eddy currents, and bearings to guide the moving parts with low friction. These changes should improve the performance of the device, especially at low frequencies. Higher speed motion tracking algorithms will also be explored, taking into consideration the residual motion of the probe with respect to the skin. One of the main issues of an active DVA design is the tradeoff between the absorber mass and its displacement. A large absorber mass results in a heavier device, but smaller displacement of the absorber mass. Conversely a small mass requires a large displacement. Yet large displacement actuators are mainly available from manufacturers in larger sizes, so choosing this approach still does not result in a very light system. Future work will continue to investigate this tradeoff with the aim of minimizing overall mass and size. The control parameters were calculated for the experiments using the nominal values of the properties of the absorber system. In practice, these properties may change over time, especially friction. An automatic tuning algorithm can be implemented. This algorithm may require addi-tional sensors (say the acceleration of the hand-held device). To cover a sufficient frequency range for the vibro-elastography, the excitation frequency is swept linearly in this work. Currently, no feed-forward scheme that uses the information about the future frequency of excitation is exploited in the control method. By introducing such a feed-forward scheme in the control algorithm, the vibration absorption performance may be improved. 6.2 Future Work 103 One of the main advantages of the vibro-elastography technique over other elastography tech-niques is its potential for calculating viscosity and density. Vibro-elastography experiments on different phantoms with different layers that exhibit different viscosity and different density values can be performed and the viscosity and density values can be calculated from the equations of motion method. Extracting viscosity and density from the transfer functions can be another area of research. The experimental results can then be compared to the experimental results obtained with the equations of motion method. This work is part of a larger effort in the lab to develop VE, and includes validation on phantoms, clinical tests, motion tracking, software and hardware acceleration, and 2D and 3D extensions of the algorithm. I will collaborate with others to combine systems, theory, and software for the development of clinical applications. Bibliography [1] R.S. Adler, J.M. Rubin, P.H. Bland, and P.L. Carson. Quantitative tissue motion analysis of digitized M-mode images: gestational differences of fetal lung. Ultrasound in Medicine and Biology, 16:561-569, 1990. [2] H. Andreas, A.H. von Flotow, A. Beard, and Bailey D. Optimum delayed feedback vibration absorber for MDOF mechanical structures. Proceedings - National Conference on Noise Control Eng., Progress in Noise Control for Industry, pages 437-54, May 1994. [3] T. Asami and O. Nishihara. H2 optimization of the three-element type dynamic vibration absorbers. Transactions of the ASME. Journal of Vibration and Acoustics, 124(4):583-92, October 2002. [4] S. A. Austin. The vibration damping effect of an electrorheological fluid. ASME journal of vibration and acoustic, 115(1):136-140, April 1993. [5] P.E. Barbone and N.H. Gokhale. Elastic modulus imaging: on the uniqueness and nonunique-ness of the elastography inverse problem in two dimensions. Inverse Problems, 20(l):183-96, January 2004. [6] G. A. Becus, C. Y. Lui, V. B. Venkayya, and V. A. Tischler. Simultaneous structural and con-trol optimization via linear quadratic regulator eigenstructure assignment. NASA Conference Publication, pages 225-232, 1987. [7] S.P. Bhattacharyya, H. Chapellat, and L.H. Keel. Robust Control: The Parametric Approach. Prentice Hall, 1995. [8] A. Bruner, W.K. Belvin, L. Horta, and J.N. Juang. Active vibration absorber for the CSI evolutionary model. Design and experimental results. Proceedings of the 32nd AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, page 2928, 1991. [9] E.J. Chen, W. Novakofski, K. Jenkins, and W.D. O'Brien. Young modulus measurements of soft tissue with application to elasticity imaging . IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 43(l):191-94, January 1996. [10] W. Chen, X.G. Wang, C. Sun, F. Devine, and C.W. De Silva. Active vibration control with state feedback in woodcutting. Journal of Vibration and Control, 9(6):645-64, June 2003. [11] G. Cheng and J.W. Zu. Dynamics of a dry friction oscillator under two-frequency excitations. Journal of Sound and Vibration, 275:591-603, 2004. 104 BIBLIOGRAPHY 105 C L . de Korte, A.F.W. van der Steen, E. Ignacio, and E.I. et al. Cespedes. Characterization of plaque components and vulnerability with inter vascular ultrasound elastography. Physics in Medicine and Biology, 45:1465-75, June 2000. R.L. Dickonson and CR. Hill. Measurement of soft tissue motion using correlation between A-scans. Ultrasound in Medicine and Biology, 8:263-71, 1982. BEI KIMCO Magnetic Division, http://www.beikimco.com/products/vca.php. World Wide Web Page. T.G. Duclos. Design devices using electrorheological fluids. SAE paper number 881134, 1988. E.M. ElBeheiry, D.C. Karnopp, M.E. ElAraby, and A.M. AbdElraaouf. Advanced grouned vehicle suspension systems-a classified bibliography. Vehicle System Dynamics, 24(3):231-258, April 1995. H. Elmali, M. Renzulli, and N. Olgac. Experimental comparison of delayed resonator and PD controlld vibration absorbers using electromagnetic actuators. Journal Dynamic Systems, Measurement and Control, 122:514-20, September 2000. S.Y. Emelianov, M.A. Lubinski, and A.R. et al. Skovoroda. Reconstructive ultrasound elas-ticity imaging for renal transparent diagnosis: kidney ex-vivo results. Ultrasonic Imaging, 22:178-94, July 2000. E. Esmailzadeh and F. Fahimi. Optimal adaptive active suspensions for a full car model. Vehicle System Dynamics, 27(2):89-107, February 1997. D. Filipovic and N. Olgac. Delayed resonator with speed feedback - deign and performance analysis . Mechatronics, 12:393-413, 2002. D. Flipovic and D. Schroder. Bandpass vibration absorber. Journal of sound and vibration, 214(3):553-66, July 1998. M. A. Franchek, M. W. Ryan, and R. J. Bernhard. Adaptive passive vibration control. Journal of Sound and Vibration, 189(5):565-585, 1995. Y . C Fung. Biomechanics, mechanical properties of living tissues. Springer-Verlag, 1993. B.S. Garra, E.I. Cespedes, and J. et al. Ophir. Elastography of breast lesions: initial clinical results. Radiology, 202:79-86, 1997. J.M. Greenleaf, M. Fatemi, and M. Insana. Selected methods for imaging elastic properties of biological tissues. Annual Review of Biomedical Engineering, 5:57-78, April 2003. E. Guglielmino, K. Edge, and R. Ghiglizza. On the control of the friction force. Meccanica, 39:395-406, 2004. B.J. Hamrock, Schmid. S.R., B.O. Jacobson, B. Hamrock, and S. Schmid. Fundumentals of machine elements. 2nd Edition, McGraw-Hill, 2004. A.P.G. Hoeks, P.J. Brands, J.M. Willigers, and R.S. Reneman. Non-invasive measurement of mechanical properties of arteries in health and disease . Proceedings of the Institution of Mechanical Engineers, Part H (Journal of Engineering in Medicine), 213:195-202, 1999. BIBLIOGRAPHY 106 [29] D.V. Horton and D.A. Crolla. Theoretical analysis of a semi-active suspension fitted to off-road. Vehicle System Dynamics, 15:351-372, February 1986. [30] D. Inman. Engineering Vibration. New York: John Wiley and Sons, Second Edition, 2001. [31] R.G. Jacquot. Optimal dynamic vibration absorbers for general beam. Journal of Sound and Vibration, 60:535-542, 1978. [32] N. Jalili. A comparative study and analysis of semi-active vibration-control systems. Journal of Vibration and Acoustics, Transactions of the ASME, 124(4):593-605, October 2002. [33] N. Jalili and N. Olgac. Optimum delayed feedback vibration absorber for MDOF mechanical structures. Proceedings of the IEEE Conference on Decision and Control, 4:4734-39, July 1998. [34] F. Kallel and M. Bertrand. Tissue elasticity reconstruction using perturbation method. IEEE Transactions on Medical Imaging, 15:299-313, June 1996. [35] B.G. Korenev and L.M. Reznikov. Dynamic Vibration Absorbers: Theory and Technical Applications. New York: John Wiley and Sons, 1993. [36] T.A. Krouskop, T.M. Wheeler, M. Kallel, B.S. Garra, and T.J. Hall. Elastic moduli of breast and prostate under compression. Ultrasonic Imaging, 20:260-74, 1998. [37] S.A. et al. Kruse. Tissue characterization using magnetic resonance elastography: preliminary results. Physics in Medicine and Biology, 45:1579-90, June 2000. [38] Y. Kuroda, T. Ura, and S. Morishita. Damping-controllable dynamic damper with neural network based adaptive control system. IEEE Int Jt Conf Neural Networks, pages 1807-12, 1991. [39] R.M. Lerner, S.R. Huang, and K.J. Parker. 'Sonoelasticity' images derived from ultrasound signals in mechanically vibrated tissues. Ultrasound in Medicine and Biology, 16(3):231-9, 1990. [40] R.M. Lerner and K.J. Parker. Sonoelasticity in ultrasonic tissue characterization and echo-graphic imaging. Eur. Commun. Worksh.), October 1987. [41] S.F. Levinson, M. Shinagawa, and T Sato. Sonoelastic determination of human skeletal muscle elasticity. Journal of Biomechanics, 28:1145-11154, October 1995. [42] G.S Lin. US Pat. 6,068,597, May 2000. [43] G.S. Lin. Vibrational resonance ultrasonic Doppler spectrometer and imager. Patent number 6068597, 2000. [44] M.A. Lubinski, S.Y. Emelianov, and M. O'Donnell. Speckle tracking methods for ultrasonic elasticity imaging using short time correlation. IEEE Transactions on Ultrasonics, Ferro-electrics and Frequency Control, 46:82-96, January 1999. [45] D.N. Manikanahalli and M.J. Crocker. Vibration absorbers for hysterically damped mass-loaded beams. Transactions of ASME Journal of Vibration and Acoustics, 113:116-122, 1991. BIBLIOGRAPHY 107 [46] H. Nishimura, K. Nonami, W. Cui, and A. Shiba. Hoo control of multi degree of freedom structures by hybrid dynamic vibration absorber (experimental consideration of robustness and control performance. Transactions of the Japan Society of Mechanical Engineers, 59:714-20, 1993. [47] K. Ogata. Modern control engineering. J^th edition, Prentice Hall, New Jersey, 2002. [48] N. Olgac, H. Elmali, M. Hosek, and M. Renzulli. Active vibration control of distributed systems using delayed resonant with acceleration feedback. Journal of Dynamic Systems, Measurement and Control, 119:380-89, September 1997. [49] N. Olgac and B.T. Holm-Hansen. Novel active vibration absorption technique: delayed res-onator. Journal of Sound and Vibration, 176(1):93-104, September 1994. [50] N. Olgac, M. Hosek, H. Elmali, and Renzulli M. Implementation of DR and DFFDR techniques on distributed systems. Active Control of Vibration and Noise, 93:83-90, 1996. [51] J. Ophir, I. Cespedes, H. Ponnekanti, Y. Yazdi, and X. Li. Elastography: a quantitative method for imaging the elasticity of biological tissues. Ultrasonic Imaging, 13:111-134, April 1991. [52] H.N. Ozguven and B. Candir. Suppressing the first and second resonance of beam by dynamic vibration absorbers. Journal of Sound and Vibration, 111:377-390, 1986. [53] M.K. Petek, D.L. Romstadt, M.B. Lizell, and T.R. Weyenberg. Demonstration of an au-tomotive semi-active suspension using electro-rheological fluid. SAE paper number 950586, 1995. [54] C. Philips and R. Harbor. Feedback control systems . Prentice Hall, 1998. [55] M.J. Reeves, p.L. Newcomb, and P.M. Marcus. Determinants of breast cancer detection among Wisconsin (United States) women. Cancer-Cause-Control, 6:103-111, 1995. [56] M. Renzulli, R. Ghosh-Roy, and N. Olgac. Robust control of the delayed resonator vibration absorber. IEEE Transactions on Control Systems Technology, 7(6):683-91, November 1999. [57] R.S. Sharp. Variable geometry active suspension for cars. Computing and Control Engineering Journal, pages 217-222, October 1998. [58] T.T. Soong. Active structural control: Theory and practice. New York: John Wiley, 1990. [59] B.F. Jr. Spencer, S.J. Dyke, and H.S. Deoskar. Benchmark problems in structural control-part I: Active mass driver system. Earthquake Engineering and Structural Dynamics, 27(11):1127-39, November 1998. [60] A. Stribersky, H. Muller, and B. Rath. The development of an integrated suspension control technology for passenger trains. Journal of Rail and Rapid Transit, 212(Fl):33-42, 1998. [61] J.Q. Sun, M.R. Jolly, and M.A. Norris. passive, adaptive and active tuned vibration absorbers-a survey. . Journal of Vibration and Acoustics, 117:234-42, June 1995. BIBLIOGRAPHY 108 [62] W.T. Thomson and M.D. Dahleh. Theory of vibration with applications. 5th edition, Prentice Hall, New Jersey, 1998. [63] E. Turgay. Imaging Visco-Elastic Properties of Soft Tissue with Ultrasound. PhD thesis, University of British Columbia, 2004. [64] E. Turgay, S. Macintosh, R. Rohling, and S. Salcudean. Parameter identification of tissue lumped models based on sequences of ultrasonic strain images. Second International Confer-ence on Ultrasonic Measurement and Imaging of Tissue Elasticity., page 9, 2003. [65] M. Turgay, Salcudean S.E., and R.N. Rohling. Method for imaging the mechanical properties of tissue. Patent number 2457376, 2004. [66] Y. Yamakoshi, J. Sato, and T. Sato. Ultrasonic imaging of internal vibration of soft tissue under forced vibration. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 37(2):45-53, March 1990. [67] I. Youn and A. Hac. Semi-active suspension with adaptive capability. Journal of Sound and Vibration, 180(3) :475-492, 1995. [68] R. Zahiri and S. Salcudean. Time domain cross correlation with prior estimates. Third Inter-national Conference on Ultrasonic Measurement and Imaging of Tissue Elasticity., page 63, 2004. [69] W. Zhang, H. Matsuhisha, and S. Sato. Dynamic absorber for the general multi-degree-of-freedom vibration system. Transactions of Japan Society of Mechanical Engineers, Part C, 56:921-27, 1990. Appendix A Voice Coil Actuator Selection A . l Theory of Voice Coil Actuators Voice coil actuators utilize a permanent magnet field and coil winding to produce a force that is proportional to the current applied to the coil. These non-commutated electromagnetic devices can be used in linear and rotary motion. Originally used in radio laud speakers, voice coil actuators axe gaining popularity in applications where proportional or tight servo control is a necessity. The force generated by a voice coil actuator, F, is given by the Lorentz Force law [14]: F = kBLIN (A- l ) where k is a constant, B is the magnetic flux density, L is the length of the conductor, I is the current and N is the total number of the conductors (number of windings). Furthermore, a conductor moving through a magnetic field will have a voltage induced across the conductor. The magnitude of the voltage, E, can be found as [14] E = kBLvN (A-2) where v is the speed of the conductor and k, B , L and N were defined for equation ( A - l ) . 109 A.2 Sizing Voice Coil Actuators 110 A.2 Sizing Voice Coil Actuators We first select the first actuator for vibrating the tissue. Three parameters will determine actuator selection: 1. Peak force requirement (Fp). 2. RMS (root mean square) force requirement (FRMS)-3. Total stroke or move distance (£>). Peak force, is the sum of the force due to load, Fi, friction, Ff, and the acceleration of the mass, Fm [14]: FP = Fl+Ff + Fm (A-3) Looking at the separate components, the force due to the load is the force acting directly against the actuator all the times. The force due to friction is determined by the mechanical configuration of the complete motion assembly and includes such factors as bearings, linkages and surface to surface contacts. The force due to the acceleration of the mass is the product of load (including actuator coil) mass, mtotai, and the acceleration of the load, a [14]: Fm = mtotaia (A-4) mtotal =mi+mc (A-5) where m/ is the load mass and mc is the mass of the coil. The maximum force requirement happens at the maximum operation frequency. Since it may be required later to vibrate tissue at higher frequencies, the maximum frequency of excitation is assumed to be 40 Hz. The required peak force is 16.2 N for a ±0.002m stroke, assuming FL = 0, Ff = 1 N, mi = 0.1 kg and mc = 0.02 kg, using equation A-3. The peak force of the selected actuator is 22.2 N . RMS force is used to approximate the average continuous force requirement of an application and is described as [14] V h + t2 + t3 + u A.2 Sizing Voice Coil Actuators 111 where t± is the acceleration time, t2 is the run time, £3 is the deceleration time, and £4 is the dwell time is a move profile (Figure A.l). v . v max fc> l2 S t4 t Figure A . l : Voice coil move profile. The RMS force is 11.46 N, calculated by the definition. Stroke may be specified as the total displacement from one end of travel to the other end, or as a plus/minus (±) displacement from a mid-stroke reference point. A ±0.002 m stroke is usually required to vibrate the tissue. The stroke of the selected actuator is ±0.005 m. The peak force of the absorber actuator is different for different frequencies of excitation. Equa-tion 3.4 gives the relation between the control force and the frequency of excitation, depicted in Figure 3.4 (a) and (b). The maximum value of the control force is required at 5 Hz and 25 Hz excitation and is 2.9 N according to Figure 3.4 (b). This value is much less than the peak force that the selected actuator can apply, 22.2 N. The required stroke is also a function of the frequency of excitation, as stated by equation 3.2. The maximum required stroke is 0.005 m, according to this equation (also depicted in Figure 3.4 (c). The stroke of the selected actuator is ±0.005 m.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- A hand-held probe for vibro-elastography
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
A hand-held probe for vibro-elastography Rivaz, Hassan 2005
pdf
Page Metadata
Item Metadata
Title | A hand-held probe for vibro-elastography |
Creator |
Rivaz, Hassan |
Date Issued | 2005 |
Description | Vibro-elastography is a new medical imaging method that identifies the mechanical properties of tissue by measuring tissue motion in response to a multi-frequency external vibration source. Previous research on vibro-elastography used ultrasound to measure the tissue motion and system identification techniques to identify the tissue properties. This thesis describes a hand-held probe with a combined vibration source and ultrasound transducer which is operational in the 5-25 Hz range to cover a significant bandwidth of the tissue response. The design uses a vibration absorption system to counter-balance the reaction forces from contact with the tissue. Four different types of vibration absorption are briefly described and among them, active dynamic vibration absorption is selected for the hand-held device. A proportional integrator control is designed for the active vibration absorption and is compared with a well-established active vibration absorption method, the delayed resonator. An electromagnetic actuator is selected for active vibration absorption with a single accelerometer providing the feedback. The dynamics of the electromagnetic actuator as well as the response of the different niters are considered in the control law. The stability of the both absorber system and the combined system are considered to find the operational frequency range. The concept of tuning speed for different vibration absorbers is elaborated and is related to the vibration absorption speed. The design of the hand-held device, which includes the active vibration absorber, is presented next. The design utilizes another electromagnetic actuator which is used to vibrate the tissue and is operated with displacement feedback. The design of the vibration absorber is elaborated next. In this design the effect of different absorber parameters on the actuator stroke, control force, stability of the combined system and the tuning speed of the absorber is studied. After manufacturing the device, the absorber parameters are identified in order to optimize the vibration absorption performance. Simulation results are provided next that verify the theories presented on the stability and tuning speed. The results show that 100% vibration absorption can be achieved in steady state for both the proportional integrator and the delayed resonator controllers. They also show that nonlinearities in the absorber system decrease the amount of vibration absorption. Experimental results are presented next, showing approximately 85% vibration absorption in the operational range. The existing parameter identification methods that identify the mechanical properties of the tissue using the ultrasound are modified for the hand-held device. Simulation results are provided to validate the revised parameter identification methods. The first elastograms obtained with the hand-held device are finally presented. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2009-12-15 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
IsShownAt | 10.14288/1.0080706 |
URI | http://hdl.handle.net/2429/16693 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2005-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2005-0616.pdf [ 7.17MB ]
- Metadata
- JSON: 831-1.0080706.json
- JSON-LD: 831-1.0080706-ld.json
- RDF/XML (Pretty): 831-1.0080706-rdf.xml
- RDF/JSON: 831-1.0080706-rdf.json
- Turtle: 831-1.0080706-turtle.txt
- N-Triples: 831-1.0080706-rdf-ntriples.txt
- Original Record: 831-1.0080706-source.json
- Full Text
- 831-1.0080706-fulltext.txt
- Citation
- 831-1.0080706.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0080706/manifest