Abstract

Groundwater resources along the coast of Cameroon (Kribi–Campo Sub-Basin) are under siege from point and non-point pollution sources, climate change, urbanization and infrastructure development. This situation is made worse by the absence of a water management and development strategy. Managing and monitoring the area's water resources requires an understanding of the groundwater systems, and thus a thorough understanding of the geology. In this study, a 3D geological model was built from electro-seismic data and the structure of the area's aquifer system developed. The aquifer system structure was transferred into Visual MODFLOW Flex and then used to develop a typical hydrogeological model, which will help the management and monitoring of the area's groundwater resources. As more geological data become available, the current model can be updated easily by editing and recomputing. This work is expected to have a positive impact quite quickly on the provision of potable water and on public health.

INTRODUCTION

Geological exploration for water resources, geothermal energy, oil, gas and minerals, and infrastructure projects have not been particularly successful in most African countries. Since such resources are needed for some degree of self-sufficiency and prosperity (growth), it is necessary to understand the geometry of the subsurface using 3D geological modelling (Hassen et al. 2016). 3D modelling is useful and has been used increasingly to integrate available datasets, leading to more realistic representations of given geological settings (Wu et al. 2005).

Water is vital and usually in limited supply for its innumerable uses, the most important of which are domestic purposes. Of earth's water resources, 97.2% are saline (mainly in oceans), and only 2.8% fresh. Groundwater accounts for about 26% of global renewable water resources (Elbeih 2015). Groundwater is the principal source of domestic water for many communities and the most extracted raw material globally. Withdrawal rates are currently estimated at about 982 km3/a (Senthilkumar & Rajakumar 2014). For groundwater resources in any area to be managed and monitored effectively, a good understanding of the geological setting is required. This will enable quantification of the resource and determination of its quality, which is important because of the potential effects of drought, floods, agriculture, urbanization, saline intrusion, and climate change.

A geological model, a spatial representation of the subsurface distribution of sediments and rocks, is traditionally presented in 2D cross-sections; for example, geological maps and profiles, and isopach maps showing the horizontal extent and thickness of subsurface structures (Andersen et al. 2018). 3D models are used increasingly (Sausse et al. 2010; Høyer et al. 2015; Hassen et al. 2016; Malquaire et al. 2017), and are constructed by interpolation from sparse geological and/or geophysical data. Geological and hydrogeological studies in urban areas are also compounded by the presence of anthropogenic deposits, making the subsurface difficult to predict (Andersen et al. 2018).

Collecting the geophysical datasets applied to characterize urban areas is time consuming, and expensive in terms of the survey team and resources mobilized. The methods used include GPR (ground penetrating radar), EM (electromagnetism), DC (direct current) and seismic reflection profiles (Andersen et al. 2018). Simple and robust technologies are now available that save time and money in data collection, while yielding a high-density geophysical dataset. One of these is the electro-seismic (ES) method (Mikhailov et al. 1997). Electro-kinetic prospecting methods are especially promising for detecting materials whose electrical resistivity varies markedly from the background. EM methods may also be used for this, but their spatial resolution is usually less than desirable because of EM's long wavelengths in the earth. Seismic methods have better spatial resolution but respond to smaller contrasts in material properties (White 2005). ES combines the virtues of both approaches – electro-kinetic and seismic – in the context of this study. The datasets used to develop the geological model in this study were collected using the ATS GeoSuite Application installed on a mobile handset and operating under the premise of the ES method (ATSGeoConsultants 2019).

The Kribi-Campo Sedimentary Sub-Basin aquifer system presents many challenges to geological interpretation, mainly because of the limited availability of data in most areas. No geological model has ever been developed for the area, geared towards groundwater resource management. Extensive geological data exist on the area, however, but are not available for academic research (Ntamak-Nida et al. 2010), because they relate to both on- and off- shore petroleum exploration. The lack of sufficient knowledge about the geological and hydrological settings introduces significant uncertainties in projections of the fates of contaminants and public health (Andersen et al. 2018).

Some 250 ES sounding points were surveyed, with an average inter-point profile grid spacing of 200 m inland and variable along the shore. This low-resolution grid spacing was used to provide representative coverage of the area, to extract the hydrological structures and parameters. Some survey points were difficult to access due to the terrain and, in general, it was not possible to use regular, straight grid lines because soundings had to be made along existing tracks and open fields.

