Integrating data from remote sensing, geology and gravity for geological investigation in the Tarhunah area, Northwest Libya

Abstract

The present work deals with the integration of remote-sensing, surface-geology and gravity-survey data to improve the structural knowledge of the Tarhunah area, northwest Libya. Geological information and remote-sensing data provided information about the surface structure. A gravity survey was conducted to decipher the subsurface structure. The results revealed that a basin having a width of 39 to 48 km trends NE. A two-dimensional (2-D) schematic model shows that the basin gradually deepens toward the southwest. Faults determined from a horizontal gradient, tilt derivative, and Euler deconvolution show a depth range of 2.5 to 7.5 km. The integration and interpretation of the results indicate that volcanic activity was related to the tectonic activity of an anticlinal structure called the Jabal Uplift.

Keywords:

• Libya
• Tarhunah area
• remote sensing
• gravity
• integration

1. Introduction

The integration of geological, remote-sensing, and geophysical data has the potential of constraining quantitative details to hitherto unknown areas and reduces the ambiguity of geological interpretation (Lunden et al. 2001, Lamontagne et al. 2003, Chen and Zhou 2005, Yassaghi 2006). The collection of new data and the analysis of a previously little-known area, the Tarhunah area, can improve scientific fields for real data with the aim of formulating new geophysical models.

The Tarhunah area mostly has a cover of basalt sheets with scattered basalt cones and phonolite intrusions, which in some ways affect the accessibility of the area. The rest of the area is complex, with a fragmented distribution of different land-cover types. Several authors studied the geological setting of the Tarhunah area in the 1960s and 1970s (Burollet 1963, Jordi and Lonfat 1963, El Hinnawy and Cheshitev 1975, Mann 1975, Antonovic 1977, Zivanovic 1977). This study emphasises the importance of the synergistic use of traditional geological investigations, remote-sensing techniques, geophysical data and different analytical and numerical modelling methodologies for further geological investigations.

The application of remote sensing for the mapping of regional structure has a long tradition worldwide (Bonham Carter 1989, Harris 1991, Fraser et al. 1997, Zeinalov 2000, Leech et al. 2003, Cengiz et al. 2006, Raharimahefa and Kusky 2006, Pereira et al. 2008). In investigating this area an enhanced thematic mapper plus (ETM + ) image was processed and analysed to identify surface structural features. Filtering operations and principal component analysis (PCA) were implemented to enhance the visual interpretation for revealing geological lineaments, thus increasing the possibilities for extracting useful geological information. Four geological maps were prepared for showing the lithology, stratigraphy and structure of this study area. Geological data were combined with remote-sensing imagery to improve our knowledge of the surface structure.

Geophysical data record the physical characteristics of subsurface geological features, and the subsurface structures can be retrieved and interpreted (Ayala et al. 2003, Yegorova et al. 2004, Bilim 2007, Tašárová 2007). Gravity data were used in this study to give a general picture of the subsurface structure of the study area. Power-spectrum analysis provides information about basement depths. Other techniques were used in an attempt to delineate the borders and the depth of the basement. Horizontal gradient and tilt derivative (TDR) techniques were used to locate the main geological boundaries of the explored basin. The analysis and interpretation of gravity data revealed a basin trending NE. The estimated faults vary in depth from 2.5 to 7.5 km. The geological map indicates that some of these faults might be interpreted as vertical extensions of surface faults. Two-dimensional (2-D) forward modelling showed the basin gradually to deepen toward the SW. The assessment of lineament length, density and trends provides useful information about the tectonic evolution of the study area.

The integration of results shows a close relationship between volcanic activity and an anticlinal structure called the Jabal Uplift.

2. The study area

The study area (Figure 1) lies within the southernmost region of the city of Tarhunah in the northwestern part of Libya, bounded by longitudes 13°04′ E to 14°00′ E and latitudes 31°45′ N to 32°23′ N. It covers a surface area of approximately 6200 km².

Figure 1.  The study area, overlain on the tectonic map of Libya, modified from Rusk (2001).

The study area has an elevation that ranges approximately between 300 and 900 m. Low topographic areas are covered mostly by fluvioeolian deposits, whereas most of the area is covered by basalt sheets, forming flat, broad plains.

The study area contains exposures of sedimentary rocks ranging in age from Middle and Late Triassic to Quaternary. The Aziziya Formation (Middle–Upper Triassic) is partially exposed in the northwestern part of the study area. The Abu Shaybah Formation (Upper Triassic) is present in the northwestern part of the study area and consists of sandstone alternating with bands of clays and scattered limy bands (Mann 1975). The Abu Ghaylan Formation (Upper Triassic–Middle Jurassic) appears as a narrow limestone ridge in the northwestern part of the study area (El Hinnawy and Cheshitev 1975).

