base module

This module is the core of the ArchPy and contains most of the central classes and functions.

class base.Arch_table(name, working_directory='ArchPy_workspace', seed=197540, write_results=False, fill_flag=False, verbose=1, ncpu=-1)

Bases: object

Major class of ArchPy. Arch_table is the central object that can be assimilated as a “project”. Practically every operations are done using an Arch_table

Parameters:

  • name (str) – Name of the project

  • working_directory (str) – Path to the working directory

  • seed (int) – Seed for random number generation

  • write_results (bool) – If True, results will be written to disk

  • fill_flag (bool) – If True, the top unit will be filled with the most common facies

  • verbose (int) – Verbosity level

  • ncpu (int) – Number of cpus to use for parallel operations, if -1, all cpus will be used

add_bh(bhs)

Method to add borehole, list of boreholes if multiples

Parameters:

bhs (borehole object or list of borehole objects) – borehole(s) to add to the Arch_Table

add_fake_bh(bhs)

Method to add a fake borehole, list if multiples Use for inversion purposes

Parameters:

bhs (list of borehole or borehole) – boreholes to add

add_geological_map(raster)

Add a geological map to Arch_table

Parameters:

raster (2D ndarray of size (ny, nx)) – geological map to add. Values are units IDs

add_grid(dimensions, spacing, origin, top=None, bot=None, rspl_method='nearest', polygon=None, mask=None)

Method to add/change simulation grid, regular grid.

Parameters:

  • dimensions (sequence of size 3,) – number of cells in x, y and direction (nx, ny, nz)

  • spacing (sequence of size 3,) – spacing of the cells in x, y and direction (sx, sy, sz)

  • origin (sequence of size 3,) – origin of the simulation grid of the cells in x, y and direction (ox, oy, oz)

  • top, bot (2D ndarray of dimensions (ny, nx) or float or raster file,) – top and bottom of the simulation domain

  • rspl_method (string) – scipy resampling method (nearest, linear and cubic –> nearest is generally sufficient)

  • polygon (2D ndarray of dimensions (ny, nx)) – boolean array to indicate where the simulation is active (1) or inactive (0). Polygon can also be a Shapely (Multi) - polygon.

  • mask (3D ndarray of dimensions (nz, ny, nx)) – 3D boolean array to indicate where the simulation is active (1) or inactive (0). If given, top, bot and polygon are ignored.

add_prop(prop)

Add a property to the Arch_Table

Parameters:

prop (Prop object) – property to add to the Arch_Table

Return type:

None

add_prop_hd(prop, x, v)

Add Hard data to the property “prop”

Parameters:

  • prop (string) – property name of the ~ArchPy.base.Prop object

  • x (ndarray of size (ndata, 3)) – x, y and z coordinates of hd points

  • v (array of size (ndata)) – HD property values at x position

celltoindex(cell_x, cell_y, cell_z)

cell index to cell position

Parameters:

  • cell_x (int) – cell index in x direction

  • cell_y (int) – cell index in y direction

  • cell_z (int) – cell index in z direction

Returns:

(x,y,z) cell position

Return type:

tuple

check_bh_inside(bh)

Check if a borehole is inside the simulation domain and cut boreholes if needed

Parameters:

bh (Borehole) – borehole object

Returns:

  • 0 if borehole is outside the simulation domain

  • 1 if borehole is inside the simulation domain

check_piles_name()

check the names of piles in order to be sure that they are different

Returns:

1 if all is ok

Return type:

int

check_units_ID()

check the IDs of units in order to be sure that they are different

Returns:

1 if all is ok

Return type:

int

compute_distribution()

Compute the probability distribution of the hard data for each unit For each unit, the distribution is computed using the Normal Score Transform

compute_domain(s1, s2)

Return a bool 2D array that define the domain where the units exist (between two surfaces, s1 and s2)

Parameters:

s1, s2 (2D ndarrays) – two surfaces over the simulation domain size: (ny, nx)

Returns:

3D ndarray of bools

Return type:

domain

compute_facies(nreal=1, verbose_methods=0)

Performs the computation of the facies

Parameters:

  • nreal (int) – number of realizations

  • verbose_methods (int) – verbose for the facies methods, 0 by default

compute_geol_map(iu=0, color=False)

Compute and return the geological map for given unit realization

Parameters:

  • iu (int) – unit realization index

  • color (bool) – if True return a 3D array with RGBA values

Returns:

geological map

Return type:

2D ndarray

compute_mean_prop(property, type='arithmetic', inside_units=None, inside_facies=None)
compute_prop(nreal=1)

Performs the computation of the properties added to the ArchTable

Parameters:

nreal (int) – number of realizations

compute_surf(nreal=1, fl_top=True, rm_res_files=True)

Performs the computation of the surfaces

Parameters:

  • nreal (int) – number of realizations

  • fl_top (bool) – assign first layer to top of the domain (True by default)

  • rm_res_files (bool) – flag to remove previous existing resulting files in working directory

coord2cell(x, y, z=None)

Method that returns the cell in which are the given coordinates

Parameters:

  • x (float) – x coordinate

  • y (float) – y coordinate

  • z (float, optional) – z coordinate, if z is None, only x and y indexes are returned

Returns:

cell – cell indexes (ix, iy, iz) or (ix, iy) if z is None

Return type:

tuple

cross_section(arr_to_plot, p_list, esp=None)

Return a cross section along the points pass in p_list

