ESTIMATION OF TEMPORAL CHANGES IN SOIL MOISTURE USING RESISTIVITY METHOD

The temporal variation in a soil moisture profile can be studied using resistivity sounding data acquired at different times. The layered earth model based estimation of soil moisture from apparent resistivity data is a two-step non-linear inversion. Firstly, the apparent resistivity data are inverted to derive the layer resistivity variations and thicknesses and, secondly, the moisture content is estimated from these layer resistivity variations using a calibration equation. The soil moisture–resistivity problem was studied using the one-dimensional formulation of resistivity problem. A generalized geoelectric earth model was considered to simulate the soil moisture distribution and its temporal variation in the unsaturated zone. An algorithm (RESMOS) for the interpretation of the apparent resistivity data in terms of soil moisture variations through this two-step inversion process is reported.     

Direct interpretation of magnetotelluric sounding data based on the frequency-normalized impedance function

An important aspect of any non-linear inversion method is the generation of a suitable or good initial model as this controls the rate of convergence and accuracy of the result. To overcome the problem, a numerical method is presented for direct interpretation of magnetotelluric sounding data based on the frequency-normalized impedance (FNI) function. The expressions used to calculate the parameters are developed, first for a two-layer case under the assumption that deeper layers do not contribute to the early part of the FNI curve, and they are then generalized for an n -layer situation. The parameters of the first layer are computed by using successive sample values and the final estimate is obtained from the arithmetic mean of selected values by excluding unacceptable results in the logarithmic space. The top layer is then removed using a reduction equation. The repetition of the procedure on successive branches of the FNI function gives successive layer parameters, the resistivity of the substratum being obtained at the final step, when the reduction equation becomes equal to the square root of that resistivity.   The proposed method can be used as a complementary method for iterative inversion as it creates an initial guess which is close to the optimal solution. The solution produced by the direct interpretation may also be modified by the interpreter to incorporate prior geological information before being input to iterative interpretation schemes.     

Joint parameter estimation from magnetic resonance and vertical electric soundings using a multi‐objective genetic algorithm

Magnetic resonance sounding (MRS) has increasingly become an important method in hydrogeophysics because it allows for estimations of essential hydraulic properties such as porosity and hydraulic conductivity. A resistivity model is required for magnetic resonance sounding modelling and inversion. Therefore, joint interpretation or inversion is favourable to reduce the ambiguities that arise in separate magnetic resonance sounding and vertical electrical sounding (VES) inversions. A new method is suggested for the joint inversion of magnetic resonance sounding and vertical electrical sounding data. A one‐dimensional blocky model with varying layer thicknesses is used for the subsurface discretization. Instead of conventional derivative‐based inversion schemes that are strongly dependent on initial models, a global multi‐objective optimization scheme (a genetic algorithm [GA] in this case) is preferred to examine a set of possible solutions in a predefined search space. Multi‐objective joint optimization avoids the domination of one objective over the other without applying a weighting scheme. The outcome is a group of non‐dominated optimal solutions referred to as the Pareto‐optimal set. Tests conducted using synthetic data show that the multi‐objective joint optimization approximates the joint model parameters within the experimental error level and illustrates the range of trade‐off solutions, which is useful for understanding the consistency and conflicts between two models and objectives. Overall, the Levenberg‐Marquardt inversion of field data measured during a survey on a North Sea island presents similar solutions. However, the multi‐objective genetic algorithm method presents an efficient method for exploring the search space by producing a set of non‐dominated solutions. Borehole data were used to provide a verification of the inversion outcomes and indicate that the suggested genetic algorithm method is complementary for derivative‐based inversions.     

Combined straightforward inversion of resistivity and induced polarization (time-domain) sounding data