The Sidi as Sid Formation (Cretaceous) is exposed through the northwestern part of the study area. It comprises two successive rock units, the Ain Tobi Limestone and Yafrin Marls (Mann 1975). The Nalut Formation (Cretaceous) is one of the most widely distributed rock units through the northern and northeastern parts of the study area. The Qasr Tigrinnah Formation (Cretaceous) is partly exposed in the western part of the study area. Jordi and Lonfat (1963) included the Qasr Tigrinnah as a member of the Mizdah Formation (Upper Cretaceous), which is widely exposed in the southern part of the study area. As recognised by Jordi and Lonfat (1963), the Mizdah Formation comprises three members, the Tigrinnah Marl, Mazuzah Limestone and Thala Member. El Hinnawy and Cheshitev (1975) restricted the Mizdah Formation to two members, the Mazuzah and Thala, and raised the Tigrinnah Marl to the status of a formation.

Volcanic rocks (Tertiary–Quaternary) are common in the central part of the study area. These volcanic rocks consist of three types: phonolite and trachyte intrusions, basalt cones and basalt flows. The phonolite and trachyte intrusions, partly observed in the northern part of the study area, exhibit topographic features of elevated, black, conical and some crater-shaped hills (Zivanovic 1977). Basalt cones are scattered over the basalt sheet in this area. They have a circular shape, with a raised rim, and invariably are intersected by lava extrusions (Antonovic 1977). The age of the basalt sheet was interpreted by Piccoli (1971) as early Eocene to Pliocene, whereas the basalt at Wadi Ghan was interpreted by Christie (1966) as early Quaternary.

Quaternary sediments were observed throughout the study area, especially in the north. These Quaternary sediments were divided according to lithostratigraphic and genetic classification into four types: the Qasr Al Haj Formation (Pleistocene), fluvioeolian deposits (Holocene), eolian deposits (Holocene), and recent wadi deposits (Holocene) (Figure 2).

Figure 2.  Geological map of the study area and a rose diagram of surface faults (I.R.C., 1975a, b, 1977a, b) (Available in colour online).

At the end of the Cretaceous or in the early Tertiary the Gharyan-Tarhunah strip was uplifted and then formed an anticlinal swelling oriented NW–SE. According to Goudarzi (1970), the main structural feature is a vast anticlinal swelling, the Gharyan-Yafran Arch, which trends NW–SE. In reference to Klitzsch (1971), the uplift may have taken place in the third period of structural development in North Africa, the time (late Paleozoic and Mesozoic) when grabens and mountain ranges, oriented mostly NE–SW, began to form, and block faulting appeared along marginal parts of the uplifts (during the Jurassic or at the Triassic–Jurassic boundary). Burollet (1963) assumed that Jabal Nefusa was uplifted during the Hercynian orogeny. The commonly accepted hypothesis is that the structural development of Jabal Nefusa was terminated toward the end of the Cretaceous (Mann 1975). Lipparini (1940) believed that the uplift was accompanied by NW–SE faulting. All authors are in agreement as to the diminishing tectonic activity in the eastern part of Jabal Nefusa. The uplifted area that wedges out in the north could be regarded as disappearing uplifts buried by extensive volcanic flows (Gray 1971).

The E–W faults are recorded mainly from the subsurface and are represented by the Al Aziziyah and Coastal faults. The Aziziyah fault was active as a large right-lateral shear during the Atlas compressional phase and probably represents a basement feature (Burollet 1991, Tawadros 2001).

3. Data processing

A Landsat ETM+ image was processed and interpreted to identify geological features and produce a regional distribution of geological formations in the study area. The gravity data were used to validate the structures on the geological map and remotely sensed imagery from a point of potential field. Furthermore, the subsurface information for the geological structures can be retrieved from the gravity data. The three different data types were integrated for the distribution and depth information of the subsurface geological structures, using the geographic information system (GIS) technique.

3.1 Remote-sensing data

The main objective of the remote-sensing analysis was to extract and map lineaments. The ETM+ program was tested for identifying and mapping geological lineaments in the entire area. Filtering operations and PCA were used to improve the visual appearance of the satellite image and prepare the image for measurement of the features and structures. The extracting criteria for lineaments and visual interpretation were based mainly on the image characteristics (tone and texture), lithological boundaries (rock units) and the geomorphological features (drainage patterns).