Parameters:

  • arr_to_plot (3D or 4D array of dimension nz, ny, nx(, 4)) – array of which we want a cross section. This array will be considered being part of the ArchPy simulation domain.

  • p_list (sequence of tuple) – list or array of tuple containing x and y coordinates (e.g. p_list=[(100, 200), (300, 200)] –> draw a cross section between these two points)

  • esp (float) – spacing to use when sampling the array along the cross section

Returns:

cross section ready to plot

Return type:

2D array

define_domains(surfaces, fl_top=True)

Performs the computation of the units domains when surfaces are provided

Parameters:

  • surfaces (dictionary of surfaces as values and pile name as key)

  • fl_top (bool) – assign first layer to top of the domain (True by default)

draw_cross_section(background='units', iu=0, ifa=0)
erase_hd()

Erase the hard data from all surfaces and facies

estimate_surf_params(default_covmodel=None, auto=False, **kwargs)

Alias for infer surface in ArchPy.infer

Parameters:

  • default_covmodel (str, default None) – default covariance model to use for surface estimation

  • auto (bool) – to automatically infer parameter (True) or not (False)

  • kwargs – various kwargs and parameters that can be passed to ArchPy.infer.infer_surface or ArchPy.infer.fit_surfaces

extract_log_facies_bh(facies, bhx, bhy, bhz, depth)

Extract the log facies in the ArchPy format at the specific location

Parameters:

  • facies (3D ndarray) – facies realization

  • bhx, bhy, bhz (float) – borehole location

  • depth (float) – depth of investigation

Returns:

log_facies – a log facies (list of facies object with elevations)

Return type:

list of Facies object

fill_ID(arr, ID=0)

Fill ID values in an 3D array given surroundings values using nearest neighbors

Parameters:

  • arr (ndarray of size (nz, ny, nx)) – simulation grid size

  • ID (ID) – ID to replace

Returns:

arr – simulation grid size

Return type:

ndarray of size (nz, ny, nx)

fill_top_unit(method='nearest_neighbors')

Function to fill each cells simulated top unit given their nearest neighbour.

Parameters:

method (str) – method to fill top unit. Default is “nearest_neighbors”

geol_contours(step=5)

This function extract information at the boundaries between units from the given geological map.(see add_geological_map()) Results are returned as a list of borehole objects

Parameters:

step (int) – step between each cells where to add a contour information

Return type:

list of borehole

get_all_facies(recompute=True)

Return a list of all facies

Parameters:

recompute (bool) – if False, the list all facies attribute will be simply retrieve. Even if changes have been made on the project.

Return type:

list of Facies objects

get_all_units(recompute=True)

Return a list of all units, even sub-units

Parameters:

recompute (bool) – if False, the list all units attribute will be simply retrieve. Even if changes have been made.

Return type:

list of Unit objects

get_bh(ID)

Return the borehole object given its ID

Parameters:

ID (string) – borehole ID

Return type:

Borehole object

get_bounds()

Return bounds of the simulation domain

Returns:

(xmin, xmax, ymin, ymax, zmin, zmax)

Return type:

tuple

get_entropy(typ='units', h_level=None, recompute=False)

Compute the Shannon entropy for units or facies and return it.

Parameters:

  • typ (string) – type of models to use to compute the entropy, valid values are “units” and “facies”

  • h_level (int) – hiearchical level to use to compute the entropy, only used if typ is “units”

  • recompute (bool) – if True, recompute the entropy

Return type:

3D ndarray of size (nz, ny, nx)

get_facies(iu=0, ifa=0, all_data=True)

Return a numpy array of 1 or all facies realization(s).

Parameters:

  • iu (int) – unit index

  • ifa (int) – facies index

  • all_data (bool) – return all the units simulations

Returns:

facies domains

Return type:

ndarray

get_facies_obj(name='A', ID=1, type='name', vb=1)

Return the facies in the Pile with the associated name or ID.

Parameters:

  • name (string) – name of the strati to retrieve

  • ID (int) – ID of the strati to retrieve

  • type (str, (name or ID)) – retrieving method

  • vb (int) – verbosity level (0 or 1)

Return type:

Facies object

get_nx()

Returns the number of cells in x direction

get_ny()

Returns the number of cells in y direction

get_nz()

Returns the number of cells in z direction

get_ox()

Returns the origin of the grid in x direction

get_oy()

Returns the origin of the grid in y direction

get_oz()

Returns the origin of the grid in z direction

get_pile_master()

Returns the Pile_master object

get_piles()

Return a list of all the subpiles

Return type:

list of Pile objects

get_prop(name, iu=None, ifa=None, ip=None, all_data=True)

Return a numpy array of 1 or all facies realization(s).

Parameters:

  • iu (int) – unit index

  • ifa (int) – facies index

  • ip (int) – property index

  • all_data (bool) – return all the units simulations

Returns:

1 or all facies realization(s)

Return type:

numpy array

get_proportions(type='units', depth_min=0, depth_max=inf, ignore_units=[], mask=None)

Function that returns the proportions of the units in the boreholes

Parameters:

  • type (str) – type of proportions to return “units” : proportions of units “facies” : proportions of facies

  • depth_max (float, default np.inf) – maximum depth of investigation in the boreholes

  • depth_min (float, default 0) – minimum depth of investigation in the boreholes

  • ignore_units (list of str, default []) – units name to ignore during the analysis

  • mask (2D ndarray of size (ny, nx), default None) – mask where to analyse the boreholes

Return type:

dictionnary

get_surface(h_level='all')

Return a 4D array of multiple surfaces according to the hierarchical level desired, by default ArchPy try to return the highest hierarchical unit’s surface

Parameters:

h_level (int or string,) – maximum level of hierarchy desired to return, “all” indicates to return the highest hierarchical level possible for each unit

