# Merge-Optimization Method of Combined Tomography of Seismic Refraction and Resistivity Data

Authors: Andy A. Bery
Journal: ISRN Geophysics
 (2012)

---

## Abstract

This paper discussed a novel application called merge-optimization method that combines resistivity and seismic refraction data to provide a detailed knowledge of the studied site. This method is interesting because it is able to show strong accuracy of two geophysical imaging methods based on many of data points collected from the conducted geophysical surveys of disparate data sets based strictly on geophysical models as an aid for model integration for two-dimensional environments. The geophysical methods used are high resolution methods. The resistivity imaging used in this survey is able to resolve the subsurface condition of the studied site with low RMS error (less than 2.0%) and 0.5 metre electrodes interval. For seismic refraction method, high resolution of seismic is used for correlation with resistivity results. Geophones spacing is 1.0 metre and the total number of shot-points is 15, which provides very dense data point. The algorithms of merge-optimization have been applied to two data sets collected at the studied site. The resulting images have been proven to be successful because they satisfy the data and are geometrically similar. The regression coefficient found for conductivity-resistivity correlation is 95.2%.


---

## Body

## 1. Introduction

The characterization of the subsurface requires a detailed knowledge of several properties of the composing rocks and fluids. Whereas some of these properties can be measured directly (seismic and borehole methods), other properties have to be estimated by indirect measurement methods such as resistivity, TEM, and magnetic. However, it is not uncommon that the geophysical data yield models of limited accuracy which may not contribute significantly to our understanding of the subsurface condition or may show incompatibilities. Thus, a new technique needs to be produced not only for better interpretation by geophysics but also for nongeophysical background people such as engineers and architects. The distribution of uncorrelated physical properties seems to be controlled by common subsurface attributes, when taken into account, able to improve and resolve the accuracy of the geophysical imaging results. An outstanding feature of the subsurface that is common to the geophysical data is the geometrical distribution of the physical properties which can be measured by the physical property changes. This condition of commonality can be incorporated in the process of estimation to obtain meaningful and more reliable subsurface imaging results.


## 2. Methodology for Merge-Optimization Method

