This file is meant to demonstrate the use of gstlearn for Faulting. It is recalled that Faulting is defined in 2-D.
import numpy as np
import pandas as pd
import sys
import os
import gstlearn as gl
import gstlearn.plot as gp
import matplotlib.pyplot as plt
Setting some global variables
ndim = 2
gl.defineDefaultSpace(gl.ESpaceType.RN, ndim)
Defining the output grid
nx = [100, 100]
grid = gl.DbGrid.create(nx)
Defining the data points
nech = 30
data = gl.Db.createFromDbGrid(nech, grid)
Defining a Model which will serve in simulating the value at the data points.
range = 40
model = gl.Model.createFromParam(gl.ECov.CUBIC, range)
Simulate the information at the data points
err = gl.simtub(None, data, model, namconv = gl.NamingConvention("Data"))
Define the set of Faults.
faults = gl.Faults()
faults.addFault(gl.PolyLine2D.create(x = [5., 10., 20.], y = [5., 25., 50.1]))
faults.addFault(gl.PolyLine2D.create(x = [50., 90.], y = [60., 70.]))
faults.addFault(gl.PolyLine2D.create(x = [20., 80., 80., 20.], y = [10., 40., 10., 10.]))
Define the Moving Neighborhood (with a maximum of 'nmaxi' samples per neighborhood) an attach the Faults.
nmaxi = 10;
neighM = gl.NeighMoving.create(nmaxi = nmaxi)
bipt = gl.BiTargetCheckFaults(faults)
neighM.addBiTargetCheck(bipt)
Plot the information
ax = data.plot()
ax = faults.plot()
ax.geometry(dims=[8,8])
ax.decoration(title="Information")
We define the parameters for the calculation of an omni-directional variogram. In a first trial, we compute the variogram without taking faults into account.
varioparam = gl.VarioParam.createOmniDirection(npas=10, dpas=5)
vario = gl.Vario.create(varioparam, data)
err = vario.compute()
vario.display()
Variogram characteristics ========================= Number of variable(s) = 1 Number of direction(s) = 1 Space dimension = 2 Variance-Covariance Matrix 0.773 Direction #1 ------------ Number of lags = 10 Direction coefficients = 1.000 0.000 Direction angles (degrees) = 0.000 0.000 Tolerance on direction = 90.000 (degrees) Calculation lag = 5.000 Tolerance on distance = 50.000 (Percent of the lag value) For variable 1 Rank Npairs Distance Value 1 5.000 5.298 0.023 2 7.000 9.910 0.166 3 16.000 15.056 1.255 4 13.000 19.965 0.687 5 18.000 24.846 0.806 6 21.000 29.884 0.785 7 19.000 34.891 0.816 8 15.000 40.427 0.772 9 19.000 44.683 1.114
Represent the variogram with the number of pairs attached to each lag.
ax = vario.plot(showPairs=True)
Performing the same variogram calculation, taking the fault into account
varioparam = gl.VarioParam.createOmniDirection(npas=10, dpas=5)
varioparam.addFaults(faults)
vario = gl.Vario.create(varioparam, data)
err = vario.compute()
vario.display()
Variogram characteristics ========================= Number of variable(s) = 1 Number of direction(s) = 1 Space dimension = 2 Variance-Covariance Matrix 0.773 Direction #1 ------------ Number of lags = 10 Direction coefficients = 1.000 0.000 Direction angles (degrees) = 0.000 0.000 Tolerance on direction = 90.000 (degrees) Calculation lag = 5.000 Tolerance on distance = 50.000 (Percent of the lag value) For variable 1 Rank Npairs Distance Value 1 5.000 5.298 0.023 2 5.000 9.895 0.212 3 9.000 14.755 1.136 4 8.000 19.703 0.301 5 10.000 24.800 0.495 6 12.000 29.845 0.597 7 11.000 34.600 0.874 8 9.000 40.303 0.693 9 11.000 44.188 0.843
Represent the variogram again with the number of pairs attached to each lag.
ax = vario.plot(showPairs=True)
err = gl.kriging(data, grid, model, neighM)
Plotting the result
gp.setDefaultGeographic(dims=[8,8])
ax = grid.plot(nameRaster="*estim")
ax = data.plot()
ax = faults.plot()
ax.decoration(title="Kriging Estimate with Faults")
ax = grid.plot(nameRaster="*stdev")
ax = data.plot()
ax = faults.plot()
ax.decoration(title="St. Dev. with Faults")
grid.display()
Data Base Grid Characteristics ============================== Data Base Summary ----------------- File is organized as a regular grid Space dimension = 2 Number of Columns = 5 Total number of samples = 10000 Grid characteristics: --------------------- Origin : 0.000 0.000 Mesh : 1.000 1.000 Number : 100 100 Variables --------- Column = 0 - Name = rank - Locator = NA Column = 1 - Name = x1 - Locator = x1 Column = 2 - Name = x2 - Locator = x2 Column = 3 - Name = Kriging.Data.estim - Locator = z1 Column = 4 - Name = Kriging.Data.stdev - Locator = NA