It is proposed that the Straightforward Inversion Scheme (SIS) developed by the authors for 1D inversion of resistivity sounding and magneto-telluric sounding data can also be used in similar fashion for time-domain induced polarization sounding data. The necessary formulations based on dynamic dipole theory are presented. It is shown that by using induced polarization potential, measured at the instant when steady state current is switched off, an equation can be developed for apparent ‘chargeability–resistivity’ which is similar to the one for apparent resistivity. The two data sets of apparent resistivity and apparent chargeability–resistivity can be inverted in a combined manner, using SIS for a common uniform thickness layer earth model to estimate the respective subsurface distributions of resistivity and chargeability–resistivity. The quotient of the two profiles will give the sought after chargeability profile. A brief outline of SIS is provided for completeness. Three theoretical models are included to confirm the efficacy of SIS software by inverting only the synthetic resistivity sounding data. Then one synthetic data set based on a geological model and three field data sets (combination of resistivity and IP soundings) from diverse geological and geographical regions are included as validation of the proposal. It is hoped that the proposed scheme would complement the resistivity interpretation with special reference to shaly sand formations.     

DAMPED LEAST-SQUARES INVERSION OF TIME-DOMAIN AIRBORNE EM DATA BASED ON SINGULAR VALUE DECOMPOSITION1

Airborne electromagnetic (AEM) methods are increasingly being used as tools of geological mapping, groundwater exploration and prospecting for coal and lignite. In such applications, quantitative interpretation is commonly based on the layered-earth model. A new approach, a damped least-squares inversion with singular value decomposition, is proposed for interpretation of time-domain, towed-bird AEM data. Studies using theoretical and field AEM data indicate that inversion techniques are dependable and provide fast converging solutions. An analysis has been made of the accuracy of model parameter determination, which depends on resistivity and thickness distribution. In the common case of conductive overburden, upper-layer resistivity and thickness are usually well determined, although situations exist where their separation becomes difficult. In the case of a resistive layer overlying a conductive basement, the layer thickness is the best-determined parameter. In both cases, estimates of basement resistivity are the least reliable. Field data obtained with the Chinese-made M-l AEM system in Dongling, Anhui Province, China, were processed using the described inversion algorithm. The survey area comprised fluvial Cenozoic clays and weathered Mesozoic sediments. Inversion of AEM data resulted in accurate depth-to-bedrock sections and realistic estimates of the resistivities of overburden and bedrock which agree with the results of drilling and resistivity sounding.     

APPLICATION OF THE MULTIFREQUENCY HORIZONTAL-LOOP EM METHOD IN OVERBURDEN INVESTIGATIONS1

In recent years, geophysical methods (shallow seismic, electromagnetic, resistivity, ground penetrating radar) have been increasingly applied to overburden investigations. Their effectiveness has been found to depend significantly on local geological conditions. Compared with advanced seismic techniques, EM methods are faster and hence more cost-effective, but they have not been considered sufficiently accurate. Analysis is carried out of data obtained with the multifrequency horizontal-loop method (HLEM) in northeastern Ontario, where the overburden consists of Quaternary glacial and glaciolacustrine sediments. Surveying along 1-6 km long profiles permitted recognition of bedrock inhomogeneities and selection of sites suitable for HLEM data interpretation using the layered model. Phasor diagrams and computer inversion based on the ridge regression technique were used to interpret HLEM soundings obtained at eight frequencies. Interpreted layer resistivities and thicknesses were correlated with the results of Rotasonic drilling at 70 sites. Relatively accurate estimates of overburden thickness (within 10%) could be obtained in about 80% of the cases. Nine examples of HLEM soundings are given and discussed: three each of one-, two- and three-layer situations. An appropriate interpretation model cannot be selected simply by minimizing the rms error or by analysing the parameter resolution matrix. Frequently, the most effective way of evaluating a solution is to consider whether resistivity values determined by inversion fit any of the ranges determined by statistical analyses of sediment resistivities. A previously published study of electrical properties of Quaternary sediments indicated that resistivities of clay, till and sand are stable within a fairly large area, such as the one under investigation. While the application of HLEM methods to mapping of Quaternary sediments can be considered a success, interpretation of EM data in regions covered by glacial sediments is more difficult than in weathered terrains, where near-surface layering is more predictable. The problem of equivalence causes non-uniqueness in interpretation. Thickness equivalence, which results in poor resistivity estimates, was found to affect areas convered by sand and till. Conductance equivalence caused poor resolution of thickness and resistivity for thin clay layers (less than 10 m).     