In this paper, seismic refraction data using the SeisOpt@2D software were developed by Pullammanappallil and Louie [1], to obtain velocity tomograms. The tomography inversion is performed by a generalized simulated-annealing method of optimization including a controlled Monte-Carlo inversion. The first arrival times and survey geometry are used as input and no initial velocity model is needed. The apparent resistivity data were inverted using the software package RES2DINV [2]. This program solves the tomography inversion problem using a smoothness-constrained least squares method. The logarithmic average of the measured apparent resistivity values is input as the homogeneous starting model [2].This paper adopts a merge problem formulation with resistivity-velocity cross gradients function in order to provide the required effective link between the resistivity model and the seismic velocity model. The cross product of the gradient is defined as(1)t→(x,y,z)=∇mr(x,y,z)×∇ms(x,y,z),
where mr and ms refer to the logarithm of resistivity and P-wave slowness, respectively. In the two-dimensional case, t→(x,y,z) always points in the strike direction, that is, it will be treated as a scalar t. It is incorporated as part of the objective function:
(2)Φ(m→r,m→s)=∥d→r-f→r(m→r)d→s-f→s(m→s)∥Cdd-12+∥αrD→m→rαsD→m→s∥2+∥m→r-m→Rrm→s-m→Rs∥CRR-12Subjecttot→(m→)=0.
Here, d→ represents the vector of observed data (logarithm of apparent resistivity, d→r and seismic travel-times, d→s), m→=[m→r∶m→s]T is the vector of the model parameters, f→ is the theoretical model respond, D→ is the discrete version of the smoothing operator acting on m→,m→R=[m→Rr∶m→Rs]T is an a priori model, C→dd  is the covariance of the field data (assumed diagonal, i.e., fully uncorrelated data), C→RR is the covariance of the priori model (assumed diagonal), and αr and αs are weighting factors that control the level of smoothing of the resistivity and seismic imaging models.The cross gradients criterion requires the problem to satisfy the conditiont→(m→)=0, where any spatial changes occurring in both resistivity and velocity must point to the same point in the same or opposite direction irrespective of amplitude. In a geological sense, this implies that if a boundary exists, then it must be sensed by both geophysical methods in a common orientation regardless of the amplitude of the physical property changes. An additional flexibility of the technique is that the cross gradients constraint is also satisfied when either ∇m→r or ∇ms vanishes in some part of the model, thus giving the models the possibility of admitting a geological boundary which has a significant change only in the electrical resistivity or seismic velocity of the adjoining rocks.In iteration 2D optimization approach, the subsurface model is discretised into rectangular cells of variable sizes optimized according to the natural sensitivity of each particular set of resistivity and seismic velocity measurements. We define the discrete version of (1) and the corresponding derivatives using the elements of the 3 cells scheme depicted in Figure 1. For the along strike components this simplifies to
(3)t≅4ΔxΔz[mrc(msb-msr)+mrr(msc-msb)4ΔxΔz+mrb(msr-msc)],
where the quantities Δx,  Δz,  mrc,  mrb, and so forth are defined in Figure 1 above.Figure 1
Definition of the resistivity-velocity cross gradients function and its derivatives on a rectangular grid domain. For 2D grid extending in thex and z directions, with each grid element characterised by logarithm of resistivity mr and seismic slowness ms, the function t is defined at the centre of a given element (marked with x) considering the parameters from two elements it is in contact with its right (subscripted r) and bottom (subscripted b).Here, the first order of Taylor series expression is used, and (2) is equivalent to
(4)min{2n→2Tm→+m→TN1m→}subjecttot(m→0)+B→m→-B→m→0=0,
where
(5)N1=[A→rTC→ddr-1A→r+αr2D→TD→+CRRr-100A→sTC→dds-1A→s+αs2D→TD→+C→RRs-1],n2=[A→rTC→ddr-1{d→r-f→r(m→0r)+A→rm→0r}+C→RRr-1m→RrA→sTC→dds-1{d→s-f→s(m→0s)+A→sm→0s}+C→RRs-1m→Rs];A→r,  A→s, and B→ are the respective partial derivatives of f→r,  f→s, and t→ evaluated at the initial model,  m→0=[m→0r∶m→0s]T. The Jacobian matrix for seismic data, A→s is computed using ray tracing as suggested by Vidale [3] and Zelt and Barton [4]. The solution to (4) used in our iterative scheme is given by
(6)m→=N→1-1n→2-N→1-1B→T(B→N→1-1BT)-1[B→N→1-1n→2-B→m→0m→=N→1-1n→2-N→1-1B→T(B→N→1-1BT)-1+t→(m→0)].
In the regularised solution process, the weighting factors are initially assigned large values which are then gradually reduced in subsequent iterations until the data are fitted to the required level. The merge optimization is initiated using a half-space model in the absence of reliable a priori information. The implemented cross gradients criterion and regularisation measures ensure that the resolution characteristics of the individual data sets are fully exploited the search for structurally linked models. Note that with the cross gradients criterion, there is no need to define or assume ab initio any interdependence resistivity and seismic velocity which could bias the inverse solution.


## 3. Seismic and Resistivity Characterization by Merge-Optimization Method

