We are IntechOpen,

the world’s leading publisher of

Open Access books
Built by scientists, for scientists

4,200
Open access books available
116,000
International authors and editors
125M
Downloads

Our authors are among the

154
Countries delivered to
TOP 1%
most cited scientists
12.2%
Contributors from top 500 universities

Selection of our books indexed in the Book Citation Index

in Web of Science™ Core Collection (BKCI)

Interested in publishing with us?

Contact [email protected]
Numbers displayed above are based on latest data collected.
For more information visit www.intechopen.com
Chapter

Electrical Resistivity Tomography:

A Subsurface-Imaging Technique
Bing Zhou

Abstract

Electrical resistivity tomography (ERT) is a popular geophysical subsurface-

imaging technique and widely applied to mineral prospecting, hydrological explo-
ration, environmental investigation and civil engineering, as well as archaeological
mapping. This chapter offers an overall review of technical aspects of ERT, which
includes the fundamental theory of direct-current (DC) resistivity exploration,
electrode arrays for data acquisition, numerical modelling methods and tomo-
graphic inversion algorithms. The section of fundamental theory shows basic for-
mulae and principle of DC resistivity exploration. The section of electrode arrays
summarises the previous study on all traditional-electrode arrays and recommends
4 electrode arrays for data acquisition of surface ERT and 3 electrode arrays for
cross-hole ERT. The section of numerical modelling demonstrates an advanced
version of finite-element method, called Gaussian quadrature grid approach, which
is advantageous to a numerical simulation of ERT for complex geological models.
The section of tomographic inversion presents the generalised standard conjugate
gradient algorithms for both the l1- and l2-normed inversions. After that, some
synthetic and real imaging examples are given to show the near-surface imaging
capabilities of ERT.

Keywords: resistivity, electrical current, geotomography, numerical modelling,

subsurface imaging

1. Introduction

Direct-current (DC) resistivity exploration is a traditional geophysical method.

It employs two electrodes to inject electric current into the ground and other two
electrodes to measure the electric potential difference. The measurements are often
carried out along a line or in an area on the earth surface, and then the observed
potential differences are converted into sounding curves or pseudo-sections of
apparent resistivities, which indicate the resistivity changes of subsurface rocks.
Analyses of these data enable us to find the underground resistivity anomalies or
outline the subsurface geological structure. With development of computer tech-
nology and numerical computational techniques, accurate numerical simulations of
subsurface electrical field and acquiring a large amount of data in fields become
possible [1–3], so that the traditional DC resistivity exploration was developed to a
computerised geotomography technique, called electrical resistivity tomography
(ERT), which employs a multielectrode equipment or system to automatically
acquire a large number of data [4, 5] and applies a computer software to the

1
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

reconstruction of subsurface resistivity structure with the observed data [6–10]. Due
to its conceptual simplicity, low equipment cost and ease of use, ERT is now widely
applied in mineral exploration, civil engineering, hydrological prospecting and envi-
ronmental investigations, as well as archaeological mapping [11]. This chapter pro-
vides an overall review of ERT techniques, which consists of four sections: (1)
fundamental theory, (2) electrode arrays, (3) numerical modelling and (4) tomo-
graphic inversion. In each section, diagrams and formulations are used to illustrate
basic concepts and principles of ERT techniques. Some synthetic experiments and
practical imaging applications are also given to show the imaging capability of ERT.

2. Fundamental theory

According to the continuity of electrical current, the following integral equation

is satisfied at any point in a conductive medium:

∯ Γ J
dΓ ¼ I, (1)

where Γ is a full or half spherical surface that encloses an electrode that injects
electric current J with a magnitude of I (Figure 1). According to Ohm’s law:

J ¼ σE ¼ σ∇U, (2)

and the property of the delta-function δ(x-xs), Eq. (1) gives

ððð
∯ Γ σ∇U
dΓ ¼ I δðx xs ÞdΩ (3)
Ω

Here σ is the conductivity of medium, U is the electric potential, Ω represents

the volume of the spherical surface Γ, and xs = (xs,ys,zs) is the location of the current
electrode. Applying the divergence theorem to Eq. (3) leads to
ððð
½∇
ðσ∇U Þ þ Iδðx xs ފ ¼ 0: (4)
Ω

Figure 1.
Electric current density J around an electrode in (a) a full space and (b) a half-space.

2
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

Note that Eq. (4) satisfies everywhere in a medium, so that the following
governing equation of electric field is obtained:

∇:ðσ∇U Þ ¼ Iδðx xs Þ, :x ∈ Ω: (5)

