=Paper=
{{Paper
|id=Vol-2783/paper05
|storemode=property
|title=Medical Ultrasound Tomography Problem: Experimental Data Processing with High-Performance Computing
|pdfUrl=https://ceur-ws.org/Vol-2783/paper05.pdf
|volume=Vol-2783
|authors=Victoria Filatova,Leonid Pestov
|dblpUrl=https://dblp.org/rec/conf/birthday/FilatovaP20
}}
==Medical Ultrasound Tomography Problem: Experimental Data Processing with High-Performance Computing==
https://ceur-ws.org/Vol-2783/paper05.pdf
Medical Ultrasound Tomography Problem:
Experimental Data Processing with
High-Performance Computing ? ??
Victoria Filatova[0000−0003−3029−1550] and Leonid Pestov[0000−0001−6831−3669]
Immanuel Kant Baltic Federal University, Kaliningrad, Russia
ViFilatova@kantiana.ru, LPestov@kantiana.ru
Abstract. The article is devoted to imaging in Ultrasound Tomography.
We describe processing of experimental data obtained from an ultrasonic
tomographic prototype in order to reconstruct an image of an acoustical
medium. For visualization, we use the well-known geophysical method
Reverse time migration (RTM). This work is the continuation of the re-
search, earlier the method was tested on synthetic data, the stability of
the method was studied and numerical experiments were carried out for
2-dimensional and 3-dimensional models of acoustic medium. The article
presents the results of reconstruction based on the data measured for a
phantom object using an experimental sample of the ultrasound tomo-
graph developed at the Faculty of Physics, Lomonosov Moscow State
University. The numerical solution of the problem is a resource consum-
ing computational task. We use the standard parallel computing library
(OpenMP) and a cluster system to reduce the processing time.
Keywords: Ultrasound tomography · Visualization · Reverse time mi-
gration · Real data processing · Breast imaging.
1 Introduction
The main purpose of the medical ultrasound tomography (UST) is the breast
cancer diagnostics. Despite the existing known diagnostic methods (mammog-
raphy, computed tomography, MRI), scientists are developing a safe, cheap and
simple methods for detecting breast neoplasms. Currently, both the ultrasonic
tomographic devices and the methods of processing the data obtained from them
are under active research. It is worth noticing several scientific groups provid-
ing modern investigation in ultrasound tomography: from Germany under the
leadership of N. Ruiter [11], from the USA under the leadership of N. Duric
[5] and from Russia under the direction of O. D. Rumiantseva [3]. Existing two-
dimensional and three-dimensional ultrasound tomographs use 256 transducers
?
This work has been supported by the Russian Science Foundation (grant No. 16-11-
10027).
??
Copyright c 2020 for this paper by its authors. Use permitted under Creative Com-
mons License Attribution 4.0 International (CC BY 4.0).
58 V. Filatova, L. Pestov
or more. In the two-dimensional case the transducers are located on a circle of
radius from 0.1 m. The dominant impulse frequency is from 1 to 8 MHz.
From a mathematical point of view, the UST problem is an inverse problem
for the wave equation associated with determining the parameters of acoustic
medium (speed of sound, attenuation, density) using boundary measurements.
The object is probed by an acoustic impulse emitted from the transducer. The
inverse data is the measured acoustic field for different positions of the transduc-
ers. For medical diagnostics, it is important to reconstruct both the image of the
acoustical medium and the acoustical parameters values. The main methods for
solving inverse problem for the wave equation are the ray method, the lineariza-
tion method (Born approximation) and the optimization method. Typical UST
systems use ray-based tomographic algorithms to reconstruct this information
from the first-arrival times [10]. A two-step algorithm based on the linearized
inverse kinematic problem and the Born approximation was applied in [2]. An-
other algorithm, which is considered by O. D. Rumiantseva and her group, is
the Grinevich-Novikov algorithm [2], where scattering data of the plane waves
in a monochromatic mode are used. N. Duric’s group are developing their own
ultrasound tomograph, they uses the Born approximation [20] too. The group of
methods from geophysics adapted for ultrasound tomography is used. The ap-
proach based on the Kirchhoff migration method is applied in [19]. An iterative
method full-waveform inversion (FWI) has been applied to 2D and 3D ultra-
sound datasets [1, 17]. The research presented in [13, 14] use the optimization
method.
We continue our research conducted before. In the papers [6, 7], we proposed
a method based on visualizing the acoustical medium using the RTM procedure
and determining the values in the inclusions using the kinematic approach. Pre-
viously, numerical research (stability study, 2D and 3D experiments) was carried
out for simulated data [6–9]. In the current work, the first results of experimental
data processing are presented. We consider experimental dataset measured in a
real tomographic experiment with phantom object (Department of Acoustics of
Moscow State University, [16]). We thank the group O. D. Rumyantseva for the
opportunity to work with the real data. The paper [3] presents the results of the
same data processing.
2 Medical Ultrasound Imaging Problem
2.1 The Problem Statement
Consider the acoustic dynamical system with zero Cauchy data in R2 × (0, T ):
1
c2 pt = div v + δ(x − xs ) r(t),
vt = ∇p, (1)
p|t=0 = 0, v|t=0 = 0.
Here p(x, t; xs ) is the acoustic pressure; v(x, t; xs ) is the particle velocity vector
field; δ(x − xs ) is Dirac delta function, which models a point source located at
Medical ultrasound tomography problem 59
xs ∈ Γ = {x : |x| = R}; r(t) is an impulse (see below Fig. 2); c = c(x) is the
speed of sound. Outside the disc Ω = {x : |x| 6 R} the function c is equal to a
constant. The waves generated by the boundary sources are reflected/scattered
from some objects and recorded on the circle Γ .
Let T ∗ is the “acoustical” radius of disc Ω, i. e. minimal time that is required
to fill disc Ω with waves initiated by all boundary sources. The registration
time T must be greater than T ∗ .
The imaging problem consists in the reconstruction of sound-speed image
based on the pressure measurements on the boundary for different positions of
the sources (inverse data):
p0 (x, t; xs ) = p(x, t; xs ), x, xs ∈ Γ, t ∈ [0, T ].
Absorption reconstruction is important for diagnosis as well, but we do not
consider this problem here.
Remark 1. We consider forward problem (1) without boundary conditions. In
practice, an absorbing tire is installed in the ultrasound tomograph, which ab-
sorbs the waves. For modeling, we use perfectly matched layers to constrain the
computational domain and attenuate waves reflected from the boundary of the
inner square (details see in Section 4).
2.2 Visualization
The RTM-image is a result of the following procedures:
i. We calculate the pair pf (x, t; xs ), v f (x, t; xs ) for known approximation of the
sound-speed c0 (x) (forward problem (1)).
ii. For a given p0 (inverse data), we solve the reversal time problem:
1
2 pt = div v,
c0
(2)
vt = ∇p + p0 δ(|x| − R)νΓ ,
p|t=T = 0, v|t=T = 0,
where νΓ is the outward unit normal vector to Γ , t ∈ (0, T ). Let the pair
pb , v b is the solution of the problem (2).
iii. For the imaging process, we calculate the so-called “Energy Imaging condi-
tion” IE (x) [18]:
T
XZ
IE (x) = [pf pb + c20 (v f , v b )](x, t; xs ) dt, (3)
xs 0
where (·, ·) is the inner product in Rn .
60 V. Filatova, L. Pestov
We omit the mathematical base of the formula (3). Notice only, that reverse
time migration is related to linearized inversion [15]. To explain this formula
(naively) let us assume that we deal with only one diffractor. Using kinematic
argument it easy to see that waves in the direct and reversal time meet each
other right over this diffractor (in the space-time domain, see Fig. 1). Notice
also, that the conventional formula for RTM image does not contain integrand
c20 (v f , v b ).
Fig. 1. The waves in the direct and reversal time in the space-time domain, see for-
mula (3).
3 Experiment Description
The data were measured in a real tomographic experiment with a phantom
object in the form of a boiled egg with two wires threaded (perpendicular to the
tomography plane) through it. The phantom object is placed in water (cwater =
1.494 m/ms).
Here is a brief description of the tomograph prototype (see details in [3]).
There are 256 quasi-point transducers (each works as both a source and a re-
ceiver) with a central operating frequency of about 1.25 MHz and located on
a circle with a radius of 0.1536 m. When one transducer emits, all the others
receive. The emitting impulse r(t) is presented in Fig. 2. The frequency of dig-
itizing the signal is 5 MHz, i.e. it is strictly equal to four times the dominant
frequency (∆t = 0.0002 ms). Sampling of each signal starts from the 341st sam-
ple from the moment of signal emission (the moment of emission is considered
the first sample), in order not to accept an uninformative initial part. An ex-
ample of the raw signal received by transducers induced by a single source is
illustrated in Fig. 3. Figures 3 and 4 have a qualitative meaning, showing the
direct and reflected waves with noise.
Medical ultrasound tomography problem 61
Fig. 2. Impulse r(t).
Fig. 3. Raw data (pressure, single source).
62 V. Filatova, L. Pestov
4 Numerical Experiment
The spatial domain is a square 0.32 m ×0.32 m with ∆x = 0.0001 m grid step.
The number of nodes is 10240000. The transducers are uniformly located on a
circle of radius R = 0.1536 m centered at the origin. The number of transducers
is 256. We processed the raw signal using muting the direct waves, band-pass
filtering and interpolating in time from 0.0002 ms to 0.00001 ms (see Fig. 4).
Time step (∆t) is equal to 0.00001 ms.
Fig. 4. Processed data (pressure, single source).
The numerical solutions of both the forward problem (1) and the reversal
time problem (2) were based on the explicit conditionally convergent Finite-
Difference Time-Domain method (FDTD) with a staggered grid [12]. The 12th-
order accuracy approximation with respect to space was used and the second-
order to time. Spatial and time steps were determined according to the Courant
condition,
∆x
∆t < √ ,
kcmax 2
where the coefficient k depends on the order of the approximation (for the 12th
order, k ≈ 1.34).
The function c is a constant outside the disc Ω and, therefore, no waves reflect
from R2 \ Ω. Thus, we restrict the computational domain by using the well-
known technique of perfectly matched layers (PML; see [4]). The computational
domain (a square) contains an absorption layer with parameters specially chosen
to attenuate waves reflected from the boundary Γ0 (see Fig. 5). We use the
Dirichlet condition p|Γ0 =0.
Medical ultrasound tomography problem 63
Fig. 5. Scheme of the computational domain.
The background speed of sound is known (c0 = cwater ). We applied the
RTM procedure and obtained the image (Fig. 6). The RTM image has a lot of
oscillation and this is due to the emitting impulse is long-periodic. The image
does not contain the wire puncture site and the egg.
Fig. 6. RTM image.
64 V. Filatova, L. Pestov
Then we applied empirical formula for the imaging process:
τ (x,x
Z s )+δ
X
I(x) = [(pb )2 + c20 |v b |2 ](x, t) dt, (4)
xs
τ (x,xs )−δ
where τ (x, xs ) is the travel time between the points x and xs . This formula is
explained by the following reasoning. If the impulse brought by waves in the
direct and reversal time is long then the integral (3) is nonzero in big enough
neighborhood of the diffractor. The shorter impulse, the smaller neighborhood
is. But we deal with long-periodic impulse. So we need to use cut-off function in
(3) which is non-zero only near the travel time τ (x, xs ).
We calculated (4) for δ = 50∆t (see Fig. 7a), δ = 500∆t (see Fig. 7b), and
δ = 700∆t (see Fig. 7c). Wire punctures are clearly visible in the image, with the
left puncture displayed clearly and the right puncture weaker. This is explained
by the fact that one puncture was filled with water, and in the other there was
a wire. As δ increases, an egg is detected in the image (see Fig. 7b, c).
5 Conclusions
The authors conducted an experimental data processing. A numerical experi-
ment was conducted on the cluster system using OpenMP technology. The RTM
image contains a lot of noise and it is impossible to identify the egg and wire.
The oscillations on the image are explained by the principle of RTM imaging
and the long-periodic impulse emitted by the source. We proposed an empirical
formula for the imaging process. As a result, the egg and punctures of the wire
were detected on the image. In [3], first large-scale inhomogeneities with a low
spatial resolution of 0.5 – 1 cm (egg) are reconstructed, then the fine structure
of the scatterers with a spatial resolution of 0.25 mm (wire). In the experiment
described in this work, both the egg and the wire are reconstructed with a spa-
tial resolution corresponding to a grid step of 0.1 mm. The sizes and locations
of the scatterers (eggs and wire) are in accordance with values indicated in the
work [3].
Acknowledgements
The work was supported by the Russian Science Foundation (grant No. 16-11-
10027). We are deeply grateful to O. D. Rumyantseva and D. I. Zotov (Depart-
ment of Acoustics of Moscow State University) for sharing the 2D data.
Medical ultrasound tomography problem 65
Fig. 7. Resulting images: δ = 50∆t (a), δ = 500∆t (b), δ = 700∆t (c).
66 V. Filatova, L. Pestov
References
1. Agudo, O.C., Guasch, L., Huthwaite, P., Warner, M.: 3D imaging of the breast
using full-waveform inversion. In: SPIE Proceedings, pp. 99–110. KIT Scientific
Publishing (2017)
2. Burov, V.A., Grishina, I.M., Lapshenkina, O.I., Morozov, S.A., Rumyant-
seva, O.D., Suhov, E.G: Recovery of the fine structure of the acoustic scatterer
in the background of the distorting influence of its large-scale components. Acous-
tical Physics 49(6), 738–750 (2003)
3. Burov, V.A., Zotov, D.I., Rumyantseva, O.D.: Reconstruction of the sound ve-
locity and absorption spatial distributions in soft biological tissue phantoms from
experimental ultrasound tomography data. Acoustical Physics 61, 231–248 (2015)
4. Colino, F., Tsogka, C.: Application of PML absorbing layer model to the linear
elastodynamic problem in anisotropic heterogeneous media. Geophysics 66(1),
294–307 (2001)
5. Duric, N., Littrup, P.: Breast ultrasound tomography. In: Kuzmiak, C. (ed.) Breast
Imaging, vol. 6. IntechOpen (2017). https://www.intechopen.com/books/breast-
imaging/breast-ultrasound-tomography
6. Filatova, V.M., Nosikova, V.V.: Determining sound speed in weak inclusions in the
ultrasound tomography problem. Eurasian Journal of Mathematical and Computer
Applications 6(1), 11–20 (2018)
7. Filatova, V.M., Nosikova, V.V., Pestov, L.N.: Application of Reverse Time Migra-
tion (RTM) procedure in ultrasound tomography, numerical modeling. Eurasian
Journal of Mathematical and Computer Applications 4(4), 5–13 (2016)
8. Filatova, V.M., Pestov, L.N., Nosikova, V.V., Rudnitskii, A.G.: Breast ultrasound
tomography problem, simulation with noisy model. In: 2018 Days on Diffraction
(DD), pp. 106–111. St. Petersburg, Russia (2018)
9. Filatova, V.M., Pestov, L.N., Sedaikina, V.A., Danilin, A.N.: Breast ultrasound
tomography problem: 3D simulation. In: 2019 Days on Diffraction (DD), pp. 42–
45. St. Petersburg, Russia (2019)
10. Gemmeke, H., Hopp, T., Zapf, M., Kaiser, C., Ruiter, N.V.: 3D ultrasound com-
puter tomography: Hardware setup, reconstruction methods and first clinical re-
sults. Nuclear Instruments and Methods in Physics Research A: Accelerators, Spec-
trometers, Detectors and Associated Equipment 873, 59–65 (2017)
11. Gemmeke, H., Ruiter, N.V.: 3D ultrasound computer tomography for medical
imaging. Nuclear Instruments and Methods in Physics Research A 580(2),
1057–1065 (2007)
12. Ke, B., Zhao, B., Liu, C., Fang, Y.: Wave-equation datuming based on a single shot
gather. In: SEG Technical Program Expanded Abstracts, pp. 2156–2159 (2004)
13. Natterer, F., Wubbeling, F.: A propagation-backpropagation method for ultra-
sound tomography. Inverse Problems 11(6), 1225–1232 (1995)
14. Natterer, F.: Reflectors in wave equation imaging. Wave Motion 45(6), 776–784
(2008)
15. Op’t Root, T., Stolk, C., de Hoop, M.: Linearized inverse scattering based on
seismic reverse time migration. Journal de Mathématiques Pures et Appliquées
98(2), 211–238 (2012)
16. Parkhomenko, P.P., Karavaj, M.F., Sukhov, E.G., Faleev, B.A., Dmitriev, O.V.,
Drozdov, S.A., Komarov, O.V., Babin, L.V., Popov, A.S., Burov, V.A., Rat-
tehl, M.I., Bobov, K.N., Konjushkin, A.L., Rumjantseva, O.D.: Ultrasound tomo-
graph and ring antenna array for ultrasound tomograph. Patent RU No. 2145797
C1 . Date of registration: 23.06.1999. Moscow
Medical ultrasound tomography problem 67
17. Pratt, R.G., Huang, L., Duric, N., Littrup, P.: Sound-speed and attenuation imag-
ing of breast tissue using waveform tomography of transmission and ultrasound
data. In: SPIE Proceedings, pp. 1–12. KIT Scientific Publishing (2007)
18. Rocha, D., Tanushev, N., Sava, P.: Acoustic wavefield imaging using the energy
norm. In: 2015 SEG Annual Meeting, pp. 49–68. New Orleans, Louisiana (2015)
19. Schmidt, S., Duric, N., Li, C., Roy, O., Huang, Z.F.: Modification of Kirchhoff
migration with variable sound speed and attenuation for acoustic imaging of media
and application to tomographic imaging of the breast. Medical Physics 38(2),
998–1007 (2011)
20. Simonetti, F., Huang, L., Duric, N., Rama, O.: Imaging beyond the Born ap-
proximation: An experimental investigation with an ultrasonic ring array. Physical
Review E 76(3) 036601 (2007)