Returns:

  • 4D ndarray of size (nlayer, nsim, ny, nx) – All computed surfaces in a nd.array of size (nlayer, nsim, ny, nx)

  • list – a list (nlayer) of units name corresponding to the surfaces to distinguish wich surface correspond to which unit

get_surfaces_unit(unit, typ='top')

Return a 3D array of computed surfaces for a specific unit

Parameters:

  • unit (Unit object) – a Unit object contained inside the master pile or another subunit

  • typ (string, (top, bot or original),) – specify which type of surface to return, top, bot or original (surfaces before applying erosion and stratigrapic rules)

Returns:

All computed surfaces in a nd.array of size (nreal_units, ny, nx)

Return type:

3D ndarray of size (nreal_units, ny, nx)

get_sx()

Returns the size of the cells in x direction

get_sy()

Returns the size of the cells in y direction

get_sz()

Returns the size of the cells in z direction

get_unit(name='A', ID=1, type='name', all_strats=True, vb=1)

Return the unit in the Pile with the associated name or ID.

Parameters:

  • name (string) – name of the strati to retrieve

  • ID (int) – ID of the strati to retrieve

  • type (str, (name or ID)) – retrieving method

  • all_strats (bool) – flag to indicate to also search in sub-units, if false research will be restricted to units directly in the Pile master

  • vb (int) – verbosity level

Return type:

Unit object

get_units_domains_realizations(iu=None, all_data=True, fill='ID', h_level='all')

Return a numpy array of 1 or all units realization(s).

iu: int

simulation to return

all_data: bool

return all the units simulations, in that case, iu is ignored

fill: string

ID or color are possible, to return realizations with unit ID or RGBA color (for plotting purpose e.g. with plt.imshow)

h_level: string or int

hierarchical level to plot. A value of 1 indicates that only unit of the master pile will be plotted. “all” to plot all possible units

get_xg()

Returns the edges of the grid in x direction

get_xgc()

Returns the centers of the grid in x direction

get_yg()

Returns the edges of the grid in y direction

get_ygc()

Returns the centers of the grid in y direction

get_zg()

Returns the edges of the grid in z direction

get_zgc()

Returns the centers of the grid in z direction

getbhindex(ID)

Return the index corresponding to a certain borehole ID

Parameters:

ID (string) – borehole ID

Returns:

index of the borehole in the pile

Return type:

int

getpropindex(name)

Return the index corresponding to a certain prop name

Parameters:

name (str) – name of the property to retrieve

Returns:

index of the property

Return type:

int

hd_fa_in_unit(unit, iu=0)

Extract facies hard data for a unit and send warning if a hard data should not be in the unit

Parameters:

  • unit (unit object) – unit to extract hard data

  • iu (int) – unit realization index to extract hard data

hd_un_in_unit(unit, iu=0)

Extract sub-units hard data for a unit

Parameters:

  • unit (unit object) – unit to extract hard data

  • iu (int) – index of the unit realization to extract hard data

Returns:

  • hd (list of tuples) – list of hard data coordinates (x,y,z)

  • sub_units (list of int) – list of sub-units ID

hierarchy_relations(vb=1)

Method that sets the hierarchical relations between units and sub-units

make_fake_bh(positions_x, positions_y, units=None, facies=None, stratIndex=0, faciesIndex=0, extractUnits=True, extractFacies=True, vb=1)

Create fake boreholes from realization of the Arch_table.

Parameters:

  • positions_x (sequence of numbers) – indicate the x positions of the borehole to create bhs

  • positions_y (sequence of numbers) – indicate the y positions of the borehole to create bhs

  • units (3D array, optional) – unit array to use to create bhs

  • facies (3D array, optional) – facies array to use to create bhs

  • stratIndex (int or sequence of int, optional) – unit index to sample

  • faciesIndex (int or sequence of int, optional) – facies index to sample

  • extractUnits (bool, optional) – flag to indicate to sample units or not

  • extractFacies (bool, optional) – flag to indicate to sample facies or not

  • vb (int, optional) – verbose level

Return type:

list of borehole objects

order_Piles()

Order all the units in all the piles according to order attribute

orientation_map(unit, azi_top='gradient', dip_top='gradient', azi_bot='gradient', dip_bot='gradient', iu=0, smooth=2)

Compute orientation maps for a given unit (azimuth and dip)

Parameters:

  • unit (Unit object) – unit object from which we want to have the orientation map

  • azi_top (string or float) – method to use to infer azimuth of the top of the unit

  • azi_bot (string or float) – method to use to infer azi of the bottom of the unit The angles are then interpolated between top and bottom (linear interpolation)

  • dip_top (string or float) – same as azimuth but for dip

  • dip_bot (string or float) – same as azimuth but for dip

  • iu (int) – unit realization index (0, 1, …, Nu)

  • smooth (int) – half-size windows for rolling mean

Returns:

orientation map

Return type:

2D ndarray

physicsforward(method, positions, stratIndex=None, faciesIndex=None, propIndex=None, idx=0, cpuID=0)
plot_arr(arr, var_name='V0', v_ex=1, plotter=None, slicex=None, slicey=None, slicez=None, filtering_interval=None, cmin=None, cmax=None, scalar_bar_kwargs=None, **kwargs)

This function plot a 3D array with the same size of the simulation domain