The aim of this work was to develop a straightforward procedure for collecting geophysical datasets, and developing geological and hydrogeological models, to enable assessment of aquifer geometries and volumes, and defining their connectivity. This should enable answers to be found for questions concerning water availability, management and monitoring, especially in an urban setting.

STUDY AREA

Site description

Kribi town is in Kribi 2 Sub-Division, of the Ocean Division of the South Region of Cameroon. It is between latitudes 2.940 and 3.007° N, and longitudes 9.90 and 9.95° E, and covers about 13 km2. The Atlantic Ocean bounds it to the west – Figure 1. The area is urban and has several beautiful beaches along approximately 187 km of coastline, which attract many tourists (Nfomou et al. 2004). The local aquifer system in the coastal area is very shallow, as indicated by the depth to the water table in some hand-dug wells. The area has an equatorial climate, so annual rainfall is sufficient. The top soil varies little.

Figure 1

The study area, showing the potential for freshwater resources in the area to be affected by saline intrusion as well as anthropogenic sources.

Figure 1

The study area, showing the potential for freshwater resources in the area to be affected by saline intrusion as well as anthropogenic sources.

Tectonic background

The Kribi-Campo sub-basin comprises an early to mid-Cretaceous series from West Africa's Atlantic coast and is in the Central African equatorial rain forest in Southern Cameroon. It is the smallest coastal basin in Cameroon and forms the southern part of the Douala/Kribi-Campo basin (Ntamak-Nida et al. 2010). The region is tectonically active, and several earthquakes of moderate magnitude (3.4–5.7) have been felt and reported, indicating that the region could be subject to seismic risk and also that the Kribi-Campo Fault (KCF) that crosses it might be active (Owona Angue et al. 2011). The KCF seems to trend in continuation of the offshore Kribi Fracture Zone, which has been inferred as the eastern end of the Ascension Fracture Zone (Ntamak-Nida et al. 2010; Malquaire et al. 2017; Nfor et al. 2018).

Geological and hydrogeological background

The Kribi-Campo sub-basin is sedimentary (Figure 2), and lies across the Cameroon coast between latitudes 2°20′ and 3°20′N, and longitudes 9°15′ and 10°00′E. It covers about 6,195 km2 (Malquaire et al. 2017). The Ntem Complex (Archaean) basement is overlain by Paleoproterozoic rocks of the Nyong Unit. The basement is composed mostly of greenstone belt rocks, Congo Craton charnockites and potassic granitoids (Owona Angue et al. 2011). The sub-basin is bordered to the south by the Campo high, which separates it and the Rio Muni Basin (Iboum Kissaaka et al. 2011). The sub-basin shows three classical stages in its structure: syn-rift (late Jurassic – early Cretaceous); transition (Aptian-Albian) and post-rift (Albian-Recent) (Ntamak-Nida et al. 2010; Iboum Kissaaka et al. 2011).

Figure 2

Geological map of the southern part of the Pan-African North Equatorial Fold Belt in Cameroon, showing the study area within the Sedimentary Basin. SF = Sanaga Fault, CC = Congo Craton, CVL = Cameroon Volcanic Line (modified from (Owona Angue et al. 2011)).

Figure 2

Geological map of the southern part of the Pan-African North Equatorial Fold Belt in Cameroon, showing the study area within the Sedimentary Basin. SF = Sanaga Fault, CC = Congo Craton, CVL = Cameroon Volcanic Line (modified from (Owona Angue et al. 2011)).

ELECTRO-SEISMIC THEORY

Electromagnetic (EM) and mechanical waves are coupled through electro-kinetics, in a porous medium such as the earth's subsurface (White 2005). The pore fluid ions are attracted to ions of the opposite sign in the solid at the pore walls, so that there is an electrical double layer (a thin charged layer coating the mineral grain surfaces) at the pore boundary – the Debye Layer. An electric field acting on this layer moves the ions relative to one another, causing movement of both the fluid and the solid. Conversely, a mechanical wave that moves the fluid and solid relative to each other will create an EM wave (White 2005). Techniques based on electro-kinetics have been proposed for use in sub-surface exploration and field experiments conducted. The experiments were variously seismo-electric and electro-seismic. In seismo-electric experiments, a mechanical (seismic) source is used and an EM wave detected, and vice versa in an electro-seismic experiment (White 2005). While such phenomena can be caused by various different physical mechanisms, they are referred to collectively as ‘electro-seismic’ (Mikhailov et al. 1997).