波形反演遗传算法及其在地震测深数据解释中的应用

为了更好地利用地震测深波形数据,提出了地震体波波形反演的遗传算法. 正演使用能精确快速计算互层结构响应的广义反、透射系数理论地震图算法;反演采用遗传算法,实现了地震体波波形反演的遗传算法. 数值试验表明,该算法具有分辨壳内低速层、高低速薄互层结构和一定的抗噪能力. 青藏高原东北缘泽库、夏河、临洮3炮地震测深P波波形反演,得到了上地壳底部低速层和中、下地壳,以及上地幔顶部薄互层的细结构图象.      

A semi-empirical velocity-porosity-clay model for petrophysical interpretation of P- and S-velocities

We design a velocity–porosity model for sand-shale environments with the emphasis on its application to petrophysical interpretation of compressional and shear velocities. In order to achieve this objective, we extend the velocity–porosity model proposed by Krief et al., to account for the effect of clay content in sandstones, using the published laboratory experiments on rocks and well log data in a wide range of porosities and clay contents. The model of Krief et al. works well for clean compacted rocks. It assumes that compressional and shear velocities in a porous fluid-saturated rock obey Gassmann formulae with the Biot compliance coefficient. In order to use this model for clay-rich rocks, we assume that the bulk and shear moduli of the grain material, and the dependence of the compliance on porosity, are functions of the clay content. Statistical analysis of published laboratory data shows that the moduli of the matrix grain material are best defined by low Hashin–Shtrikman bounds. The parameters of the model include the bulk and shear moduli of the sand and clay mineral components as well as coefficients which define the dependence of the bulk and shear compliance on porosity and clay content. The constants of the model are determined by a multivariate non-linear regression fit for P- and S-velocities as functions of porosity and clay content using the data acquired in the area of interest. In order to demonstrate the potential application of the proposed model to petrophysical interpretation, we design an inversion procedure, which allows us to estimate porosity, saturation and/or clay content from compressional and shear velocities. Testing of the model on laboratory data and a set of well logs from Carnarvon Basin, Australia, shows good agreement between predictions and measurements. This simple velocity-porosity-clay semi-empirical model could be used for more reliable petrophysical interpretation of compressional and shear velocities obtained from well logs or surface seismic data.     

REMARQUES SUR L'INTERPRETATION DES SONDAGES DE POLARISATION PROVOQUEE ESSAI DE THEORIE QUANTITATIVE *

The transient phenomena associated with induced polarization are very complex and it is necessary to have both a proper understanding of the physical processes and a good supply of mathematical aids in order to make any advance in the quantitative interpretation of I.P. surveys. The analogy between the classical electrical sounding method and the induced polarization sounding method has led a number of workers to interpret the results of the latter like the results of the former. To clarify this idea, it seems very interesting to calculate a mathematical model in order to know the form of the anomaly that can be probably expected. The case of a thin horizontal polarizable layer, with the same uniform and isotropic conductibility as the homogeneous background, has been studied in this paper. The result of the calculation for the case of a single current electrode and a single potential electrode is given. From these data, the theoretical curve of IP sounding using the Wenner electrode array is computed. The master curve is compared with experimental curves obtained on scale models. A good correlation between these curves is noticed. In addition, it can be concluded that the mathematical model is proving the validity of the experiments on interpretation scale models.     

Joint inversion of Rayleigh‐wave dispersion data and vertical electric sounding data: synthetic tests on characteristic sub‐surface models