In this paper, we used data collected during recent electrical tomography and seismic refraction surveys [5]. The same data sets were used for simultaneous subsurface imaging using algorithm for improved near-surface characterization. We expect that the merge-optimization method of the data will define more accurately the main geological and features of the subsurface. Figure 2 shows the geometry used for both infield tests. The resistivity data are collected using a Wenner-Schlumberger array using with 0.5 metre electrode intervals. As for seismic refraction, 15 shot points and 1.0 metre geophone intervals are used. The purpose of this geometry is to increase the resolution of the imaging results (less noise effect and moderate sensitivity to both vertical and horizontal structures) and hence both separate data sets are more to correlate and increase the reliability of the results. The study site is located in a slope area, thus it is more challenging than the flat ground.Figure 2
Geometry of the infield tests for both geophysical methods.The seismic refraction velocity and electrical resistivity imaging derived by merge-optimization of infield tests data sets (Figures3(a) and 4(a)) is less than 2.0% of root mean-squares (RMS) for electrical resistivity imaging results. The reconstructed distributions of the model parameters show structural similarities and hence good spatial correlation of velocity with resistivity data sets. Note that the cross gradients criterion serves for geological structural control but does not force the two models into conformity. Figures 3(b) and 4(b) show the seismic refraction results for the same area and the survey geometry. Figures 3(c) and 4(c) show the arcs formation (seismic refraction) at the subsurface in the studied site.Figure 3
March infield’s test results for optimal 2D merge resistivity and velocity models (a) resistivity model, (b) velocity model, and depth section with arcs.Figure 4
April infield’s test results for optimal 2D merge resistivity and velocity models (a) resistivity model, (b) velocity model, and depth section with arcs.The resistivity images show that the subsurface of the studied sit consists of two main zones. Resistivity values lower than 900 Ω.m are indicative of residual (clayey sand soils) while values higher than 1,100 Ω.m are indicative of a weathered layer. The presence of moist zones and dry zones can be associated with loose zones and compacted soil, respectively. The seismic refraction images showed that the subsurface consists of three layers. Velocity values of 370–500 m/s are associated with loose soil mixed with boulders (high resistivity value near surface). The second layer has velocity values of 600–800 m/s associated to hard layer (unsaturated) and the third layer has velocity values greater than 1,000 m/s which are associated with a saturated layer. Comparing the results of the two geophysical methods, we can summarize that the resistivity method has the limitation of a lower depth of investigation with respect to seismic refraction. On the other hand, seismic refraction is unable to resolve well the subsurface features and has much less resolution compare to electrical resistivity. However, both data sets showed their validity and reliability when correlated together to determine the subsurface features of studied site. In this paper, 838 of data points correlated. This can be showed by relationship resistivity-velocity (Figure5).Figure 5
The results show an L-shaped (A-B) trend that may be a consequence of a water table or natural divide between consolidated and unconsolidated materials.To ascertain the resistivity-velocity relationship for the reconstructed models, we have plotted in Figure5 these two physical parameters for all the coincident sampling positions (838 data points) in the merge-optimization models. This plot shows a structural feature that did not come out in the results of optimization. The merge-optimization results suggest distinct subgroupings in the near-surface volume. From the correlated data points, we found an L-shaped (A-B) trend that may be a consequence of a water table or natural divide between consolidated and unconsolidated materials.In Figure6, notice that the seismic rays are critically refracted at the top of the structural units mapped as layer 2, defining the possible boundary upper dry layer and moist layer where water table exist through it effect of water seepage from the top of slope. Merge-optimization thus appears to offer a tool for improved characterization of complex near-surface environments. To address the question: if the contention that electrical resistivity and seismic refraction imaging are both influenced by the same geological structures is correct, and there is some structural control on the distribution of physical properties in the subsurface, can accurately deduce any physical feature from merge structural imaging? It is possible that Figures 5 and 6 have physical implications from the external factors. Thus, before this merge-optimization method was applied, we have did a study on the physical characterizations of the soils (moisture content W, cohesion C′, friction angle Ø′, porosity n, void ratio e, and saturation degree S) at various locations along the survey lines at different periods [5, 6].(a) Interpretative model showing empirical correlation between conductivity (S/m) and resistivity (Ω.m) is found asσ=0.0037e-20.001(ρ) and strong relationship of data sets regression coefficient which is 95.2%. (b) Near-surface structure of the studied site can be recognized as a two-layer case where the upper layer is unconsolidated (weathered surface materials) and the lower layer is made of consolidated materials.
(a)
(b)


## 4. Conclusion

The incorporation of the cross-gradients criterion in 2D optimization leads to a geologically meaningful solution by improving the near-surface structural conformity between the velocity and resistivity imaging, without forcing or assuming the form of the relationship between two geophysical methods which have their own advantage and limitation. The cross-gradients criterion also allows detecting subsurface features to which only one of the geophysical techniques is sensitive, leading to a better structural characterization. The application of these geophysical techniques combination of 2D optimization with cross-gradients to the data collected from seismic refraction and electrical resistivity field surveys has led to an improved characterization of the near-surface material and features. This study suggests that unconsolidated (possibly unsaturated and saturated) materials may be subclassified on the basis of the resistivity-velocity relationship revealed from the application of merge-optimization method. The cross-gradients approach adopted in this paper can also be used for 3D problems and for any combination of independent geophysical methods.



---

*Source: 293132-2012-12-31.xml*