Parameters:

  • arr (3D array) – array to plot. Size of the simulation domain (nx, ny, nz)

  • var_name (str) – variable names to plot and to show in the pyvista plot

  • v_ex (float) – vertical exaggeration

  • plotter (pyvista.Plotter) – pyvista plotter to plot on

  • slicex (float or sequence of floats) – fraction of the x axis to plot

  • slicey (float or sequence of floats) – fraction of the y axis to plot

  • slicez (float or sequence of floats) – fraction of the z axis to plot

  • filtering_interval (list of two floats) – interval of values to keep

  • cmin (float) – minimum value for the colorbar

  • cmax (float) – maximum value for the colorbar

  • scalar_bar_kwargs (dict) – dictionary of arguments to pass to the pyvista scalar bar

plot_bhs(log='strati', plotter=None, v_ex=1, plot_top=False, plot_bot=False)

Plot the boreholes of the Arch_table project.

Parameters:

  • log (string) – which log to plot –> strati or facies

  • plotter (pyvista plotter)

  • v_ex (float) – vertical exaggeration

  • plot_top (bool) – if the top of the simulation domain must be plotted

  • plot_bot (bool) – if the bot of the simulation domain must be plotted

plot_cross_section(p_list, typ='units', arr=None, iu=0, ifa=0, ip=0, property=None, esp=None, ax=None, colorbar=False, ratio_aspect=2, i=0, dist_max=100, width=0.5, vmax=None, vmin=None)

Plot a cross section along the points given in p_list with a spacing defined (esp)

Parameters:

  • p_list (list or array of tuple containing) – x and y coordinates (e.g. p_list=[(100, 200), (300, 200)] –> draw a cross section between these two points)

  • typ (string) – units, facies or a property name

  • iu (int) – units index realization

  • ifa (int) – facies index realization

  • ip (int) – property index realization

  • property (str) – property name that have been computed

  • esp (float) – spacing to use when sampling the array along the cross section

  • ax (matplotlib axes on which to plot)

  • ratio_aspect (float) – ratio between y and x axis to adjust vertical exaggeration

plot_facies(iu=0, ifa=0, v_ex=1, inside_units=None, excludedVal=None, plotter=None, slicex=None, slicey=None, slicez=None, scalar_bar_kwargs=None, show_scalar_bar=True, **kwargs)

Plot the facies realizations over the domain with the colors attributed to facies

Parameters:

  • iu (int) – indice of units realization

  • ifa (int) – indice of facies realziation

  • v_ex (int or float) – vertical exageration

  • inside_units (array-like of Unit objects) – list of units inside which we want to have the plot. if None –> all

  • plotter (pyvista plotter if wanted)

  • slicex (float or sequence of floats) – fraction of x axis where to slice

  • slicey (float or sequence of floats) – fraction of y axis where to slice

  • slicez (float or sequence of floats) – fraction of z axis where to slice

  • scalar_bar_kwargs (dict) – kwargs for the scalar bar

  • show_scalar_bar (bool) – if True, show the scalar bar

  • kwargs (dict) – kwargs to pass to geone.imgplot.drawImage3D_slice or geone.imgplot.drawImage3D_surface

plot_geol_map(plotter=None, v_ex=1, up=0)
plot_grid(v_ex=1)

Function that plots the grid of the simulation domain

Parameters:

v_ex (float) – vertical exaggeration

plot_lines(list_lines, names=None, ax=None, legend=True)

Plot cross-sections lines given in the “list_lines” argument on a 2D top view of the domain

Parameters:

  • list_lines (seqence of points) – list of points defining the cross-section lines. Each point is a list of 2 floats [x, y]

  • names (sequence of string) – list of names to give to the cross-section lines

  • ax (matplotlib axis) – axis on which to plot the cross-section

  • legend (bool) – if True, plot a legend

plot_mean_prop(property, type='arithmetic', v_ex=1, inside_units=None, inside_facies=None, filtering_interval=None, plotter=None, slicex=None, slicey=None, slicez=None, cmin=None, cmax=None, scalar_bar_kwargs=None, **kwargs)

Function that plots the arithmetic mean of a property at every cells of the simulation domain given all the property simulations

Parameters:

  • property (string) – name of the property to plot

  • type (string) – type of mean to plot: - “arithmetic” for arithmetic mean - “std” for standard deviation - “median” for median value

  • v_ex (float) – vertical exaggeration

  • inside_units (list of string or Unit objects) – list of units to consider for the mean “std” for standard deviation, “median” for median value

  • inside_facies (list of string or Facies objects) – list of facies to consider for the mean

  • filtering_interval (list of two floats) – interval of values to keep

  • plotter (pyvista.Plotter) – pyvista plotter to plot on

  • slicex (float or sequence of floats) – fraction of the x axis to plot

  • slicey (float or sequence of floats) – fraction of the y axis to plot

  • slicez (float or sequence of floats) – fraction of the z axis to plot

  • cmin (float) – minimum value for the colorbar

  • cmax (float) – maximum value for the colorbar

  • scalar_bar_kwargs (dict) – dictionary of arguments to pass to the pyvista scalar bar

plot_proba(obj, v_ex=1, plotter=None, filtering_interval=[0.01, 1.0], slicex=None, slicey=None, slicez=None, scalar_bar_kwargs=None, excludedVal=None, **kwargs)

Plot the probability of occurence of a specific unit or facies (can be passed by a name or the object directly)