In the traditional inversion of the Rayleigh dispersion curve, layer thickness, which is the second most sensitive parameter of modelling the Rayleigh dispersion curve, is usually assumed as correct and is used as fixed a priori information. Because the knowledge of the layer thickness is typically not precise, the use of such a priori information may result in the traditional Rayleigh dispersion curve inversions getting trapped in some local minima and may show results that are far from the real solution. In this study, we try to avoid this issue by using a joint inversion of the Rayleigh dispersion curve data with vertical electric sounding data, where we use the common‐layer thickness to couple the two methods. The key idea of the proposed joint inversion scheme is to combine methods in one joint Jacobian matrix and to invert for layer S‐wave velocity, resistivity, and layer thickness as an additional parameter, in contrast with a traditional Rayleigh dispersion curve inversion. The proposed joint inversion approach is tested with noise‐free and Gaussian noise data on six characteristic, synthetic sub‐surface models: a model with a typical dispersion; a low‐velocity, half‐space model; a model with particularly stiff and soft layers, respectively; and a model reproduced from the stiff and soft layers for different layer‐resistivity propagation. In the joint inversion process, the non‐linear damped least squares method is used together with the singular value decomposition approach to find a proper damping value for each iteration. The proposed joint inversion scheme tests many damping values, and it chooses the one that best approximates the observed data in the current iteration. The quality of the joint inversion is checked with the relative distance measure. In addition, a sensitivity analysis is performed for the typical dispersive sub‐surface model to illustrate the benefits of the proposed joint scheme. The results of synthetic models revealed that the combination of the Rayleigh dispersion curve and vertical electric sounding methods in a joint scheme allows to provide reliable sub‐surface models even in complex and challenging situations and without using any a priori information.     

Mapping of quick clay by electrical resistivity tomography under structural constraint

Geotechnical projects usually rely on traditional sounding and drilling investigations. Drilling only provides point information and the geology needs to be interpolated between these points. Near surface geophysical methods can provide information to fill those gaps. Norwegian case studies are presented to illustrate how two-dimensional electrical resistivity tomography (ERT) can be used to accurately map the extent of quick clay deposits. Quick clay may be described as highly sensitive marine clay that changes from a relatively stiff condition to a liquid mass when disturbed. Quick clay slides present a geo-hazard and therefore layers of sensitive clay need to be mapped in detail. They are usually characterized by higher resistivity than non-sensitive clay and ERT is therefore a suitable approach to identify their occurrence. However, our experience shows that ERT cannot resolve this small resistivity contrast near large anomalies such as a bedrock interface. For this reason, a constrained inversion of ERT data was applied to delineate quick clay extent both vertically and laterally. As compared to the conventional unconstrained inversions, the constrained inversion models exhibit sharper resistivity contrasts and their resistivity values agree better with in situ measurements.     

随县-西安剖面地壳结构的初步研究

本文对随县-西安测线的地震测深资料进行了初步分析,讨论了该剖面地壳内主要界面的震相特征。通过对观测走时的分层反演及综合反演,得到三层的地壳模型。地壳的平均速度为6.33km/s,厚度约34km,上地幔顶部的速度为8.10km/s。下地壳是由高速及低速层交替构成,高速层的速度为7.39km/s。构制了主要界面的剖面起伏图,并进行了讨论。     

Difficulties in determining electrical anisotropy in subsurface investigations

Surface electrical and electromagnetic methods have a limited resolution capability for determining the conductivity structure of the earth. In one-dimensional modelling a collection of many thin layers is frequently considered as one composite layer, which is then macro-anisotropic. Neither galvanic methods nor inductive methods alone can resolve the anisotropy of the ground, but a joint inversion of galvanic and inductive data may do so. The necessity of including the coefficient of anisotropy in the joint inversion of galvanic and inductive sounding data is demonstrated. An analysis is made of the combined use of geoelectrical and transient soundings to resolve the coefficient of anisotropy of a subsurface layer for varying thickness, resistivity and coefficient of anisotropy. It is found that the coefficient of anisotropy is well resolved only for layers that are many times thicker than the overburden and for coefficients of anisotropy that are not too small. The ability of the joint inversion of geoelectrical and transient sounding data to resolve macro-anisotropic layers is tested using realistic earth models determined from electrical logs.     