ES phenomena in porous media have the potential to enable detection of subsurface zones of high fluid mobility and fluid chemistry contrasts, supported by observation, with explicit comparison to full waveform modelling results (Mikhailov et al. 1997). A seismic wave propagating in a medium can induce an electric field or generate an electromagnetic wave.

The ES effects employed (ATSGeoConsultants 2019) look at results of a seismic wave crossing or travelling along the interface between two porous media. When a seismic source crosses an interface, there is ES conversion (Figure 3(a)). When a spherical P-wave crosses an interface, however, the imbalance of the streaming currents induced by the seismic wave on opposite sides of the interface creates a dipole charge separation. The induced dipole radiates an electromagnetic wave, which can be detected by remote antennae (Mikhailov et al. 1997; ATSGeoConsultants 2019).

Figure 3

(a) ES conversion at an interface; (b) generation of an electrical field by a head wave crossing an interface (Mikhailov et al. 1997).

Figure 3

(a) ES conversion at an interface; (b) generation of an electrical field by a head wave crossing an interface (Mikhailov et al. 1997).

When a seismic head wave travels along an interface between media (Figure 3(b)), it causes charge separation across the interface, which induces an electric field that moves along the interface with the head wave and can be detected by ground-dipole antennae when the head wave passes under them (Mikhailov et al. 1997; ATSGeoConsultants 2019). The electric fields are transformed into a set of time-varying potential differences, which are transmitted through the signal cables to the mobile device running the GeoSuite Application, and recorded and stored for processing (ATSGeoConsultants 2019).

METHODOLOGY

ES equipment consists of a single mobile seismic source – any form of seismic impulse generator, including a low energy hammer and plate, a medium energy drop-weight, or a high-power explosive source such as a buffalo gun. The seismic source is placed over the field point and an operating area marked around it, and two stainless steel pins inserted into the ground as the ground–dipole antenna. A signal cable connected to the pins enables recording and is calibrated to the mobile device's audio port. The seismic waves are produced by striking the ground, a process repeated from 5 to 20 times at each location (referred to as data stacking).

Data processing includes filtering, extraction, correlation, stacking, and quality evaluation, all of which are completed automatically to produce up to 25 different datasets. The datasets include both hydraulic conductivity tomography (ESKT) and interface tomography (ESIT), which enable aquifer characterization, as well as others. Amongst the parameters used to determine the quality of the datasets are the correlation indexes. The spatial distribution of data points with their correlation indexes, which describe their reliability, is shown in Figure 4. In a study involving drilling, points with correlation index values below 50% would probably be re-drilled. As there was no drilling in this study, no repeated data collections were carried out at data points with correlation index below 50%. Instead, only those with correlation index values exceeding 50% were analyzed further. Water table depths were obtained from 61 hand-dug wells and shallow boreholes in the survey area, and the resulting dataset used to construct a ‘layer’ comprising the probable phreatic aquifer.

Figure 4

ES survey points with their correlation index values.

Figure 4

ES survey points with their correlation index values.

DEVELOPMENT OF THE AQUIFER SYSTEM AND GEOLOGICAL MODEL FROM THE DATASETS

The ESKT dataset effectively describes the aquifer systems under a site in 1-, 2-, or 3- D. The ESKT data are normalized to give relative hydraulic conductivity values ranging from 0 to 100%. In other words, the data are not calibrated to an absolute field reference point, at which the hydraulic conductivity with depth is known. The maximum normalized value is obtained by assigning 100% to the maximum recorded hydraulic conductivity value and normalizing all other values against this (ATSGeoConsultants 2019), thus, identifying aquifers' locations. The resulting normalized hydraulic conductivity values were classified and each class was then assigned a given lithology that could possibly correspond to the hydraulic conductivity values of that class. This classification is shown in Table 1, together with the probable lithologies used to construct a lithologic model. The classified dataset was then imported into Leapfrog Geo software (ARANZ 2018), and the aquifer model developed.

Table 1

ESKT value classes and anticipated lithologies used to construct a lithological model for the site