Parameters:

  • obj (Unit or Facies or str) – unit or facies object or string name of the unit/facies to plot the probability of occurence

  • v_ex (float) – vertical exageration, default is 1

  • plotter (pyvista.Plotter) – pyvista.Plotter object to plot on

  • slicex (float or sequence of floats) – fraction of x axis where to slice

  • slicey (float or sequence of floats) – fraction of y axis where to slice

  • slicez (float or sequence of floats) – fraction of z axis where to slice

  • scalar_bar_kwargs (dict) – kwargs to pass to the scalar bar

  • excludedVal (float or sequence of floats) – values to exclude from the plot

  • kwargs (dict) – kwargs to pass to geone.imgplot.drawImage3D_slice or geone.imgplot.drawImage3D_surface

plot_prop(property, iu=0, ifa=0, ip=0, v_ex=1, inside_units=None, inside_facies=None, filtering_interval=None, plotter=None, slicex=None, slicey=None, slicez=None, cmin=None, cmax=None, scalar_bar_kwargs=None, **kwargs)

Plot the facies realizations over the domain with the colors attributed to facies

Parameters:

  • iu (int) – indice of units realization

  • ifa (int) – indice of facies realization

  • ip (int) – indice of property realization

  • v_ex (int or float) – vertical exageration

  • inside_units (array-like of Unit objects or unit names (string)) – list of units inside which we want to have the plot. By default all

  • inside_facies (array-like of Facies objects or facies names (string)) – list of facies inside which we want to have the plot. By default all

  • plotter (pyvista plotter if wanted)

  • slicex, slicey, slicez (float or sequence of floats) – fraction in x, y or z direction where the slice is done

  • cmin, cmax (float) – min and max value of the colorbar

  • scalar_bar_kwargs (dict) – kwargs for the scalar bar

  • kwargs (dict) – kwargs to pass to geone.imgplot.drawImage3D_slice or geone.imgplot.drawImage3D_surface

plot_units(iu=0, v_ex=1, plotter=None, h_level='all', slicex=None, slicey=None, slicez=None, excludedVal=None, scalar_bar_kwargs=None, show_scalar_bar=True, **kwargs)

plot units domain for a specific realization iu

Parameters:

  • iu (int) – unit index to plot

  • v_ex (float) – vertical exageration

  • plotter (pyvista.Plotter) – pyvista.Plotter object to plot on

  • h_level (string or int) – hierarchical level to plot. A value of 1 indicates that only unit of the master pile will be plotted. “all” to plot all possible units

  • slicex (float or sequence of floats) – fraction of x axis where to slice

  • slicey (float or sequence of floats) – fraction of y axis where to slice

  • slicez (float or sequence of floats) – fraction of z axis where to slice

  • scalar_bar_kwargs (dict) – kwargs to pass to the scalar bar

  • show_scalar_bar (bool) – show scalar bar or not

  • kwargs (dict) – kwargs to pass to geone.imgplot.drawImage3D_slice or geone.imgplot.drawImage3D_surface

pointToIndex(x, y, z)

Return the index of the cell containing the point (x,y,z)

Parameters:

  • x (float) – x coordinate of the point

  • y (float) – y coordinate of the point

  • z (float) – z coordinate of the point

Returns:

index of the cell containing the point (x,y,z)

Return type:

3 int

process_bhs(step=None, facies=True, stop_condition=False, reprocess=False)

ArchPy Pre-processing algorithm Extract hard data from boreholes for all units given the Piles defined in the Arch_table object.

Parameters:

  • step (float) – vertical interval for extracting borehole

  • facies (bool) – flag to indicate to process facies data or not

  • stop_condition (bool) – flag to indicate if the process must be aborted if a inconsistency in bhs is found (False) or bhs will be simply ignored (True)

process_geological_map(typ='all', step=5)

Process the geological map attributed to ArchTable model. This function creates fake boreholes from a given geological map (raster) and them to list_map_bhs

Parameters:

  • typ (str) – flag to indicate what information to take: - “uniform” for only superficial information (no contact or boundaries) - “boundaries” for only the contact between the units - “all” for both

  • step (int) – step for sampling the geological map, small values implies that much more data are sampled from the raster but this increases the computational burden. Default is 5 (every 5th cell is sampled)

realizations_aggregation(method='basic', depth=100, ignore_units=None, units_to_fill=[], n_iter=50)

Method to aggregate multiple ArchPy realizations into one for units and facies (TO DO)

Parameters:

  • method (str) – method to use to aggregate the realizations valid method are:

    • basic, return a model with the most probable units/facies in each cells

    • probas_prop, return a model constructed sequentially by ensuring that proportions are respected (at best…)

    • mean_surfs, return a model created by meaning the surface elevation if units were simulated with categorical method, basic method is used for these units

  • depth (float) – probas_prop parameter, maximum depth of investigation to compute probas and proportions. Should be around the median depth of the boreholes

  • ignore_units (list) – probas_prop parameter, units name to ignore. These units will not be aggregated in the final model.

  • units_to_fill (list) – mean_surfs parameter, units name to fill with NN at the end. Should not be used except to fill the top unit.

Return type:

3D ndarray of size (nz, ny, nx)

rem_all_bhs(standard_bh=True, fake_bh=True, geol_map_bh=True)

Remove all boreholes from the list

Parameters:

  • fake_only (bool) – if True, only fake boreholes are removed

  • geol_map_only (bool) – if True, only boreholes from geological map are removed

rem_all_facies_from_units()

To remove facies from all units

rem_bh(bh)

Remove a given bh from the list of boreholes

Parameters:

bh (borehole object) – borehole to remove

rem_fake_bh(bh)

Remove a given bh from the list of fake boreholes

Parameters:

bh (borehole object) – borehole to remove

reprocess()

Reprocess the boreholes and erase the previous hard data

resample2grid(raster_path, band=None, rspl_method='nearest')

Resample a raster to the size of the simulation grid.