Edge enhancement was performed by one of the two types of spatial filtering called a high-pass filter, which is designed to emphasise the detailed high-frequency components of an image and de-emphasise the more general low-frequency information (Lillesand and Kiefer 2000). Edge enhancement is performed by detecting edges and then either by adding these back into the original image to increase contrast in the vicinity of an edge or by highlighting edges using saturated (black, white or colour) overlays on borders (Richards 1993). An edge represents discontinuities or a sharp change in the grey-scale value of a particular pixel at a point that might have some interpretation in terms of geological structure or relief (Mather 1993). Directional gradient Prewitt and Sobel filters (Suzen and Toprak 1998, Gonzalez and Woods 2002, Mavrantza and Argialas 2002, Chanda and Majumdar 2003) were applied to the image in the four main directions—NW–SE, NE–SW, N–S and E–W—to emphasise spatial frequency and image contrast. The lineaments extracted and mapped from the Sobel and Prewitt filters are considerably different in number. From the analysis 207 lineaments were detected by the Sobel filter, whereas 173 lineaments were detected by the Prewitt filter. A comparison of the two maps indicates that most of the additional lineaments from the Sobel filter are short and inhomogeneous over the study area, with a mean lineament length of 1.4 km. The maximum lineament length is 4.8 km. However, the lengths of these lineaments are less than the expected values.

PCA allows the redistribution of data in the original channels by means of a linear transformation of channel variables between the same number of new channels in such a manner that the images in the first three principal component images must contain the most information in the N-band. The effectiveness of the PCA approach in identifying and analysing lineaments was investigated by Canas and Barnett (1985), Qari (1991), Leech et al. (2003) and Won-In and Charusiri (2003). Nama (2004) argues that the lineaments can be easily identified by using PCA for the ETM+ image, which removes redundant information from visible and NIR multispectral data. The three first PCA bands were produced using all bands except the thermal and panchromatic bands (Figure 3). The first three principal components contain 96.03% of the total variance of the ETM+ image. In all, 144 lineaments were extracted and mapped using PCA images. The total number of lineaments that were extracted and mapped by using the two remote-sensing techniques was 283. The mean length of all lineaments is 1.7 km, and the maximum length is approximately 6 km. All lineaments were compiled on a reinterpreted geological map on the basis of the ages of the geological formations (Figure 4).

Figure 3.  False-colour composite of PC1 (red), PC2 (green) and PC3 (blue). Dashed black lines indicate formation boundaries. Solid black lines indicate basalt cones. Solid yellow lines indicate phonolite intrusions.

Figure 4.  Geological formations by age, showing lineaments extracted from remote-sensing data. A rose diagram is also shown. (Available in colour online).

3.2 Gravity data

A gravity survey of the Tarhunah area was conducted in an attempt to delineate its subsurface structure. A density of 2,670 kg/m³ was used to compute a Bouguer anomaly map of this area. Figure 5 shows that the area is characterised by negative gravity values ranging between –36 and –8 mGal, increasing in value in the north-northwestern, northeastern, and eastern parts of the map area. These areas are composed mainly of limestone, dolomite and sandstone. It can be seen that the gravity map shows a graben system trending NE in association with the tectonic features of the area. In this study we have applied an upward-continuation filter with a distance of 500 m after examining the spectral analysis of the data and the results of the matched filter so as to reduce the effect of noise.

Figure 5.  Gravity map of the study area (after upward continuation). (Available in colour online).

3.2.1 Energy spectrum analysis

A power spectral analysis yielded information on the depth of a significant density contrast in the crust and on the crustal structure (Kadirov 2000, Gotze and Krause 2002). Spector and Bhattacharyya (1966) studied the energy spectrum, which was calculated from different 3-D model configurations. Then, Spector and Grant (1970) studied the statistical ensemble of 3-D maps and concluded a general form of the spectrum

1

where

E: total energy (frequency Hz).

R: radial wave number (cycles/m).

θ: azimuth of radial wave number (degree).

〈〉: expresses ensemble average.

h: depth (meter).

H: depth factor.

S: horizontal size (width) factor.

C: vertical size (thickness or depth extent) factor.

a, b: parameters related to the horizontal dimensions of the sources.

t, Φ: parameters related to the vertical depth extent of the source.

It is clear that only three factors (H, S and C) are functions of the radial frequency r, thus in the case of profile form, equation (7) can be written as

2