Obtaining basic simulation data for a heterogeneous field with stratified marine soils

Thinly stratified sedimentary deposits in a heterogeneous field were investigated to obtain basic physical data for the simulation of water flow. A procedure is described which translates a thinly stratified soil profile into a number of functional layers using functional hydrological properties. A functional layer is defined as a combination of one or more soil horizons and should (i) be recognizable during a soil survey using an auger and (ii) show significantly different functional hydrological properties when compared with another functional layer. This procedure gave three easily recognizable functional layers. Sets of hydrological characteristics of these three functional layers were obtained by physical measurements of the soil and by estimation, using textural data for classification into a standard Dutch series. The performance of several combinations of these sets was tested by comparing simulated and measured soil matric potentials for seven plots during one year. The best simulation results were obtained if measured soil hydraulic characteristics were used for relatively homogeneous functional layers and if the soil hydraulic characteristics were estimated at each location for the most heterogeneous layer.     

Use of VFSA for resolution, sensitivity and uncertainty analysis in 1D DC resistivity and IP inversion

We present results from the resolution and sensitivity analysis of 1D DC resistivity and IP sounding data using a non-linear inversion. The inversion scheme uses a theoretically correct Metropolis–Gibbs' sampling technique and an approximate method using numerous models sampled by a global optimization algorithm called very fast simulated annealing (VFSA). VFSA has recently been found to be computationally efficient in several geophysical parameter estimation problems. Unlike conventional simulated annealing (SA), in VFSA the perturbations are generated from the model parameters according to a Cauchy-like distribution whose shape changes with each iteration. This results in an algorithm that converges much faster than a standard SA. In the course of finding the optimal solution, VFSA samples several models from the search space. All these models can be used to obtain estimates of uncertainty in the derived solution. This method makes no assumptions about the shape of an a posteriori probability density function in the model space. Here, we carry out a VFSA-based sensitivity analysis with several synthetic and field sounding data sets for resistivity and IP. The resolution capability of the VFSA algorithm as seen from the sensitivity analysis is satisfactory. The interpretation of VES and IP sounding data by VFSA, incorporating resolution, sensitivity and uncertainty of layer parameters, would generally be more useful than the conventional best-fit techniques.     

Surface‐wave inversion for a P‐velocity profile with a constant depth gradient of the squared slowness

Surface waves are often used to estimate a near‐surface shear‐velocity profile. The inverse problem is solved for the locally one‐dimensional problem of a set of homogeneous horizontal elastic layers. The result is a set of shear velocities, one for each layer. To obtain a P‐wave velocity profile, the P‐guided waves should be included in the inversion scheme. As an alternative to a multi‐layered model, we consider a simple smooth acoustic constant‐density velocity model, which has a negative constant vertical depth gradient of the squared P‐wave slowness and is bounded by a free surface at the top and a homogeneous half‐space at the bottom. The exact solution involves Airy functions and provides an analytical expression for the dispersion equation. If the ratio is sufficiently small, the dispersion curves can be picked from the seismic data and inverted for the continuous P‐wave velocity profile. The potential advantages of our model are its low computational cost and the fact that the result can serve as a smooth starting model for full‐waveform inversion. For the latter, a smooth initial model is often preferred over a rough one. We test the inversion approach on synthetic elastic data computed for a single‐layer P‐wave model and on field data, both with a small ratio. We find that a single‐layer model can recover either the shallow or deeper part of the profile but not both, when compared with the result of a multi‐layer inversion that we use as a reference. An extension of our analytic model to two layers above a homogeneous half‐space, each with a constant vertical gradient of the squared P‐wave slowness and connected in a continuous manner, improves the fit of the picked dispersion curves. The resulting profile resembles a smooth approximation of the multi‐layered one but contains, of course, less detail. As it turns out, our method does not degrade as gracefully as, for instance, diving‐wave tomography, and we can only hope to fit a subset of the dispersion curves. Therefore, the applicability of the method is limited to cases where the ratio is small and the profile is sufficiently simple. A further extension of the two‐layer model to more layers, each with a constant depth gradient of the squared slowness, might improve the fit of the modal structure but at an increased cost.     