Parameters:

  • raster_path (str) – path to the raster file

  • band (int, optional) – raster band to use, if None 0 is used. The default is None.

  • rspl_method (str, optional) – resampling method to use. availables are : nearest, linear, cubic. The default is “nearest”.

Return type:

2D array of size (self.ny, self.nx)

set_Pile_master(Pile_master)

Define a pile object as the main pile of the project

Parameters:

Pile_master (Pile object) – pile object to set as the main pile of the project

unit_mask(unit_name, iu=0, all_real=False)

Return the mask of the given unit for a realization iu

Parameters:

  • unit_name (string) – unit name defined when creating the unit object for more details: Unit

  • iu (int) – unit realization index (0, 1, …, Nu)

  • all_real (bool) – flag to know if a mask of all realizations must be returned

Returns:

mask array

Return type:

3D (or 4D) ndarray

base.Arr_replace(arr, dic)

Replace value in an array using a dictionnary linking actual values to new values

Parameters:

  • arr (np.ndarray) – Any numpy array that contains values

  • dic (dict) – A dictionnary with arr values as keys and new values as dic values

Returns:

A new array with values replaced

Return type:

nd.array

class base.Facies(ID, name, color)

Bases: object

class for facies (2nd level of hierarchy)

Parameters:

  • ID (int) – facies ID that is used in the results

  • name (string) – facies name, used as an identifier

  • color (string) – color of the facies, used for plotting

class base.Geol

Bases: object

ArchPy output class which contain the results for the geological simulations

class base.Pile(name, nature='surfaces', verbose=1, seed=1)

Bases: object

This class is a major object the Stratigraphic Pile (SP). It contains all the units objects and allow to know the stratigraphical relations between the units. One Pile object must be defined for each subpile + 1 “master” pile.

Parameters:

  • name (str) – name of the pile

  • nature (str) – units type of interpolation, can be “surfaces” or “3d_categorical” if “surfaces” is chosen (default), 2D surfaces interpolations of the surfaces are performed. The surfaces are then used to delimit unit domains. if “3d_categorical”, a facies method is used to simulate position of the units. The available methods are the same to simulate the facies.

  • verbose (int) – level of verbosity

  • seed (int) – seed for random number generation

add_unit(unit)

Add a unit to the pile

Parameters:

unit (Unit object or list of Unit objects) – unit(s) to add to the pile

compute_surf(ArchTable, nreal=1, fl_top=False, subpile=False, tops=None, bots=None, vb=1)

Compute the elevation of the surfaces units (1st hierarchic order) contained in the Pile object.

Parameters:

  • ArchTable (Arch_table object) – ArchTable object containing the architecture of the pile

  • nreal (int) – number of realization

  • fl_top (bool) – to not interpolate first layer and assign it=top

  • subpile (bool) – if the pile is a subpile

  • tops, bots (sequence of 2D arrays of size (ny, nx)) – of top/bot for subpile surface

  • vb (int) – level of verbosity, 0 for no print, 1 for print

define_domains(ArchTable, surfaces, tops=None, bots=None, subpile=False, vb=0, fl_top=True)

Define units domains from precomputed surfaces. This method is used when the surfaces are precomputed (e.g. from a previous simulation) in place to the compute_surfaces method.

Parameters:

  • ArchTable (Arch_table) – ArchTable object

  • surfaces (array) – array of surfaces (nreal, nlay, ny, nx)

  • tops (array) – array of top elevations (ny, nx)

  • bots (array) – array of bottom elevations (ny, nx)

  • subpile (bool) – if the pile is a subpile, i.e. if True, the top and bottom of the pile are not used

  • vb (int) – verbosity level, 0 is silent, 1 is verbose

  • fl_top (bool) – if True, the top of the pile is used to define the top of the first unit

get_pile_ref()
order_units(vb=1)

Order list_liths according the order attributes of each lithologies

Parameters:

vb (int) – level of verbosity, 0 for no print, 1 for print

remove_unit(unit_to_rem)

Remove the given unit object from Pile object

Parameters:

unit_to_rem (Unit object) – unit to remove from the pile

class base.Prop(name, facies, covmodels, means, int_method='sgs', x=None, v=None, def_mean=1, vmin=None, vmax=None)

Bases: object

Class for defining a propriety to simulate (3rd level)

Parameters:

  • name (string) – Property name

  • facies (list of Facies objects) – facies in which we want to simulate the property (in the others the propriety will be set homogenous with a default value)

  • covmodels (list of geone.covModel.CovModel2D objects) – covmodels for the simulation (same size of facies or only 1)

  • means (list of floats) – mean of the property in each facies (same size of facies)

  • int_method (string) – method of interpolation, possible methods are: - sgs, default - fft - homogenous

  • def_mean (float) – default mean to used if none is passed in means array

  • vmin, vmax (float) – min resp. max value of the property

  • x (ndarray of size (n, 2 –> x, y)) – position of hard data

  • v (ndarray of size (n, 1 –> value)) – value of hard data

add_hd(x, v)

add hard data to the property

Parameters:

  • x (ndarray of size (n, 2 –> x, y)) – position of hard data

  • v (ndarray of size (n, 1 –> value)) – value of hard data

class base.Surface(name='Surface_1', dic_surf={'N_transfo': False, 'covmodel': None, 'int_method': 'nearest'}, contact='onlap')

Bases: object

Class Surface, must be linked to a Unit object

