A variational formulation to image inclusions in 2D travel time cross hole tomography

Abstracts

It was developed a simple and high speed method with the purpose to model inclusions of irregular forms (simply connected and convex domain) in the context of cross hole travel time 2D tomography. The extension to 3D can be easily implemented and is not covered in this paper. The following assumptions were made: a) Homogeneus host medium (therefore straight ray propagation assumed) b) Homogeneus inclusion. c) Not very high velocity contrast between host medium (V1) and the inclusion (V2). The algorithm calculate the host velocity V1 considering the (straight) ray between source 1 and receiver 1, and considering that this ray not intersect the anomaly. Then the velocity V2 of the inclusion is estimated upon determining the difference deltat among measured times source i - receiver j and corresponding calculated times. If the global difference is negative signifies that V2 will be greater than V1 and viceversa. The small values of deltat are treated as outliers and therefore eliminated. The last initial parameters needed by the algorithm are the coordinates of the center of the inclusion, which are estimated by a graphical procedure. Afterwards the program approximates the contour of the inclusion using the Nelder-Mead 'simplex' algorithm. At last the perfomance of the algorithm was applied with real data.

Tomography; Cross Hole; Simplex Downhill; Variational; Nelder Mead

Se desarrolló un método simple y veloz para modelar inclusiones de formas irregulares (dominio convexo, simple conexo) en el contexto del método de tomografía sísmica en cross hole 2D. La implementación a 3D es fácil de llevar a cabo y no es tratada en este trabajo. Las siguientes premisas se tuvieron en cuenta: a) Homogeneidad del sustrato (por consiguiente se consideró la teoría de propagación de rayo recto). b) Inclusión homogénea. c) Un bajo contraste de velocidades entre el sustrato (V1) y la inclusión (V2). El algoritmo calcula la velocidad del sustrato V1 considerando el rayo (recto) entre la fuente '1' y el receptor '1' (superiores) y considerando que este rayo en particular no intersecta la anomalía. Luego la velocidad V2 de la inclusión se calcula a partir de la diferencia ddeltat entre los tiempos medidos fuente 'i' receptor 'j' y los correspondientes tiempos calculados con las premisas citadas. Si la diferencia global es negativa significa que V2 es mayor que V1 y viceversa. Los valores muy pequeños de ddeltat son tratados como puntos fronterizos y no tenidos en cuenta. Los últimos parámetros iniciales necesarios por el algoritmo son las coordenadas del centro de la inclusión, las cuales se estiman a partir de un procedimiento gráfico. Posteriormente el programa aproxima el contorno de la inclusión mediante el empleo del algoritmo 'simplex' de Nelder Mead. Por último el funcionamiento del programa se puso a prueba con datos reales.

Tomografía; Cross Hole; Simplex Downhill; Variacional; Nelder Mead

A variational formulation to image inclusions in 2D travel time cross hole tomography

Armando Luis Imhof^I; Carlos Adolfo Calvo^II

^IInstituto Geofísico Sismológico Volponi / Facultad de Ciencias Exactas Físicas y Naturales. Ignacio de la Roza y Meglioli. Rivadavia. C.P. 5400, San Juan, República Argentina. Fax: [phone]. Tel.: [phone]. E-mail: [e-mail]

^IIDepartamento de Matemática/ Facultad de Ingeniería. Av Libertador Gral San Martín 1109 (O) - CP 5400, San Juan, República Argentina. E-mail: [e-mail]

ABSTRACT

It was developed a simple and high speed method with the purpose to model inclusions of irregular forms (simply connected and convex domain) in the context of cross hole travel time 2D tomography. The extension to 3D can be easily implemented and is not covered in this paper.

The following assumptions were made:

a) Homogeneus host medium (therefore straight ray propagation assumed)

b) Homogeneus inclusion.

c) Not very high velocity contrast between host medium (V₁) and the inclusion (V₂).