In general, the conductivity σ(x) = σ(x,y,z) changes with three coordinates so

that it defines a 3D geological model, and the electric potential U=U(x,y,z) becomes
a 3D function of a given current injection Iδðx xs Þ and conductivity σ(x)
(Figure 2a). If the conductivity is a constant with the y-coordinate, i.e. σ(x) = σ(x,z),
it defines a 2D geological model (Figure 2b), but the electric potential U(x,y,z) is still
a 3D function of the coordinates due to the point current injection Iδðx xs Þ. Apply-
ing cosine Fourier transform Fc{
} to the y-coordinate of Eq. (5), the governing
equation is changed into
 
e þ k2 U
∇: σ∇U e ¼ Iδðx xs Þ=2, x ∈ Ω, (6)
y

which is often named 2.5D governing equation. Here ky is the wavenumber, and
 
e
U x; ky ; z ¼ Fc fUðx; y; zÞg becomes the spectrum of electric potential in the
wavenumber domain. The coordinate vector x = (x,z) and gradient ∇ = (∂x, ∂z) are
the 2D versions. If a geological model σ(x,y,z) or σ(x,z) and current injection
Iδðx xs Þ are given, solving the governing equation (5) or (6), one obtains electrical
 
potential U(x,y,z) for a 3D conductivity model or its spectrum U e x; ky ; z for a 2D
conductivity
n  model o and then performs the inverse cosine Fourier transform U(x, y, z)

=F 1 U e x; ky ; z to obtain electric potential U(x,y,z). These computations are
c
called forward modelling. In theoretical computation or numerical forward modelling,
e
Green’s function G or G—the electric potential response to a unit current injection
(I = 1 Am)—is often applied, so as to remove the magnitude of electric current I in
Eqs. (5) and (6). One may focus on the computation of Green’s function, and then the
e is computed by the multiplication of Green’s function G with
electric potential U or U
the magnitude of practical electric current I, i.e. U = I
G.
The simplest geological model is a homogenous half-space. Applying Eq. (3)
(equivalent to Eq. (5)) to a constant medium, the surface integral is calculated by
 
∯ Γ σ∇U
dΓ ¼ σ∇U 2πr2 ¼ I, (7)

which gives the electric potential at distance r from a current source:

Figure 2.
Sketch of (a) a 2D and (b) a 3D geological model defined by conductivity σ(x,z) and σ(x,y,z), respectively.

3
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

I ρI
U¼ ¼ , (8)
2πσr 2πr

where ρ = 1/σ is resistivity. In practice, to inject electric current into the ground,
a pair of current electrodes C1 and C2 must be employed; one is positive (+I) and
another is negative (-I). Thus, the electric potential at a point P with a pair of
current electrodes is calculated by.
 
ρI 1 1
U¼ : (9)
2π rC1P rC2P

Here, rC1P and rC2P are distances of the observed point P to two current elec-
trodes C1 and C2. To measure the potential on the earth surface, a pair of potential
electrodes P1 and P2 is also required (Figure 1). According to Eq. (9), one has the
following expression for the electric potential difference between two potential
electrodes:
   
ρI 1 1 1 1
ΔU ¼ , (10)
2π rC1P1 rC2P1 rC1P2 rC2P2

from which apparent resistivity ρa is defined as follows:

ΔU
ρa ¼ K ¼ KΔG, (11)
I

where ΔG = ΔU/I and K is called geometry factor of electrode array given by

    1
1 1 1 1
K ¼ 2π , (12)
rC1P1 rC2P1 rC1P2 rC2P2

which depends on the positions of four electrodes. Different layouts of four elec-
trodes have variable geometry factors and are often called electrode arrays. In the
traditional electrode arrays, C2 and P2 may be set up very far away from C1 and P1, so
that C2 and P2 are treated as remote electrodes in theory, and distances rC1P2 , rC2P2 ,
rC2P1 , and rC2P2 are supposed to be infinite (∞) in Eq. (12). In these cases, one can find
that geometry factors are still applicable for these electrode arrays that involve one or
two remote electrodes, which are named pole-pole, pole-dipole and dipole-pole arrays.
If subsurface resistivity is homogenous (ρ0), Eq. (11) shows that no matter
which electrode array is used, apparent resistivity is constant (ρa = ρ0). Otherwise,
ρa indicates resistivity variation of the underground. For a certain range of apparent
resistivity ρa, Eq. (11) also reveals that the geometry factor is inversely proportional
to the potential difference ΔG = ΔU/I. It implies that a big value of geometry factor
K will cause a small reading of the electric potential difference ΔG in fields. Such a
small reading is easily obscured or contaminated by background noise. For a good
data quality, one should avoid very big values of geometry factors in data acquisi-
tion for ERT. Therefore, Eqs. (11) and (12) are the fundament formulae of the
traditional DC electrical resistivity exploration.