Parameters:

  • name (string) – to identify the surface for debugging purpose

  • contact (string) – onlap or erode. Onlap indicates that this surface cannot erode older surfaces, on contrary to erode surfaces

  • dic_surf (dict) –

    parameters for surface interpolation:

    • int_methodstring

      method of interpolation, possible methods are: kriging, MPS, grf, grf_ineq, linear, nearest,

    • covmodelstring

      geone covariance model (see doc Geone for more information) if a multi-gaussian method is used

    • N_transfobool

      Normal-score transform. Flag to apply or not a Normal Score on the data for the interpolation

    • Units kwargs (passed in dic_surf directly):

    for the search ellipsoid (kriging, grf, …):

    • r, relative (to covariance model) radius of research (default is 1)

    • neig, number of neighbours

    • krig_type: string, kriging method (ordinary_kriging, simple_kriging)

    for MPS:

    • TI

    • various MPS parameters, see geone documentation

copy()
get_surface_covmodel(vb=1)
set_contact(contact)
set_covmodel(covmodel, vb=0)

change or add a covmodel for surface interpolation of self unit.

Parameters:

covmodel (geone.covModel.CovModel2D) – covariance model to be used for surface interpolation if the chosen method is grf, grf ineq or kriging

set_dic_surf(dic_surf)
class base.Unit(name, order, color, surface=None, ID=None, m_unit_name=None, dic_facies={'Flag': None, 'G_cm': None, 'SubPile': None, 'TI': None, 'f_covmodel': None, 'f_method': 'homogenous'}, contact='onlap', verbose=1)

Bases: object

This class defines Unit objects, which are the building blocks of the Pile class.

Parameters:

  • name (string) – name of the unit that will be used as an identifier

  • order (int) – order of the unit in the pile (1 (top) to n (bottom), where n is the total number of units)

  • color (string) – color to use for plotting and representation. The color can be any color that is accepted by matplotlib

  • surface (Surface object) – surface object that defines the surface of the unit

  • ID (int) – ID of the unit, if None, the ID will be set to the order of the unit

  • m_unit_name (string) – identifier name if the unit is a multiple unit (can appear several times in the pile) Not yet implemented

  • dic_facies (dict) – dictionary that defines the facies of the unit. The mandatory keys for the dictionary are:

    • f_method : string, valable method are:

      • “homogenous” : homogenous facies

      • “SubPile” : facies are simulated with a subpile

      • “SIS” : facies are simulated with stochastic indicator simulation

      • “MPS” : facies are simulated with Multiple points statistics

      • “TPGs”: facies are simulated with Truncated pluri-Gaussian method

    • f_covmodel : 3D geone covariance model object or list of 3D geone covariance model objects only used with “SIS” method

    • TI : geone.image, Training image for “MPS” method

    • SubPile : Pile object, subpile for “SubPile” method

    Facies keyword arguments to pass to the facies methods (these should be pass through dic_facies !):

    • SIS:

      • “neigh” : int, number of neighbors to use for the SIS

      • r : float, radius of the neighborhood to use for the SIS

      • probability : list of float, proportion of each facies to use for the SIS

    • MPS:

      • TI : geone img, Training image(s) to use

      • mps “classic” parameters (maxscan, thresh, neig (number of neighbours))

      • npost : number of path postprocessing, default 1

      • radiusMode : radius mode to use (check deesse manual)

      • anisotropyRatioMode : Anisotropy ratio to use in reasearch of neighbours (check deesse manual)

      • rx, ry, rz : radius in x, y and z direction

      • angle1, angle2, angle3: ellipsoid of research rotation angles

      • ax, ay, az : anisotropy ratio for research ellipsoid

      • rot_usage : 0, 1 or 2 (check deesse manual)

      • rotAziLoc, rotDipLoc, rotPlungeLoc : local rotation, True or False ?

      • rotAzi, rotDip, rotPlunge : global rotation angles: values, min-max, maps, see deesse_input doc

      • xr, yr, zr : ratio for geom transformation

      • xloc, yloc, zloc : local or not transformation

      • homo_usage : homothety usage

      • probaUsage : probability constraint usage, 0 for no proba constraint, 1 for global proportion defined in globalPdf, 2 for local proportion defined in localPdf

      • globalPdf : array-like of float of length equal to the number of class, proportion for each class

      • localPdf : (nclass, nz, ny, nx) array of floats probability for each class, localPdf[i] is the “map defined on the simulation grid

      • localPdfRadius : support radius for local pdf, default is 2

      • deactivationDistance : float, distance at which localPdf are deactivated (see Deesse doc)

      • constantThreshold : float, threshold value for pdf’s comparison

      • dataImage : geone img used as data, see deesse/geone documentation

      • outputVarFlag : bool or list of bool of size nv, to output or not the variables

      • distanceType : string or list of strings, “categorical” or “continuous”

    • TPGs:

      • flag: dictionary of limits of cuboids domains in Gaussian space, check … for the right format.

      • Ik_cm: indicator covmodels

      • G_cm: list of Gaussian covmodels for the two gaussian fields, can be infered from Ik_cm

      • various parameters: (du, dv –> precision of integrals for G_cm inference),

        (dim: dimension wanted for G_cm inference, (c_reg: regularisation term for G_cm inference), (n: number of control points for G_cm inference))

add_facies(facies)

Add facies object to strati

Parameters:

facies (Facies object or list of them) – Facies object to model inside the unit

compute_facies(ArchTable, nreal=1, mode='facies', verbose=0)

Compute facies domain for the specific unit

Parameters:

  • ArchTable (Arch_table object) – ArchTable containing units, surface, facies and at least a Pile (see example on the github)

  • nreal (int) – number of realization (per unit realizations) to make

  • mode (str) – “facies” or “strati”

  • verbose (int) – verbosity level, 0 for no print, 1 for print