The algorithm calculate the host velocity V₁ considering the (straight) ray between source 1 and receiver 1, and considering that this ray not intersect the anomaly. Then the velocity V₂ of the inclusion is estimated upon determining the difference dt among measured times source_i - receiver_j and corresponding calculated times. If the global difference is negative signifies that V₂ will be greater than V₁ and viceversa. The small values of dt are treated as outliers and therefore eliminated.

The last initial parameters needed by the algorithm are the coordinates of the center of the inclusion, which are estimated by a graphical procedure.

Afterwards the program approximates the contour of the inclusion using the Nelder-Mead 'simplex' algorithm.

At last the perfomance of the algorithm was applied with real data.

Keywords: Tomography, Cross Hole, Simplex Downhill, Variational, Nelder Mead.

RESÚMEN

Se desarrolló un método simple y veloz para modelar inclusiones de formas irregulares (dominio convexo, simple conexo) en el contexto del método de tomografía sísmica en cross hole 2D. La implementación a 3D es fácil de llevar a cabo y no es tratada en este trabajo.

Las siguientes premisas se tuvieron en cuenta:

a) Homogeneidad del sustrato (por consiguiente se consideró la teoría de propagación de rayo recto).

b) Inclusión homogénea.

c) Un bajo contraste de velocidades entre el sustrato (V₁) y la inclusión (V₂).

El algoritmo calcula la velocidad del sustrato V₁ considerando el rayo (recto) entre la fuente '1' y el receptor '1' (superiores) y considerando que este rayo en particular no intersecta la anomalía. Luego la velocidad V₂ de la inclusión se calcula a partir de la diferencia ddt entre los tiempos medidos fuente 'i' receptor 'j' y los correspondientes tiempos calculados con las premisas citadas. Si la diferencia global es negativa significa que V₂ es mayor que V₁ y viceversa. Los valores muy pequeños de ddt son tratados como puntos fronterizos y no tenidos en cuenta.

Los últimos parámetros iniciales necesarios por el algoritmo son las coordenadas del centro de la inclusión, las cuales se estiman a partir de un procedimiento gráfico.

Posteriormente el programa aproxima el contorno de la inclusión mediante el empleo del algoritmo 'simplex' de Nelder Mead.

Por último el funcionamiento del programa se puso a prueba con datos reales.

Palabras Clave: Tomografía, Cross Hole, Simplex Downhill, Variacional, Nelder Mead.

INTRODUCTION

It was developed a simple and high speed algorithm with the purpose to delineate inclusions in the context of cross hole 2D tomography.

The assumptions made were:

a) Homogeneus host medium (therefore straight ray propagation assumed)

b) Homogeneus inclusion.

c) Not very high velocity contrast between host medium (V₁) and the inclusion (V₂).

The third assumption due to the fact that the algorithm consider that the acoustic impedance (z=r.V) among both media do not differ very much.

In a real case the error due to these assumptions is expected not to be large due to the fact that it was proven that a small perturbation in travel time (say 1%) represents an alteration in travel distance (due to refraction mainly) up to 7% (SANTAMARINA; CESARE, 1994).

Figure 1 gives a resemblance of the cross-hole type array.

BASIC THEORY

At each reception point (r_j) it was calculated the corresponding travel time from source 's_i' to receiver 'r_j', conforming the travel time matrix T:

T = {t_ij}

Denoting D_ij the distance from s_i to r_j:

D = {D_ij}

D_ij = ||s_i - r_j||

Then the t_ij are:

t_ij =

If there exists an anomaly (medium with velocity V₂) insert in medium 1, the matrix of travel times T is altered and it is formed a new matrix T':

T' = {t'_ij }

Inside the inclusion the ij-ray travels a d_ij distance with velocity V₂. Figure 2 represents this situation:

The computation of the travel time is:

t'_ij = t_ij^<1> + Dt_ij^<2>

where t'_ij is the ij-th perturbed travel time due to the inclusion presence; t_ij^<1> is the travel time in host medium and Dt_ij^<2> into the anomaly. Writing the equation in terms of velocity and distances:

In matrix notation:

Note that in real situations T' corresponds to the matrix of measured travel times.

Extracting the matrix 'd' of travel distances into the anomaly gives:

It is important to note that the matrix 'd' has non zero-elements corresponding to the rays that intersect the inclusion. Later when proccesing these elements d_ij ; there must be annulated those with values near zero (treated as outliers).

It was assumed that the inclusion pertains to a simply connected and convex domain whose contour 'C' is a soft curve provided that every i-j ray intercepts C in two points (except for the tangency points). Any convex domain fulfill this condition.

Let be these intersection points and (Figure 2)

Then the distances d_ij will be:

d_ij = || - ||

VARIATIONAL FORMULATION

Let y=y(x) be the contour curve equation . Denoting as 'L' to the closed curve length, it appears the functional L = L(y) conditioned by matrix 'd'.

If from the contour C it is altered 'y' to 'y+dy' (Figure 3) it'll be seen that 'L' turns to be greater; this is equivalent to say that the contour C corresponds to a minimum length solution (always conditioned by 'd')

It was accepted then, as a plausible hipothesis (variational one), that the solution contour correspond to that of minimum length:

Numerical Solution to Variational Calculus

The analitical solution to this variational problem is not possible to obtain in a general case; therefore it was faced by direct methods. These consist to replace the curve y(x) for a finite number of points.

For the inclusion it was replaced the continuous curve for a discretized one with the points and , resulting a polygonal with length L* of easy determination.

Due to the fact that the matrix of distances 'd' is known, is easily determined from (; ) and since ray path of ij-th ray is known; only remains as independent variable (see Figure 4)

Finally, it stays formulated the problem to get minimum the function L*:

L*= L*()

where is a vector whose 'p' components correspond to the rays that intersects the inclusion.

Computational Implementation

To get minimum the function L*, it was employed the SIMPLEX algorithm of Nelder-Mead (PRESS et al., 1992).

An n-simplex is a convex polyhedron of 'n' dimensions determined by 'n+1' vertex (P₀,P₁,...P_n) linearly independent, whose x points are defined by:

with

0 < ti < 1 and = 1

SIMPLEX is an algorithm that do not use derivatives; it is slow at the beginning of the iteration process, but nevertheless converges more rapidly in the final stage than other methods (which employ derivatives).

The source code was implemented using MatLab.

Experimental Study

The efficiency of the computational algorithm was demonstrated on real data, which were obtained in the laboratory using sound waves in air. Figure 5 represents schematically the complete installation and the solution parameters.

The objetive was to determine five unknowns with the proposed algorithm:

V_inc, R_inc; X_inc; Y_inc and V_back meaning 'inc' for inclusion and 'back' for background.

As it is represented in Figure 5, the case under study involved a 2D circular inclusion (helium balloon with radius = 23cm) insert in a square image region. Tomographic crosshole data were obtained by radiating the host medium (air) from seven discrete source positions (SANTAMARINA; REED, 1994; SANTAMARINA; FRATTA, 1998) and acquiring the corresponding records at the receivers on the opposite side.

First of all it is necessary to determine which rays intersect the anomaly (from the delay times). Those delay times below a pre-established margin were discarded with the purpose to discard noise.

Figure 6 represents the rays from the source 'i' to the receiver 'j' that intersects the inclusion. Then it was necessary to evaluate the center coordinates of it by inferential means, clicking with the mouse at the inferred position.

It is easy to appreciate that the inclusion could be at x = 0.5m. from the origin.

This is the beginning datum. From it the algorithm calculate the and Simplex were run calculating the C contour of minimum length. Afterwards the points were interpolated and represented in three forms:

a) Representing the raypaths inside the anomaly (Figure 7)

b) Representing the anomaly contour (Figure 8)

c) Interpolating the points with an ellipse (Figure 9)

Solution Parameters (Calculated)

Figure 10 represents a detail of the X_inc,Y_inc and approximation of R_inc.