3. Electrode arrays

In order to obtain apparent resistivity in fields, many electrode arrays were

developed in the traditional DC resistivity exploration. In principle, ERT requires a

4
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

high data density and good coverage of the earth surface for high-resolution images
of subsurface targets. Dahlin and Zhou [11] carried out synthetic experiments of
ERT using 10 electric arrays and compared their imaging results for four geological
models: a buried channel, a narrow dike, dipping blocks and waste ponds. They
demonstrated that two three-electrode arrays (pole-dipole and dipole-pole) and
three four-electrode arrays (dipole–dipole, Schlumberger and gradient arrays) pro-
duce satisfactory images of the subsurface targets. However, due to the use of
remote electrodes, the three-electrode arrays are rarely applied for ERT in practice;
thus, the four-electrode arrays become popular. Particularly, gradient array [12] is
well suited for multichannel data acquisition and can significantly increase the
speed of data acquisition in the field, and at the same time, it gives higher data
density and lower sensitivity to noise than dipole-dipole array. Figure 3 shows  three
four-electrode arrays (upper row) and their pseudo-sections of data points xρa ; zρa
(middle row) and their geometry factors (bottom row) for a layout of total 81
electrodes on the earth surface. Figure 3 shows that dipole–dipole array has the
biggest range of geometry factor among these three-electrode arrays, and gradient
array performs the smallest range, so that it is well suitable for data acquisition with
a high data density and good data quality. Zhou and Dahlin [13] based on many sets
of field data obtained with different electrode arrays and found that the measured
potential errors depend on the reading of potential difference, ΔG = ΔU/I, and the
measured datum in field may be expressed by

 α β
∗
ΔG ¼ ΔG 1 þ R =100 , (13)
ΔG

where ΔG* denotes a noise-contaminated datum, R is a random number in ( 1,1)

and α and β are constants. Due to having a big range of geometry factor (see
Figure 3), the dipole-dipole array is more sensitive to noise than Schlumberger and
gradient arrays, so that one must pay much attention to the noise contamination when
conducting dipole-dipole measurements [13]. According to the relationship between
the geometry factor K and reading of potential difference ΔG, Eq. (13) may be applied
to predications of noise-contaminated data for different electrode arrays [11].
Reviewing Figure 3, one may find that gradient array only uses the potential
electrodes between C1 and C2 to measure the potential differences. It does not do
the measurements at the potential electrodes outside of C1 and C2. It means that
gradient array actually misses the outside data for every pair of current electrodes.
To obtain a better data coverage and complete subsurface information, the outside