where : the average depth, half width and thickness of the source ensemble.

This equation demonstrates that contributions from the depths, widths and thicknesses of the source ensemble can affect the shape of the energy spectral decay curve.

Figure 6(a) shows the energy power spectrum of the gravity data of the Tarhunah area. Figure 6(b) shows the depths of the gravity sources. Two different gravity contributions can be seen: a deep-seated source at a depth of 10 km, and a shallow source at a depth of 2 km.

Figure 6.  (a) Radially averaged power spectrum of the Bouguer anomalies of the study area; (b) gravity depth estimation based on linear segments of the energy decay curve.

3.2.2 Horizontal gradient technique

The horizontal gradient method has been used intensively to locate boundaries of density contrast from gravity or magnetic data (Ma et al. 2006, Bilim 2007, Pilkington 2007). Generally, the horizontal gradient of the gravity anomaly caused by a tabular body tends to overlie the edges of the body if the edges are vertical and isolated from each other (Cordell 1979, Cordell and Grauch 1985).

The best advantage of the horizontal gradient method is that it is least susceptible to noise in the data, because it only requires the calculations of the two first-order horizontal derivatives of the field (Phillips 1998). This method is also robust in delineating both shallow and deep sources, in comparison with the vertical gradient method, which is useful only in identifying shallower structures. The amplitude of the horizontal gradient (Cordell and Grauch 1985) is expressed as

3

where (∂g/∂x) and (∂g/∂y) are the horizontal derivatives of the gravity data. Figure 7 shows the horizontal-gradient map of the study area. High gradient values were observed at locations that could be the boundaries of the basin.

Figure 7.  Horizontal-gradient map of the upward-continued gravity data. Circles indicate the maxima of the horizontal gradient at which contacts could be located. These peaks are at the contacts or faults. In some places they are at the borders of the basin. (Available in colour online).

Once the horizontal-gradient map is obtained, the maximum peaks can be determined, using the criteria of Blakely and Simpson (1986) shown as black circles in Figure 7.

3.2.3 Tilt derivative (TDR) technique

The TDR is used to enhance and sharpen the potential field anomalies (Verduzco et al. 2004, Cooper and Cowan 2006). The advantage of the TDR is that its zero contour line is located on or close to the contacts. The operating definition (Miller and Singh 1994) is

4

where

5

The results of the TDR are presented in Figure 8.

Figure 8.  Tilt derivative of gravity data. Solid black line is the zero contour line, which matches well with the basin borders and also with the horizontal gradient map (Figure 7) of the gravity data. (Available in colour online).

3.2.4 Euler deconvolution technique

The Euler deconvolution technique is used to estimate depth and location of the magnetic gravity source anomalies. Euler deconvolution has come into wide use for automatic estimates of field source location and depth (Bournas et al. 2003, Shepherd et al. 2006, Chennouf et al. 2007, Li et al. 2008). This method was established by Thompson (1982) and has been applied essentially to real magnetic data along profiles. Reid et al. (1990) followed up a suggestion in Thompson's paper and developed the equivalent method, which operates on gridded magnetic data. The application of Euler deconvolution to gravity data has been carried out by several authors—e.g., Wilsher (1987), Corner and Wilsher (1989), Klingele et al. (1991), Marson and Klingele (1993), Fairhead et al. (1994) and Huang et al. (1995).

The 3-D equation of Euler deconvolution given by Reid et al. (1990) is

6

Equation (6) can be rewritten as

7

where (x _o , y _o , z _o ) are the position of a source whose gravity anomaly is detected at (x, y, z), β is the regional value of the gravity (background) and η is the structural index (SI), which can be defined as the rate of attenuation of the anomaly with distance. The SI must be chosen according to a prior knowledge of the source geometry. For example, in the gravity case, SI = 2 for a sphere, SI = 1 for a horizontal cylinder, SI = 0 for a fault and SI = –1 for a contact (FitzGerald et al. 2004). The two horizontal gradients (∂g/∂x, ∂g/∂y) and vertical (∂g/∂z) derivatives are used to apply the Euler method.

By considering four or more neighbouring observations at a time (an operating window), source location (x _o , y _o , z _o ) and β can be computed by solving a linear system of equations generated from Equation (7). Then, by moving the operating window from one location to the next over the anomaly, multiple solutions for the same source are obtained. In our study the Euler method has been applied to the gravity data using a moving window of 0.01 km×0.01 km, and the grid cell size of the gravity data is 750 m. We have assigned several structural indices values and found that an SI of zero gives good clustering solutions. Reid et al. (1990, 2003) and Reid (2003) presented a structural index equal to zero for the gravity field for detecting faults. Results of the Euler deconvolution of gravity data are shown in Figure 9.