The V_back was calculated by the algorithm considering that the raypath among source 1 and receiver 1 (upper) does not intersect the anomaly.

R_inc= 0.26m.

X_inc= 0.5m.

Y_inc= 0.772m.

V_inc.= 600m/s

V_back.= 349m/s

CONCLUSIONS

It was proposed a new and simple method for delimitating homogeneous inclusions into a homogeneous background which is accurate and can model both high and low velocity anomalies.

The method was succesfully applied in more general domains (i.e. simply connected) without the convexity constraint.

The exact circular form is not obtained from the inversion scheme. This is caused by the lack of sources or receivers (or both) in the 'x' direction, which diminishes the resolution (SANTAMARINA; REED, 1994; SANTAMARINA; FRATTA, 1998). Inherently, the cross hole configuration is not the optimum to model anomalies accurately. It is necessary to place sources and/or receivers in the 'x' direction too.

Because of the straight ray theory assumption was considered, the errors upon the size of the inclusion will be greater when the contrast between host and inclusion velocities increases.

The initial parameters Vback, Xinc and Yinc, needed for an adequate start with simplex method were easily calculated with a simple calculus (for Vback) and inferential method (Xinc, Yinc).

Acknowledgements

This study is part of a research work on wave material interactions and applications. Tomographic data were downloaded with permission from the associated website (SANTAMARINA, FRATA, 1998).

The authors want to thank to professor Juan Carlos Santamarina from Georgia Tech (USA) for permission to use real case data from his laboratory in this paper.

Recebido em 16 dez. 2003 / Aceito em 28 maio 2004

Received dec. 16, 2003 / Accepted may 28, 2004

NOTES ABOUT THE AUTHORS

Armando Luis Imhof is geophysicist graduated in Universidad Nacional de San Juan (UNSJ, Argentina) in 1989. Specialist on Environmental Studies from Universidad Nacional de San Luis (Argentina); graduated 2002. Doctoral student at the Engineering Faculty (Universidad Nacional de Cuyo, Mendoza, Argentina). Researcher at the Instituto Geofísico Sismológico Volponi (Facultad de Cs. Exactas, Fcas y Naturales - UNSJ) and adjunct professor at the Electrical Prospecting, Geophysics and Astronomy Department (UNSJ) since 1990. Member of Asociación Argentina de Geofísicos y Geodestas (A.A.G.G.). Topics of main Research: Geophysics applied to mining, geotechnics an hidrogeological studies.

Carlos Adolfo Calvo is electromechanical engineer graduated in Universidad Nacional de San Juan (UNSJ, Argentina) in 1978. MSc in Applied Mathematics and Informatics at Moscow Physics Institute (MEPHI). Main Professor of Mathematics Analysis, Applied Mathematics and Numerical Methods; Engineering Faculty (U.N.S.J). Researcher on topics covering applied mathematics to solve engineering problems. Topics of main research: Applied mathematics.

• POTTS, B. D.; SANTAMARINA, J. C. Geotechnical tomography: the effects of diffraction. Geotechnical Testing Journal, [S.l.], v. 16, n. 4, p. 510-517, 1993.
• PRESS, W. H. et al. Numerical recipes in Fortran: the art of scientific computing. New York. Cambridge University Press, 1992.
• SANTAMARINA, J. C. 2000. Available in: <http://www.ce.gatech.edu/~carlos/laboratory/publication/sigma/ch11/ch11.htm
• ______; CESARE, M. A. Velocity inversion in the near surface: vertical heterogeneity and anisotropy. Internal report University of Wat+erloo, Canada, 1994.
• ______; REED, A. C. Ray tomography: errors and error functions. Journal of Applied Geophysics, [S.l.], n. 32, p. 347-355, 1994.
• ______; FRATTA, D. Introduction to discrete signals and inverse problems in civil engineering. Reston: ASCE Press, 1998.

Publication Dates

• Publication in this collection
  24 Jan 2007
• Date of issue
  Dec 2003

History

• Accepted
  28 May 2004
• Received
  16 Dec 2003