Figure 3.
Dipole-dipole, Schlumberger and gradient arrays (upper row) for ERT data acquisition and examples of their
pseudo-sections of data points (xρa , zρa Þ (middle row) and geometry factors (bottom row).

5
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

potential electrodes of each pair of C1 and C2 should be employed in the gradient

measurements. These additional data apparently complement to common gradient
measurements and may improve ERT imaging. Accordingly, a new electrode array
is naturally formed and called full-range gradient array shown in Figure 4a. The
difference from the common gradient array is that the new electrode array uses not
only the potential electrodes between C1 and C2 but also the outside potential
electrodes of every pair of C1 and C2, if their geometry factors fall in a reasonable
range for the data acquisition. Figure 4b and c shows the pseudo-section points and
geometry factors for a layout of total 81 electrodes. Comparing with Figure 3, the
full-range gradient array does improve the data coverage, and its geometry factors
are controlled in a reasonable range.
If there are boreholes in a field, cross-hole ERT may be carried out to image the
geological structure between the boreholes. Zhou and Greenhalgh [14] investi-
gated all possible electrode arrays for cross-hole ERT data acquisition and found
that the electrode arrays of pole-pole (A-M), pole-bipole (A-MN), bipole-pole
(AB-M), and bipole-bipole (AM-BN) with their multi-spacing cross-hole profiling
and scanning surveys are useful for cross-hole ERT. Here the capital letters A and
B stand for two current electrodes, and M and N denote two potential electrodes.
These cross-hole electrode arrays are shown in Figure 5 with their sensitivity
functions in backgrounds [15]. They also found that the electrode arrays which
have either both current electrodes or both potential electrodes in the same bore-
hole, e.g. A-MN, AB-M and AB-MN, have a singularity problem in data acquisi-
tion (geometry factor goes to infinite so that apparent resistivity and pseudo-
reaction are not applicable), namely, zero readings of the potential or potential
difference in cross-hole measurements, so that the potential data are easily
obscured by background noise and their images are inferior to those from other

Figure 4.
Full-range gradient array (a) for ERT data acquisition and examples of its pseudo-section of data points
xρa ; zρa (b) and geometry factors (c) with a layout of total 81 electrodes.

6
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

Figure 5.
Electrode arrays for cross-hole ERT data acquisition. A and B stand for two current electrodes. M and N denote
two potential electrodes. The background contours are the sensitivity functions of the electrode array.

cross-hole electrode arrays. The data having the singularity problem may be
predicated by zero values of the inverse geometry factors, which should be
avoided in cross-hole data acquisition. Therefore, A-M, AM-B, and AM-BN with
multi-spaces are recommended for cross-hole ERT.

4. Numerical modelling

To compute theoretical electric potential U(x,y,z) in a complex 2D or 3D

geological model, one has to solve the governing equation (6) or (5) using a
numerical approach and a computer. For computational efficiency, an under-
ground artificial boundary Γ0 is required to truncate an infinite geological model,
so that the boundary condition on Γ0 must be introduced to the governing equa-
tions. Applying Eq. (7) to the artificial boundary Γ0 and zero normal component
(J
n = 0) of electric current density J on the earth surface, one may find that the
electric potential U holds

σ∇U
n þ νU ¼ 0, (14)

where n is the normal vector of the artificial boundary Γ0 and ν is calculated by

ν = (r
n)/r2. Eq. (14) is the 3D boundary condition on Γ0. Similarly, one may find
the boundary condition for a 2D geological model [16]:

e
n þ νU
σ∇U e ¼ 0: (15)

Adding the artificial boundary conditions Eqs. (14) and (15) to Eq. (5) and (6),
numerical modelling becomes to solve the following definite governing equation:

∇
ðσ∇GÞ ¼ δ ðx xs Þ, x ∈ Ω,
3D : (16)
σ∇G
n þ vG ¼ 0, x ∈ Γ0

7
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

or
(  
~ þ k2 G
∇
σ∇G ~ ¼ δðx xs Þ=2, x ∈ Ω,
y
2D : (17)
~
n þ vG
σ∇G ~ ¼ 0, x ∈ Γ0 :

The numerical approach to the definite governing equations is called numerical

modelling. The most popular numerical approaches are finite-difference and
finite-element methods. The former is simple and straightforward, due to directly
applying a finite-difference formula to the derivatives of the gradient ∇, so that the
definite governing equations at all points in the model domain Ω and on the artificial
boundary Γ0 are discretised and assembled into a linear equation system [16–20]. The
latter converts the definite governing equation into an integral equation by weighting
residual principle or variational principle [21] and then carries out numerical volume
integration to produce a linear equation system [22–26]. Therefore, no matter which
approach is applied, the numerical modelling becomes to solve a linear equation
system. For an instant, the following paragraphs briefly show an advanced 2D version
of the finite-element method, called Gaussian quadrature grid (GQG) approach,
because of its advantages over finite-difference method and other schemes of finite-
element approaches, e.g. high accuracy and easy implementation of discretisation of a
geological model having arbitrary free-surface topography. For a 3D numerical
modelling, one may follow the 2D procedures.
Weighting residual principle is to calculate the following integral of the
governing equation in Eq. (17):
ðð h   i
W ∇
σ∇G e þ k2 G e þ δðx xs Þ=2 dΩ, (18)
y
Ω

where W is an arbitrary weighting function. Applying the divergence theorem to

the above and then submitting the artificial boundary condition yields
ðð   ð
e 2 e e
σ∇W
∇G þ σky WG dΩ þ νW GdΓ ¼ Wðxs Þ=2: (19)
Ω Γ0

In order to calculate the integrals, the model domain Ω is divided into a set of the
no-overlap subdomains {Ωe, e = 1,2,…,Ne} that matches the free-surface topography
and subsurface interfaces. In each subdomain, Gaussian abscissae {(xα, zβ),
α,β = 1,2,…,NG} are employed (see examples shown in Figure 6) to discretise the
subdomain and then Lagrange interpolation:

G e ðeÞ ,
e ¼ ∑N G li ðxÞ, lj ðzÞG (20)
i, j¼1 ij

and Gaussian weights {wαβ} are applied to calculation of the submain integrals.
Consequently, Eq. (19) becomes
n h   i ðeÞ
0 0 2 e
∑e ∑α, β, i, j wαβ σαβ li ðxα Þδjβ ∂x Wαβ þ δiα lj zβ ∂z Wαβ þ ky Wαβ δiα δjβ G αβ
o (21)
þ∑α να Wα G e ðeÞ ¼ Wðxs Þ=2:
α

 
Here, G e ðeÞ are the discrete values of the electric potential spectra G
e x; ky ; z in
αβ
the subdomain Ωe. Choosing the Lagrange basis polynomials as the weighting func-
tions, i.e. {Wpq ¼ lp ðxÞlq ðzÞ, ∀p, qg, Eq. (21) is changed into

8
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

Figure 6.
2D and 3D Gaussian quadrature grids for numerical modelling: (a) a 2D model, 8  10 subdomains and
7  7 Gaussian abscissae in each subdomain and (b) a 3D model, 10  10  10 subdomains and 5  5  5
Gaussian abscissae in each subdomain.

n h     i ðeÞ
∑e ∑α, β, i, j wαβ σαβ l0i ðxα Þl0p ðxα Þδjβ δqβ þ l0j zβ l0q zβ δiα δpα þ k2y δpα δqβ δiα δjβ G
e
αβ
    ðeÞ o
þ∑γ νγ lp xγ lq zγ G e ¼ lp ðxs Þlq ðzs Þ=2,
γ

(22)

which can be rewritten in the matrix form for all points (p,q) and (i,j) in Ω:

e ¼ bs ,
MG (23)

where M is a square symmetric matrix assembled by the coefficients of G e ðeÞ in

αβ
e
Eq. (22). G and bs are two vectors: the former consists of the discrete values of G e at
all points of the Gaussian quadrature grid, and the latter is a zero vector except for
the component of 1/2 at the electric current source located at xs. Eq. (22) shows that
the matrices M and G e depend on the wavenumber ky and conductivity model σαβ.
The wavenumber can be predicated by the minimum and maximum spaces of
electrodes [17]. The linear equation system may be efficiently solved for multiple
current sources B = {b1,b2,…,bN} by the Cholesky decomposition [27]. After
e the theoretical electric potentials are obtained by inverse
obtaining the spectra G,  
1 e
cosine Fourier transform Uðx; y; zÞ ¼ I
F G . c
From Eq. (22), one can find that the model conductivity σ(x,z) is discretised by
σαβ at the Gaussian abscissae instead of finite elements. The subdomain integrals are
calculated by the weights wαβ to the discrete integrands. The dense abscissae and
variable σαβ may describe the details of a complex geological model and generate
high accurate solutions, so that sizes of the subdomains are not necessarily small but
shaped to match the free-surface and subsurface interfaces. Apparently, accuracy of
the numerical modelling depends on the number of the Gaussian abscissae in the
subdomains. The more abscissae are applied, the more accurate results are yielded
but cost more computer time. For efficient and accurate modelling, the minimum
electrode space of electrode array may be chosen for the size of the subdomain, to
which five Gaussian abscissae is applied to produce satisfactory solutions of
numerical modelling. As examples, Figure 7 shows the pseudo-sections of dipole–
dipole, Schlumberger and gradient arrays for a layout of total 141 electrodes over an
anticline model. These results show that similar pseudo-sections of apparent resis-
tivity are obtained by using the three electrode arrays. From the above analysis, one

9
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

Figure 7.
Numerical modelling for (a) an anticline model and the apparent resistivity pseudo-sections of (b) dipole-
dipole, (c) gradient and (d) Schlumberger arrays. The discrete anticline model is also given in the pseudo-
sections.

can see that the finite-element approach does not calculate the high-order deriva-
tives for the governing equation, but any finite-difference scheme does. Therefore,
finite-element approach is often called a ‘weak’ solution of the governing equation
against a ‘strong’ solution obtained by finite-difference method.

5. Tomographic inversion

Tomographic inversion is to reconstruct the geological model that offers syn-

thetic data matching with observed data. Due to incompleteness and noise contam-
ination of observed data, the model reconstruction is ill posed (multiple solutions).
Therefore, tomographic inversion is often defined as an optimisation of data fit-
tingness with regularisation of the model [28], e.g. a generalised objective function
is applied:
 lp lp
ΦðmÞ ¼ dob dsyn ðmÞW þ λk m k W m
, (24)
d

where dob is a vector of observed data, which are either apparent resistivities or
potential differences measured by different electrode arrays. dsyn(m) stands for
synthetic data and is calculated by the finite-difference or finite-element method
for a guessed geological model m, which consists of discrete conductivities or
resistivities. λ is a regularisation parameter that plays a trade-off role between the
data fittingness (the first term in Eq. (24)) and model smoothness (the second term
p l
in Eq. (24)). k
kW stands for the weighted lp-norm with a weighting matrix W, e.g.
kεkW ¼ ε Wε and kεk1W ¼ ∑i W i ∣εi ∣ are the weighted l2-norm (generalised least
2 T

10
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

Figure 8.
An l2-normed (left) and an l1-normed (right) objective function. Here I is an unit matrix.

square) and the weighted absolute norm. Figure 8 illustrates the difference between
a l2-normed and l1-normed objective functions. Wd and Wm are weighting matrices
to data and model, respectively. It is common to choose a combination of finite-
difference operators for Wm, e.g. W m ¼ λ0 Iþλ1 DTx Dx þλ2 DTy Dy þλ3 DTz Dz [29].
Here I is a unit matrix, and Dx, Dy and Dz are the finite-difference operators in the
x-, y- and z-directions. λk (k = 0,1,2,3) are constants and called extensional
regularisation parameters used for searching for the smoothest model in three
directions. Therefore, tomographic inversion becomes solving the following opti-
mization problem:

m∗ ¼ minfΦðmÞg: (25)

To do so, a global or a local search may be applied to Eq. (25) [30], but the global
search is extraordinarily computer time-consuming if m has a larger dimension
[31]. Therefore, the local search of the standard conjugate gradient method is
commonly applied for tomographic inversion. Figure 9 gives a flowchart of the
conjugate gradient algorithm and shows that the gradient ∇Φ(mi) and the Hessian
matrix H(mi) of the objective function are required. Applying the linearised
∂d
approximation to the synthetic data, dsyn(mi + 1) ≈ dsyn(mi) + ( ∂msyn ÞT (mi + 1-mi),
and then substituting it for Eq. (24), one can obtain the gradient and the Hessian
matrix [28]:

Figure 9.
Flowchart of a standard conjugate gradient algorithm for tomographic inversion. Here m0 is an initial model,
ε is a small value, maxITS is the maximum iterations, and H(mi) is the Hessian matrix of object function Φ(mi).

11
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

 T
∂d 
∇Φðmi Þ ¼ Wd dob dsyn ðmi Þ þ λWm ðmi m0 Þ, (26)
∂mi
 T
∂d ∂d
Hðmi Þ ¼ Wd þ λWm , (27)
∂mi ∂mi

for the l2-normed objective function.

Note that the gradient ∇Φ(mi) of the l1-normed objective function has multiple
values at zero misfit Φ(mi) = 0 (see Figure 8), but for any non-zero misfit
Φ(mi) ≥ ε0 (ε0 is a very small value), the computations of the gradient ∇Φ(mi) and
the Hessian matrix H(mi) become simple, and they are similar to Eqs. (26) and (27)
except for the weighting matrices Wd and Wm. Therefore, as Φ(mi) < ε0, the l1-
normed gradient ∇Φ(mi) and the Hessian matrix H(mi) may be replaced with
Eqs. (26) and (27), and as Φ(mi) ≥ ε0 the l1-normed gradient and the Hessian
matrix are directly applied to the standard conjugate gradient algorithm. Therefore,
the l1-normed inversion and the l2-normed inversion are implemented with an
almost same algorithm. However, their inversion results may be quite different in
the case that the observed data have outliers. Many synthetic and practical experi-
ments have shown that the l1-normed inversion is less sensitive to the outliers [13].
If the data have no outliers or high qualities, two inversions converge [11].
Figure 10 shows synthetic experiments for imaging the anticline structure
shown in Figure 7a with dipole–dipole, Schlumberger and gradient arrays. The l2-
normed inversions were conducted with the apparent resistivities shown in
Figure 7b–d, which were computed with the GQG software. The commercial soft-
ware RES2DINV was applied to the inversions. From these results, one can see that
the dipole-dipole and gradient arrays yield competitive images of the anticline
structure. Figure 11 gives real applications of surface and cross-hole ERT conducted
at a rural site in Australia, where there are many existing drill boreholes which can
be used for examinations of both surface and cross-hole ERT experiments. From
these boreholes and logging data, the subsurface rocks and main structures were
well known (shown in Figure 11b). The main objective of this work was to use the
borehole geological and logging information and exam the imaging capabilities of
surface and cross-hole ERT in this area, particularly for mapping the base of alluvial
overburden and the base of pisolite, as well as predicting clay contamination within
pisolite. Surface data acquisition with dipole-dipole and Schlumberger arrays was

Figure 10.
Tomographic inversion of (a) dipole-dipole, (b) Schlumberger and (c) gradient arrays for the anticline model
shown in Figure 7.

12
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

Figure 11.
Applications of surface and cross-hole electrical resistivity tomography for mapping base of alluvial overburden:
(a) integrated resistivity images from surface and cross-hole ERT and (b) geological section from borehole rock
samples and logging data.

conducted along a 750 m line, and two pairs of cross-hole ERT on the same line were
also carried out for details of the geological structure between the boreholes.
Figure 11a shows the integrated cross-hole ERT results with the surface ERT imag-
ing. It gives the resistivity structure along the line from the surface and cross-hole
ERT. Figure 11b is the geological section from the existing boreholes and logging
data along the same line. Comparing the resistivity imaging results with the geolog-
ical section, one can see that the integrated ERT results well map the base of alluvial
overburden and the base of pisolite. The clay contamination within pisolite is also
shown in cross-hole ERT images from the horizontal distance 200–300 m.

6. Conclusions

ERT is a useful near-surface imaging technique, which mainly include data

acquisition, numerical modelling and tomographic inversion. For surface data
acquisition, dipole-dipole, Schlumberger and gradient arrays are applicable for
high-resolution image; particularly, the full-range gradient array may complement
to gradient array for better data coverage and completeness of data information. For
cross-hole data acquisition, pole-pole (A-M), bipole-pole (AM-N) or pole-bipole
(N-AM) and bipole-bipole (AM-BN) can be employed, and the geometry factor and
numerical modelling may be applied for designing efficient and effective arrays and
exam the imaging capability of ERT for specified targets.
GQG approach is a new version of finite-element modelling. It uses Gaussian
abscissae to discrete the model domain and Gaussian weights to compute the vol-
ume integral. Therefore, it is much easier to match arbitrary free-surface topogra-
phy and subsurface interface and computation of the subdomain integrals. It does
not require a small size of element and complex element mesh generator. The
accuracy of modelling depends on the number of Gaussian abscissae in subdomains.
The more abscissae are employed in subdomains, the more accurate modelling
result is generated but costs more computer time.
Tomographic inversion is generally implemented by a standard conjugate gradi-
ent algorithm, which requires to compute the gradient and the Hessian matrix of an
objective function. Two types of objective functions can be applied. One is the l1-
norm and the other is the l2-norm. The former is less sensitive to an outlier in data
but hardly computes the gradient at zero misfit. The latter has no problem to
compute the gradient and the Hessian matrix, and it is much easier to implement a
standard conjugate gradient algorithm. With high-quality data, the two inversions
converge. The field experiments show that surface and cross-hole ERT can be

13
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

applied to map the base of alluvial overburden and the base of pisolite, as well as the
clay contamination within pisolite.

Acknowledgements

The author thanks the Abu Dhabi Education Council for the Award for Research
Excellence (AARE17-273) to financially support this work and greatly appreciates
the ROBE company for sharing the results of ERT imaging experiments.

Author details

Bing Zhou
Department of Earth Science, Khalifa University of Science and Technology,
Abu Dhabi, UAE

*Address all correspondence to: [email protected]

© 2018 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms
of the Creative Commons Attribution License (http://creativecommons.org/licenses/
by/3.0), which permits unrestricted use, distribution, and reproduction in any medium,
provided the original work is properly cited.

14
Electrical Resistivity Tomography: A Subsurface-Imaging Technique
DOI: http://dx.doi.org/10.5772/intechopen.81511

References

[1] Smith NC, Vozoff K. Two- [10] Loke MH, Chambers JE, Rucker DF,
dimensional DC resistivity inversion for Kuras O, Wilkinson PB. Recent
dipole–dipole data. IEEE Transactions development in the direct-current
on Geoscience and Remote Sensing. geoelectrical imaging method. Journal of
1984;GE-22:21-28 Applied Geophysics. 2013;95:135-156

[2] Sasaki Y. Resolution of resistivity [11] Dahlin T, Zhou B. A numerical

tomography inferred from numerical comparison of 2D resistivity imaging
simulation. Geophysical Prospecting. with 10 electrode arrays. Geophysical
1992;40:453-464 Prospecting. 2004;52:379-398

[3] Dahlin T. 2D resistivity surveying for [12] Dahlin T, Zhou B. Multiple-gradient

environmental and engineering array measurements for multichannel
applications. First Break. 1996;14: 2D resistivity imaging. Near Surface
275-283 Geophysics. 2006;4:113-123

[4] Zhe J, Greenhalgh SA, Marescot L. [13] Zhou B, Dahlin T. Properties and
Multi-channel, full waveform and effects of measurement errors on 2D
flexible electrode combination resistivity imaging surveying. Near
resistivity imaging system. Geophysics. Surface Geophysics. 2003:105-117
2007;72:F57-F64
[14] Zhou B, Greenhalgh SA. Cross-hole
[5] LaBrecque D, Miletto M, Daily W, resistivity tomography using different
Ramirez A, Owen E. The effects of noise electrode configurations. Geophysical
on “Occam” inversion of resistivity Prospecting. 2000;48:887-912
tomography data. Geophysics. 1996;61:
538-548 [15] Zhou B, Greenhalgh SA. Rapid 2-D/
3-D crosshole resistivity imaging using
[6] Loke MH, Barker RD. Rapid least- the analytic sensitivity function.
squares inversion of apparent resistivity Geophysics. 2002;67:755-765
pseudosections by a quasi-Newton
Method. Geophysical Prospecting. 1996; [16] Mufti IR. Finite difference
44:131-152 resistivity modelling for arbitrary
shaped two-dimensional structures.
[7] Mauriello P, Monna D, Patella D. 3-D Geophysics. 1976;41:62-78
geoelectric tomography and
archaeological applications. Geophysical [17] Dey A, Morrison HF. Resistivity
Prospecting. 1998;46:543-570 modelling for arbitary shaped two-
dimensional structures. Geophysical
[8] Mauriello P, Patella D. Resistivity Prospecting. 1979;27:106-136
anomaly imaging by probability
tomography. Geophysical Prospecting. [18] Dey A, Morrison HF. Resistivity
1999;47:411-429 modelling for arbitrary shaped three-
dimensional structures. Geophysics.
[9] Zhou B, Greenhalgh SA. Explicit 1979;44:753-780
expressions and numerical calculations
for the Fréchet and second derivatives [19] Mundry E. Geoelectrical model
in 2.5D Helmholtz equation inversion. calculations for two-dimensional
Geophysical Prospecting. 1999;47: resistivity distributions. Geophysical
443-468 Prospecting. 1984;32:124-131

15
Applied Geophysics with Case Studies on Environmental, Exploration and Engineering Geophysics

[20] Spitzer K. A 3-D finite difference [29] Ellis RG, Oldenburg DW. The pole-
algorithm for DC resistivity modelling pole 3-D-resistivity inverse problem: A
using conjugate gradient methods. conjugate-gradient approach.
Geophysical Journal International. 1995; Geophysical Journal International. 1994;
123:903-914 119:187-194

[21] Zienkiewicz OC. The Finite Element [30] Kirsch A. An Introduction to the
Method in Engineering Science. Mathematical Theory of Inverse
London, New York: McGraw-Hill Book Problems. Basel: Springer; 1996
Co, 1971
[31] Schwarzbach C, Ralph-Uwe B,
[22] Coggon JH. Electromagnetic and Klaus S. Two-dimensional inversion of
electrical modelling by the finite direct current resistivity data using a
element method. Geophysics. 1971;36: parallel, multi-objective genetic
132-155 algorithm. Geophysical Journal
International. 2005;162:685-695
[23] Fox RC, Hohmann GW, Killpact TJ,
Rijo L. Topographic effects in resistivity
and induced polarization surveys.
Geophysics. 1980;45:75-93

[24] Pridmore D, Hohmann GW, Ward

SH, Sill WR. An investigation of the
finite element method for electrical and
electromagnetic modelling data in three
dimensions. Geophysics. 1981;46:
1009-1024

[25] Queralt P, Pous P, Marcuello A. 2D

resistivity modelling: An approach to
arrays parallel to the strike direction.
Geophysics. 1991;56:941-950

[26] Zhou B, Greenhalgh SA. Finite

element three-dimensional direct
current resistivity modelling: Accuracy
and efficiency considerations.
Geophysical Journal International. 2001;
145:676-688

[27] Press WH, Teukolsky SA, Vetterling

WT, Flannery BP. Numerical Recipes in
C: The Art of Scientific Computing.
Second ed. England: Cambridge
University, EPress; 1992

[28] Greenhalgh SA, Zhou B, Green A.

Solutions, algorithms and inter-relations
for local minimization search
geophysical inversion. Journal of
Geophysics and Engineering. 2006;3:
101-113

16