Figure 9.  Euler deconvolution of the upward continued gravity data (posted circles). The solutions of Euler are consistencies with the tilt derivative (TDR; zero contour line in black). At the zero contour line of the TDR, you can easily find the Euler solutions for faults/contacts. Background is image of gravity map.

The interpretation of Euler solutions indicates that the NE–SW, NW–SE, E–W, and N–S trends characterise the structural setting of the Tarhunah area. The depth of faults ranges from 2.5 to 7.5 km. These contacts were interpreted as faults and plotted and statistically analysed (Figure 9). The rose diagram (Figure 10) indicates that four major fault patterns (N–S, E–W, NE–SW, NW–SE) characterise the study area.

Figure 10.  Rose diagram of faults extracted from gravity data.

3.2.5 Two-dimensional (2-D) forward modelling

For quantitative gravity interpretation, a 2-D forward modelling method (Talwani et al. 1959) was applied to the gravity data. Two profiles were selected (A–A′ and B–B′, Figure 5). Both profiles were extracted from the Bouguer anomaly map to understand the subsurface structure and for identifying the borders of the basin and its depth. The density contrast between the two modelled layers is –0.2 g/cm³, as shown in Figure 11. The density contrast is –0.2 g/cm³.

Figure 11.  The 2-D forward gravity model, cross sections A–A′ and B–B′.

3.3 Data integration

The GIS-based method was used to integrate raster and vector results extracted from geological field work, and remote-sensing and geophysical data analysis. The integrated analysis of the different data types was applied to find the mutual relationship as a base for geological analysis. The ultimate application of this technique should result in a more detailed and accurate geological interpretation. In particular, the integration of such data makes it possible to define both the relations of geological formations with the structural context (Rowland and Duebendorfer 1994), and provide new opportunities for studying the tectonic evolution of the study area (Figure 12).

Figure 12.  Basin trending NE. The estimated depth of faults ranges from 2.5 to 7.5 km. The borders of the basin are in black solid lines, based on the results of horizontal gradient, tilt derivative and Euler deconvolution. Dashed black lines indicate the interpreted faults from gravity data. Solid blue lines are the faults from the geological map. Solid green lines are the lineaments from the satellite image. Violet spots are basalt cones, and orange spots are phonolite intrusions. (Available in colour online).

4. Results

In accordance with the geological maps, the study area is mostly overlain by a basalt flow, with scattered basalt cones and phonolite intrusions. The faults trend predominantly NW–SE, whereas the E–W faults are subordinate. The NE–SW and N–S trends are not primary. A dominant NW–SE fault trend is evident in Triassic rocks. In Cretaceous rocks the dominant fault trend is NW–SE, with a subordinate E–W trend. Faults in Tertiary rocks are oriented predominantly NW–SE. In Quaternary sediments the predominant fault trend is NW–SE. In the Triassic and Cretaceous rocks the frequency and lengths of the NW–SE faults indicate reactivated faulting.

Lineaments extracted from remote-sensing data revealed different trends, but the main trend is NW–SE, parallel to the main tectonic line of the Jabal Uplift. The E–W lineaments represent a secondary trend. Other significant directions are N–S and NE–SW. The N–S trend is more conspicuous. The remaining directions are statistically unimportant. Variations in lineament trends can be explained by gentle anticlinal or broad synclinal forms. Lineaments extracted from remote-sensing data were divided into four groups, based on the age of the geological formations. Lineaments in Triassic rocks are dominantly NW–SE, whereas NE–SW is the subordinate direction. The prevailing trend of Cretaceous rock lineaments is NW–SE, and the second dominant trend is NE–SW. Lineaments of Tertiary rocks dominate in an E–W direction, with the NW–SE direction being subordinate. Quaternary sedimentary lineaments trend dominantly NW–SE and are the longest lineaments in these sediments.