ESKT Value Classes (%)Probable Lithology
[100,60] Sand 
[60,30] Fine sand 
[30,10] Clayey and/or silty sand 
[10,0] Silt and/or clay 
ESKT Value Classes (%)Probable Lithology
[100,60] Sand 
[60,30] Fine sand 
[30,10] Clayey and/or silty sand 
[10,0] Silt and/or clay 

The ESIT dataset is obtained by processing ES data, and indicates the positions of interfaces between formations with differing electrical and elastic properties. The interfaces are detected on the basis of the amplitude of the seismic reflection coefficient, and all interfaces are detected irrespective of its amplitude. This highlights geological features that might be missed if viewed in the context of reflector strength. The data show the reflector's total gradient polarity (Table 2). Blue responses indicate negative polarity and red positive polarity reflectors (ATSGeoConsultants 2019). Visual inspection and the resulting combination of some interfaces/media led to development of the geological model.

Table 2

ESIT values and associated interpretation used in constructing the geological model

Interface tomography values (%)Reflector's polarity
−100 Blue response 
Material medium with normal wave propagation (uniform acoustic impedance) 
100 Red response 
Interface tomography values (%)Reflector's polarity
−100 Blue response 
Material medium with normal wave propagation (uniform acoustic impedance) 
100 Red response 

RESULTS

Dataset calibration

The datasets analyzed were converted to resemble drilling data and imported into the Voxler software (Golden Software 2019b). The water table surface grid was prepared in Surfer 16 software (Golden Software 2019a) and also imported into Voxler. The resulting model shows a system that fits together well, comprising a shallow unconfined aquifer along that length of coast; that is, the datasets reflect the geological and hydrogeological formations (see Figure 5).

Figure 5

ES dataset calibration in Voxler software with the water table grid, reflecting the existing shallow aquifer system. The slope falls to sea level from about 19 m.

Figure 5

ES dataset calibration in Voxler software with the water table grid, reflecting the existing shallow aquifer system. The slope falls to sea level from about 19 m.

Geological model construction

The geological model component of Leapfrog Geo software incorporates five sub-components – (model) boundary, fault system, lithologies, surface chronology and output volume. The boundary inherits the clipped boundary set at topography level, which defines the upper boundary automatically. The fault system and lithologies sub-components are self-explanatory, although it is noted that the lithology field can be over-ridden manually if the dataset contains no lithologic information. The surface chronology component deals with all contact surfaces, which are arranged in chronological order, youngest first. It is these surfaces and their chronology that determine how the volume inside the model is divided into lithologic units. The output volumes comprise all those generated in building the model, in chronological order from youngest to oldest.

Leapfrog Geo can create a number of different surface types. The stratigraphic sequence was chosen for two reasons. Firstly, it should be used when there is a series of continuous layers that can be modeled separately. Secondly, it is known that this is a sedimentary basin, so that the lithologies can be modelled as layered horizontal blocks. Figure 6 shows the datasets in Leapfrog Geo and Figure 7 the ensuing geological model.

Figure 6

ES datasets converted to pseudo-drilling data and uploaded into Leapfrog Geo software (vertical exaggeration = 10×).

Figure 6

ES datasets converted to pseudo-drilling data and uploaded into Leapfrog Geo software (vertical exaggeration = 10×).

Figure 7

Geological model generated in Leapfrog Geo from ES datasets (vertical exaggeration = 10×).

Figure 7

Geological model generated in Leapfrog Geo from ES datasets (vertical exaggeration = 10×).

Aquifer system

An aquifer system was developed from the geological model. The unconfined aquifer modelled is shown in Figure 8 as the Upper Sandy Aquifer. The latter was compared with that generated from the water table grid – see Figure 9 – and the two layers generally coincide, with minor discrepancies in regions of extreme extrapolation. The final aquifer system determined for the study area is shown in Figure 10.

Figure 8

Unconfined aquifer as modelled using Leapfrog Geo and the ES datasets (vertical exaggeration = 10×).

Figure 8

Unconfined aquifer as modelled using Leapfrog Geo and the ES datasets (vertical exaggeration = 10×).

Figure 9

Unconfined aquifer as modelled in Leapfrog Geo compared with the water table grid developed from measured values (see Figure 5) (vertical exaggeration = 10×).

Figure 9

Unconfined aquifer as modelled in Leapfrog Geo compared with the water table grid developed from measured values (see Figure 5) (vertical exaggeration = 10×).

