watex.geology.stratigraphic.GeoStrataModel#

class watex.geology.stratigraphic.GeoStrataModel(area=None, beta=5, n_epochs=100, ptol=0.1, eta=0.0001, kind='linear', degree=1, z0=0.0, max_depth=700.0, step=50.0, doi='1km', tolog10=False, verbose=False, **kwargs)[source]#

Create a stratigraphic model from 2D inversion models blocks.

The Inversion model block is two dimensional array of shape (n_vertical_nodes, n_horizontal_nodes ). Can use external packages to

build blocks and provide the 2Dblock into crm parameter.

The challenge of this class is firstly to delineate with much accuracy the existing layer boundary (top and bottom) and secondly,to predict the stratigraphy log before the drilling operations at each station. Moreover, it’s a better way to select the right drilling locationand also to estimate the thickness of existing layer such as water table layer as well as to figure out the water reservoir rock in the case of groundwater exploration [1].

Note that if the model blocks is build from external softwares. User can provided model usefull details into the keyword arguments of the fit methods. For instance the model_station_locations, model_depth, model_x_nodes etc. See more details in fit method docstring.

Parameters:
  • area (str, optional) – The name of the area where the survey is done.

  • beta (int, default=5) – Value to divide into the CRM blocks to improve the computation times.

  • n_epochs (int, default=100) – Number of iterations for new resistivity model construction (NM).

  • ptol (float, default=.1) – the error value to tolerate during NM construction between the tres values given and the calculated resistivity in crm. Higher is the error, less is the accuracy.

  • eta (float, default=1e4) – It is the learning rate for gradient descent computing to reach its convergence. If kind is set to polynomial the default value should be set`1e-8`.

  • kind (str , default='linear') – Kind of model function to compute the best fit model to replace the value in crm . Can be ‘linear’ or ‘polynomial’. if polynomial is set, specify the degree.

  • degree (int, default=1) – Polynomial function degree to implement gradient descent algorithm. If kind is set to Polynomial the default degree is 3.

  • z0 (float, default=0.) – The elevation at the first stations. Note here, we consider that no topography is included and the default is the level of the sea. If the model_depth`if provided in the `fit_params, no need to specify.

  • max_depth (float, default=700) – The maximum depth for building the block. by default, the unit is in meters. If model_depth is given, no need to provide.

  • step (float, default=50) – The step between stations. Here the unit is in meters. If the model_station_locations if given through the fit_params, no need to specified.

  • doi (default='1km') – the depth of investigation of the skin-depth. By default without any unit suffix, it should be in meters.

  • tolog10 (bool, default=False) – Convert the true values of layer resistivities to log10.

  • verbose (bool, default=False) – Output messages and warn users.

nm_#

The New resistivity model (NM) matrix with the same dimension with crm model blocks.

Type:

Ndarray of ( n_depth , n_stations )

nmSites_#

The recomputed station locations of the new resistivity model (NM)

Type:

ArrayLike

crmSites_#

The recomputed stations locations of the model-calculated resistivity (CRM)

Type:

ArrayLike,

Examples

>>> import numpy as np
>>> from watex.geology.stratigraphic import GeoStrataModel
>>>
>>> # (1):  Using the CRM without any inversion specified files.
>>>
>>> # test while files
>>> np.random.seed (42)
>>> # generate a model calculated resistivity, layers and true
>>> # resistivities values
>>> crm = np.abs( np.random.randn (215 , 70 ) *1000 )
>>> tres = np.linspace (crm.min() +1  , crm.max() +1 , 12 ) #
>>> layers = ['granites', 'gneiss', 'sedim.']
>>> gs= GeoStrataModel( to_log10 =True )
>>> gs.fit(crm, tres = tres, layers =layers).buildNM(display_infos =True )
>>> print( gs.nm_.shape)
build-NM: 18B [00:01, 12.60B/s]
----------------------------------------------------------------------
                      layers [auto=automatic](13)
----------------------------------------------------------------------
   1. hard rock (auto)                2. clay (auto)
   3. conglomerate (auto)             4. sedimentary rock (auto)
   5. gravel (auto)                   6. igneous rock (auto)
   7. fresh water (auto)              8.Sedim.
   9.Gneiss                          10. *struture not found (auto)
  11. saprolite (auto)               12. duricrust (auto)
  13.Granites
----------------------------------------------------------------------
(215, 70)
>>> gs.strataModel () # plot strata model, by default kind ='NM'
>>> gs.plotStrata ('s7')  # plot strata log at station S7
Out[7]: watex.geology.stratigraphic.PseudoStratigraphic
>>>
>>> # (2): Works with occam2d inversion files if 'pycsamt' or 'mtpy'
>>> # is installed. It will call the module Geodrill from pycsamt to
>>> make occam2d 2D resistivity block for our demo.
>>> # It presumes pycsamt is installed.
>>>
>>> from pycsamt.geodrill.geocore import Geodrill
>>> path=r'data/inversfiles/inver_res/K4' # path to inversion files
>>> inversion_files = {'model_fn':'Occam2DModel',
...                   'mesh_fn': 'Occam2DMesh',
...                    "iter_fn":'ITER27.iter',
...                   'data_fn':'OccamDataFile.dat'
...                    }
>>> tres =[10, 66, 70, 180, 1000, 2000, 3000, 7000, 15000 ]
>>> layers =['river water', 'fracture zone', 'granite']
>>> inversion_files = {key:os.path.join(path, vv) for key,
                vv in inversion_files.items()}
>>> gdrill= Geodrill (**inversion_files,
                     input_resistivities=input_resistivity_values
                     )
>>> # we can collect the 'model_res' and occam2d inversion usefull
>>> # 'attributes' from  `gdrill object` and passed to the
>>> # 'GeoStrataModel' then fit_params keyword arguments method as
>>> geosObj = GeoStrataModel(ptol =0.1).fit(
                     crm = gdrill.model_res ,
                     tres=gdrill.input_resistivities,
                     layers=gdrill.input_layer_names,
                     model_x_nodes=gdrill.model_x_nodes,
                     model_stations= gdrill.station_names,
                     model_depth= gdrill.geo_depth,
                     model_station_locations= gdrill.station_locations,
                     data_fn = gdrill.data_fn ,
                     mesh_fn=gdrill.mesh_fn,
                     iter_fn= gdrill.iter_fn
                     model_fn= gdrill.model_fn)
>>> geosObj.buildNM ()
>>> zmodel = geosObj._zmodel
>>> geosobj.nm_ # New constructed resistivity 2D model block

Notes

Modules work properly with occam2d inversion files if ‘pycsamt’ or ‘mtpy’ is installed and inherits the GeoBase class which works with geological structures and properties. Furhermore, Occam2d inversion files are also acceptables for building model blocks. However the MODEM resistivity files development is still ongoing.

References

[1]

Kouadio, L. K., Liu, R., Malory, A. O., Liu, W., Liu, C., A novel approach for water reservoir mapping using controlled source audio - frequency magnetotelluric in Xingning area , Hunan Province, China. Geophys. Prospect., https://doi.org/10.1111/1365-2478.1338