Reconnaissance investigations and study of gravity data revealed a NE–trending basin with a width of 39 to 48 km. Faults extracted from horizontal gradient, tilt derivative, and Euler deconvolution methods show a dominant N–S trend with subordinate E–W faults. The depths of these faults range from 2.5 to 7.5 km. According to the gravity data and 2-D models, the basin gradually deepens toward the southwest. The faults extracted from gravity data can be placed into four groups. The first group (N–S) is most common in the eastern part of the area and represents the dominant trend. The second group (E–W) is a subordinate fault trend, based on gravity data. These faults are far apart and are about 5 km long. The third group trends NE–SW and is mostly along the edges of the basin, particularly on the west. The fourth group (NW–SE) is common in the northwestern part of the study area. These faults are relatively closely spaced and are about 6 km long. Generally, the two fault systems that trend NE–SW and NW–SE are far less common in comparison with faults of N–S and E–W trends. The depths of faults obtained from gravity data are shallow (<2.5 to >7.5 km deep). More analytically, the NW–SE faults appear mostly to be less than 2.5 km deep, whereas the NE–SW faults vary in depth from 2.5 to 7.5 km. The depths of the N–S faults range from less than 2.5 to 5 km. The E–W faults vary in depth between 2.5 and 7.5 km. Notably, the basalt sheet covers approximately 50% of the study area, although it disappears on the gravity map, this can indicate a thin, shallow basalt sheet.

5. Discussion and conclusion

The Gharyan-Tarhuna strip was uplifted during the end of the Cretaceous or the early Tertiary. This uplift caused anticlinal swelling trending in the NW–SE direction. This movement was followed by parallel faults (NW–SE). Moreover, the study area is surrounded by two main faults and structures, the Aziziyah Fault (E–W) from the north and the Hun Graben faults (NW–SE) from the southeast. Therefore a complex structure is anticipated in the study area.

The integration of results shows that four lineament trends were recognised and mapped in the study area. The NW–SE trend probably represents the longitudinal lineaments associated with the NW–SE uplift that affected the area. The NE–SW trend was probably rejuvenated from the NNE–SSW lineaments associated with the Tibesti-Sirte Uplift, which is well developed in the Paleozoic rocks. Another idea, put forth by Klitzsch (1971), was that this trend could be related to the third period of the structural development of North Africa. The E–W trend is mostly related to the Al Aziziyah and Coastal faults. Two ideas have been advanced as to the genesis of the N–S faults: (1) this trend was caused by the same tectonic activity that was responsible for the Jabal Uplift (conjugate faults). (2) This trend was caused by some unknown tectonic movement. The first conclusion seems more logical.

Lithologically, the arrangement of phonolite hills and basalt cones in parallel lines with the dominant fault trends (NW–SE) indicates that the volcanic activity is related to the tectonic activity of the Jabal Uplift (Figure 12).

The study area was probably affected by three main movements: (1) the oldest movement (Tibesti-Sirte Uplift) began in the Late Silurian Caledonian orogeny, which initially defined the limits of the Paleozoic basins of Libya (Rusk 2001). This uplift was associated with and followed by parallel faults trending NE–SW. (2) The second movement, at the end of the Cretaceous or in the early Tertiary, was the uplift that resulted in the formation of an anticlinal swelling in a NW–SE direction. This movement was followed by parallel faulting and magmatic activities represented in the study area by intrusive and extrusive volcanic rocks. (3) The third movement was in the pre-Miocene, related to the Al Aziziyah and Coastal faults. From this movement the study area underwent E–W faulting.

Notes on contributors

Nureddin Saadi is a PhD candidate at Kyushu University (Japan). He received his BS degree in geology from the Al-fateh University in 1994, and his MS degree in remote sensing and image processing from Dundee University (UK) in 2003. His current research interests include data fusion, data integration, and geological investigations.

Essam Aboud is a doctor of engineering in exploration geophysics (Graduate School of Engineering, Kyushu University, Japan). Now he is a staff at the National Research Institute of Astronomy and Geophysics (NRIAG), Helwan, Cairo, Egypt. His research interests include potential-fields (gravity and magnetic).

Hakim Saibi is a doctor of engineering in geothermics at the Department of Earth Resources Engineering, Kyushu University, Japan and currently a JSPS Post-doctoral fellow at the same university. His research interests include geothermal engineering (geothermal reservoir simulation), hydrogeology, hydrochemistry and geodesy (gravity and GPS). His particular research interest is the monitoring of gravity change at geothermal areas using relative and absolute gravimeters. Now, he is developing gravity post-processor for a computer simulation model.

Koichiro Watanabe is a professor at the Faculty of Engineering, Kyushu University, Japan. He is the head of the Laboratory of Economic Geology. His research interest is petrological monitoring system for active volcanoes. The major aim is to progress the technologies of geological exploration, including modern geochronology, isotope geology and geological modelling.

Acknowledgements

The authors are grateful to the Editor and three anonymous referees for their comments and suggestions to improve this paper.