Figure 10

Aquifer system within Kribi town as derived from analyses of ES data.

Figure 10

Aquifer system within Kribi town as derived from analyses of ES data.

Fence sections

Two fence sections were drawn in Leapfrog Geo to show the lithologic variations within the geological model. They run along the landward side of the study area (Figure 11(a) and 11(b)) and the other across the area towards the ocean (Figure 12(a) and 12(b)).

Figure 11

(a) Landward side fence section of the study area through the geological model; (b) Fence section layout; loess (top soil) occurs sparingly but is commoner along the landward side.

Figure 11

(a) Landward side fence section of the study area through the geological model; (b) Fence section layout; loess (top soil) occurs sparingly but is commoner along the landward side.

Figure 12

(a) Fence section across the study area approximately perpendicular to the coast; (b) Fence section layout. The loess (top soil) is completely absent towards the ocean because the coastal zone is sandy.

Figure 12

(a) Fence section across the study area approximately perpendicular to the coast; (b) Fence section layout. The loess (top soil) is completely absent towards the ocean because the coastal zone is sandy.

Hydrogeological model

The conceptual system developed using Leapfrog Geo was converted into surfaces and exported as Surfer grid files to build the hydrogeological model in Visual MODFLOW Flex software (Hydrogeologic Waterloo 2015). The latter, a basic high-level representation of the hydrogeological system, will form the basis of other numerical models developed in Visual MODFLOW Flex, since it is grid and simulator independent.

The flow condition chosen for the groundwater model was saturated (constant density) flow, with no transport objective defined. After the model structure had been defined, the model's property zones were defined using the synthetic hydraulic conductivity values given in Table 3. The hydrogeological model is shown in Figure 13.

Table 3

Hydraulic conductivity values used to define property zones for aquifers and aquitards in the Visual MODFLOW Flex-based hydrogeological model

LayerKx, Ky and Kz
Surface layer 0.001 
Upper, sandy, aquifer 0.002 
Clay/silt aquitard 0.00001 
Lower, fine sand, aquifer 0.0015 
LayerKx, Ky and Kz
Surface layer 0.001 
Upper, sandy, aquifer 0.002 
Clay/silt aquitard 0.00001 
Lower, fine sand, aquifer 0.0015 
Figure 13

Hydrogeological model developed in Visual MODFLOW Flex to be used for grid generation and conversion to groundwater numerical models.

Figure 13

Hydrogeological model developed in Visual MODFLOW Flex to be used for grid generation and conversion to groundwater numerical models.

CONCLUDING REMARKS

The idea behind this study was conceptualization of the procedure for generating hydrogeological models, based on ES dataset collection. Construction of a hydrogeological model starts with site and data acquisition point definition, data collection and processing, and using software to simulate the 3D hydrogeological system; for example, Leapfrog Geo and Visual MODFLOW Flex. This method of geophysical exploration is easily implemented in mobile devices and the groundwater exploration easily handled, so good understanding of the procedure will limit the time taken for groundwater studies. Although this model was not calibrated using down-hole log data, the high-density geophysical sets gave coherence in interpolating the surfaces, which could be helpful in groundwater modelling.

In this study the attempt was made to show how easy it can be to carry out groundwater studies. Due to the maiden nature of the work in general and the paucity of data in the study area, the work suffered from shortcomings that included, particularly:

  • The model domain required further refinement by re-evaluating all sounding points with correlation index values (Figure 4) below 50% – they were excluded from analyses – to further reduce interpolation errors;

  • Borehole data should be used to validate the ES dataset in the area. In other words, the lithology obtained from the datasets should conform with that determined from drilling logs and hydraulic conductivities. If this is done, the method can be reliably exported to other areas.

ACKNOWLEDGEMENTS

The datasets generated and/or analyzed during this study are available in the 4TU.Center for Research Data repository, http://researchdata.4tu.nl/.

REFERENCES