copy()
get_baby_units(recompute=False, vb=1)

Method to return all units that are under the current unit in the hierarchy

Parameters:

  • recompute (bool) – If True, recompute the list of units if there had been modifications in the hierarchy

  • vb (int) – verbosity level, 0 for no print, 1 for print

Returns:

  • list

  • list of childs Unit objects of the current unit

get_big_mummy_unit()

Return lowest unit (main unit) in hierarchical order

get_h_level()

Return hierarchical level of the unit

goes_up_until(unit_target)

Climb the hierarchical tree of self unit until finding a specific unit and return it

Parameters:

unit_target (Unit object) – Unit object to find in the hierarchical tree

Returns:

unit – Unit object found in the hierarchical tree if not found, return None

Return type:

Unit object or None

rem_facies(facies=None, all_facies=True)

To remove facies from unit, by default all facies are removed

Parameters:

  • facies (Facies object or list of them) – Facies object to remove from the unit

  • all_facies (bool) – If True, all facies are removed

set_SubPile(SubPile)

Change or define a subpile for filling

Parameters:

SubPile (Pile object)

set_Surface(surface)

Change or add a top surface for Unit object. This will define how the top of the formation will be modelled.

Parameters:

surface (Surface object) – Surface object to use for the simulation of the top of the formation

set_dic_facies(dic_facies)

Set a dictionary facies to Unit object

Parameters:

dic_facies (dict) – dictionnary containing the parameters to pass at the facies methods

set_f_TI(TI)

Change training image for a strati unit

Parameters:

TI (geone.image) – Training image to use with the MPS

set_f_covmodels(f_covmodel)

Remove existing facies covmodels and one or more depending of what is passed (array like or only one object)

Parameters:

f_covmodel (geone.covModel object) – that will be used for the interpolation of the facies if method is SIS. Can be a list or only one object.

class base.borehole(name, ID, x, y, z, depth, log_strati, log_facies=None)

Bases: object

Class to create a borehole object. Borehole are used to define the lithology and add conditioning data into ArchPy models

Parameters:

  • name (string) – name of the borehole, not used

  • ID (int) – ID of the borehole

  • x, y, z (float) – x, y, z coordinates of the top of the borehole

  • depth (float) – depth of the borehole

  • log_strati (list of tuples) – log_strati contains the geological information about the units in the borehole information is given by intervals of the form (strati, top) where strati is a Unit object and top is the top altitude of unit interval

  • log_facies (list of tuples) – log_facies contains the facies information about the borehole information is given by intervals of the form (facies, top) where facies is a Facies object and top is the top altitude of facies interval

create_thickness_log(pile_ref)

Create the thickness log based on the pile reference and the raw log

extract(z, vb=1)

extract the units and facies information at specified altitude z

get_list_facies()
get_list_stratis()
prop_units()

Return a dictionnary of the proportion of the units in the borehole

set_log_facies(log_facies)
set_log_strati(log_strati)
set_thickness_log(log_thk)
base.get_size(obj, seen=None)

Recursively finds size of objects

Parameters:

obj (any python object) – The object, variables or anything whose size we are looking for

Returns:

Size of the object in bytes

Return type:

int

Note

Function taken from stack overflow, written by Aaron Hall

base.interp2D(litho, xg, yg, xu, verbose=0, ncpu=1, mask2D=None, seed=123456789, **kwargs)

Function to realize a 2D interpolation based on a multitude of methods (scipy.interpolate, geone (kriging, MPS, …))

Parameters:

  • litho (Surface object) – ArchPy surface object on which we want to interpolate the surface

  • xg (ndarray of size nx) – central coordinates in x direction

  • yg (ndarray of size ny) – central coordinates in y direction

  • xu (ndarray of size (n, 2)) – position at which we want to know the estimation (not used for every interp method)

  • **kwargs – Different parameters for surface interpolation

    • nit: int

      Number of iterations for gibbs sampler (depends on the number of data)

    • nmax: int

      Number of neighbours in grf with inequalities (to speed up the simulations)

    • krig_type: str

      Can be either “ordinary_kigring” or “simple kriging” method for kriging interpolation and grf with inequalities

    • mean: float or 2D array

      Mean value to use with grf methods

    • unco: bool

      Unconditional or not

    • All other MPS and GRF parameters (see geone documentation)

Returns:

array of same size as x, interpolated values

Return type:

ndarray

base.resample_to_grid(xc, yc, rxc, ryc, raster_band, method='nearest')

Function to resample the raster data to a user supplied grid of x, y coordinates.

x, y coordinate arrays should correspond to grid vertices

Parameters:

  • xc (np.ndarray or list) – an array of x-cell centers

  • yc (np.ndarray or list) – an array of y-cell centers

  • rxc (ndarray) – raster xcell centers

  • ryc (ndarray) – raster ycell centers

  • raster_band (2D ndarray) – raster band to re-sample

  • method (str) – scipy interpolation method options

    “linear” for bi-linear interpolation “nearest” for nearest neighbor “cubic” for bi-cubic interpolation

Return type:

np.array

Note

Function taken from flopy (3.3.4)

base.running_mean_2D(x, N)

Smooth a 2d surface

Parameters:

  • x (2D ndarray) – Array to smooth

  • N (int) – Window semi-size

Return type:

2D ndarray

base.split_logs(bh)

Take a raw borehole with hierarchical units mixed and sort them by group and hierarchy.

Parameters:

bh (borehole object)

Returns:

A list of new boreholes

Return type:

list