A FAST METHOD FOR DETERMINING THE LAYER DISTRIBUTION FROM THE RAISED KERNEL FUNCTION IN GEOELEGTRICAL SOUNDING*

In a previous publication (Koefoed 1968) a function called the “raised kernel function” has been introduced as an intermediate function in the interpretation of resistivity sounding data, and methods have been described both for the determination of the raised kernel function from the apparent resistivity function, and for the determination of the layer distribution from the raised kernel function. In the present paper a procedure is described by which the second step in this interpretation method–i.e. the determination of the layer distribution from the raised kernel function–is considerably accelerated. This gain in interpretation speed is attained by the use of a standard graph for a function which defines the reduction of the raised kernel function to a lower boundary plane.     

STUDYING THE INTERCONNECTION OF GROUND AND SURFACE WATERS IN REPRESENTATIVE BASINS UNDER THE PROGRAMME OF THE INTERNATIONAL HYDROLOGICAL DECADE

Abstract

Four geoelectrical soundings were measured with a combination of Schlumberger and azimuthal or equatorial dipole electrode arrays on a Carboniferous limestone basin of the Condroz area, Belgium. The measuring technique is briefly outlined as well as the interpretation procedure, which follows a closed-loop scheme with control of calculated model curves. Some general problems of interpretation of geoelectrical sounding curves are tackled, as far as they have a practical bearing on the treatment of Condroz soundings.

The problem of determining the very high resistivity of limestone is approached through ARCHIE's formula, an empirical relation between the bulk rock resistivity, the porosity and the electrolyte resistivity. An evaluation of the latter two parameters, combined with electrical horizontal conductance measurements directly made on resistivity sounding curves, offers a possibility for fast determination of the total water storage in a limestone aquifer. Such storage determinations could be applied whenever an aquifer shows up as a conductive layer interbedded between two highly resistant layers (e.g. nonsaturated limestone and compact, non-fractured limestone).     

海洋CSEM的共中心点域快速反演

海洋可控源电磁法(CSEM)对海底高阻体的反映比较灵敏,可用于天然气水合物探测资料的定性解释和反演研究.海洋CSEM资料的共中心点(CMP)域转换方式,可在横向上较好地分辨高阻储层.本文提出在CMP域实现一维频率域海洋可控源电磁资料高斯-牛顿反演算法.鉴于一维反演是解释地球物理资料的基础,较于二维和三维反演方法有着更高的计算效率和更低的硬件要求,将二维模型的响应在CMP域单元内表达为一维模型的响应,进而运用一维高斯-牛顿反演解释二维海洋CSEM资料.模型数据试算表明,海洋CSEM的CMP域反演速度较快,能够实现二维CSEM资料的反演解释.     

倾子资料三维共轭梯度反演研究

在对倾子响应和共轭梯度算法深入分析的基础上,我们实现了倾子资料三维共轭梯度反演算法.基于倾子资料的三维共轭梯度反演研究,探讨了利用倾子资料进行三维反演定量解释的方法.通过对理论模型合成数据 进行反演试算,验证了所实现的倾子资料三维共轭梯度反演算法的有效性和稳定性.该反演算法可用于对大地电磁测深和地磁测深(地震地磁台站进...