REFERENCES
Andersen
T. R.
Poulsen
S. E.
Thomsen
P.
Havas
K.
2018
Geological characterization in urban areas based on geophysical mapping: a case study from Horsens, Denmark
.
Journal of Applied Geophysics
150
,
338
349
.
ARANZ Geo Leapfrog
2018
Leapfrog geo (Version 4.0.0)
.
ATSGeoConsultants
2019
Aquatronic Solutions Geo-Consultants and Developers of Geophysical Software and Tools
.
Elbeih
S. F.
2015
An overview of integrated remote sensing and GIS for groundwater mapping in Egypt
.
Ain Shams Engineering Journal
6
(
1
),
1
15
.
Golden Software
2019a
Surfer® (Version 16)
.
Available from: www.goldensoftware.com
Golden Software
2019b
Voxler® (Version 4.3.771)
.
Available from: www.goldensoftware.com
Høyer
A.-S.
Jørgensen
F.
Sandersen
P.
Viezzoli
A.
Møller
I.
2015
3D geological modelling of a complex buried-valley network delineated from borehole and AEM data
.
Journal of Applied Geophysics
122
,
94
102
.
Hydrogeologic Waterloo
2015
Visual MODFLOW Flex (Version 5.0)
.
Iboum Kissaaka
J. B.
Ntamak-Nida
M. J.
Ngueutchoua
G.
Mvondo Owono
F.
Djomeni Nitcheu
A. L.
Fowe Kuetche
P. G.
Bourquin
S.
2011
Syn-rift tectonic markers from outcrops and offshore seismic data of the southern part of Kribi-Campo sub-basin (Cameroon, West African margin)
.
Unidade de Geologia, Hidrogeologia e Geologia Costeira
99
,
35
42
.
Malquaire
K. P. R.
Louise
O. A. M.
Nfor
N.
Eliezer
M.-D.
2017
3D modelling from new and existing gravity data of an intrusive body in the Northern Part of Kribi-Campo Sub-Basin in Cameroon
.
International Journal of Geosciences
8
(
08
),
984
.
Mikhailov
O. V.
Haartsen
M. W.
Toksöz
M. N.
1997
Electroseismic investigation of the shallow subsurface: field measurements and numerical modeling
.
Geophysics
62
(
1
),
97
105
.
Nfomou
N.
Tongwa
A. F.
Ubangoh
U. R.
Bekoa
A.
Metuk
N. J.
Victor
H. J.
2004
The July 2002 earthquake in the Kribi region: geological context and a preliminary evaluation of seismic risk in southwestern Cameroon
.
Journal of African Earth Sciences
40
(
3–4
),
163
172
.
Nfor
N.
Clotilde
O. A. M. L.
Malquaire
K. P. R.
Eliezer
M. D.
Achille
L. T. J.
Jöel
P. K.
2018
A comparative geophysical study of the sedimentary-metamorphic contacts in the Douala and Kribi-Campo Sub-Basins of Cameroon
.
Trends Journal of Sciences Research
3
(
1
),
33
51
.
Ntamak-Nida
M. J.
Bourquin
S.
Makong
J.-C.
Baudin
F.
Mpesse
J. E.
Ngouem
C. I.
Abolo
G. M.
2010
Sedimentology and sequence stratigraphy from outcrops of the Kribi-Campo sub-basin: Lower Mundeck Formation (Lower Cretaceous, southern Cameroon)
.
Journal of African Earth Sciences
58
(
1
),
1
18
.
Owona Angue
M.
Nguiya
S.
Nouayou
R.
Tokam Kamga
A.
Manguelle-Dicoum
E.
2011
Geophysical investigation of the transition zone between the Congo Craton and the Kribi-Campo sedimentary basin (Southwestern Cameroon)
.
South African Journal of Geology
114
(
1
),
145
.
Sausse
J.
Dezayes
C.
Dorbath
L.
Genter
A.
Place
J.
2010
3D model of fracture zones at Soultz-sous-Forêts based on geological data, image logs, induced microseismicity and vertical seismic profiles
.
Comptes Rendus Geoscience
342
(
7–8
),
531
545
.
Senthilkumar
M.
Rajakumar
T. J.
2014
Groundwater quality assessment in Lower Tamirabharani River Basin Tamil Nadu, India
.
International Research Journal of Earth Sciences
2
(
10
),
8
15
.
White
B. S.
2005
Asymptotic theory of electroseismic prospecting
.
SIAM Journal on Applied Mathematics
65
(
4
),
1443
1462
.
Wu
Q.
Xu
H.
Zou
X.
2005
An effective method for 3D geological modeling with multi-source data integration
.
Computers & Geosciences
31
(
1
),
35
43
.