Code

Module with classes and methods to perform implicit regional modelling based on the potential field method. Tested on Ubuntu 14

Created on 10/10 /2016

@author: Miguel de la Varga

gempy_front.compute_model(interp_data, output=’geology’, u_grade=None, get_potential_at_interfaces=False)[source]

Computes the geological model.

Parameters:
  • interp_data (gempy.DataManagement.InterpolatorInput) – Rescaled data.
  • u_grade (list) – grade of the polynomial for the universal part of the Kriging interpolations. The value has to
  • either 0, 3 or 9 (be) –
  • on the number of points given as input to try to avoid singular matrix. NOTE (depends) – if during the computation
  • the model a singular matrix is returned try to reduce the u_grade of the series. (of) –
  • get_potential_at_interfaces (bool) – Get potential at interfaces
Returns:

if compute_all was chosen in gempy.DataManagement.InterpolatorInput, the first row will be the lithology block model, the second the potential field and the third the fault network block. if compute_all was False only the lithology block model will be computed. In addition if get_potential_at_interfaces is True, the value of the potential field at each of the interfaces is given as well

Return type:

numpy.array

gempy_front.create_data(extent, resolution=[50, 50, 50], **kwargs)[source]

Method to create a InputData object. It is analogous to gempy.InputData()

Parameters:
  • extent (list or array) – [x_min, x_max, y_min, y_max, z_min, z_max]. Extent for the visualization of data and default of for the grid class.
  • resolution (list or array) – [nx, ny, nz]. Resolution for the visualization of data and default of for the grid class.
Keyword Arguments:
 
  • Resolution ((Optional[list])) – [nx, ny, nz]. Defaults to 50
  • path_i – Path to the data bases of interfaces. Default os.getcwd(),
  • path_f – Path to the data bases of orientations. Default os.getcwd()
Returns:

Object that encapsulate all raw data of the project

dep: self.Plot(GeMpy_core.PlotData): Object to visualize data and results

Return type:

GeMpy.DataManagement

gempy_front.data_to_pickle(geo_data, path=False)[source]

Save InputData object to a python pickle (serialization of python). Be aware that if the dependencies versions used to export and import the pickle differ it may give problems

Parameters:path (str) – path where save the pickle
Returns:None
gempy_front.export_to_vtk(geo_data, path=None, name=None, lith_block=None, vertices=None, simplices=None)[source]

Export data to a vtk file for posterior visualizations

Parameters:
  • geo_data (gempy.InputData) – All values of a DataManagement object
  • block (numpy.array) – 3D array containing the lithology block
  • path (str) – path to the location of the vtk
Returns:

None

gempy_front.get_data(geo_data, dtype=’all’, numeric=False, verbosity=0)[source]

Method that returns the interfaces and orientations pandas Dataframes. Can return both at the same time or only one of the two

Parameters:
  • dtype (str) – input data type, either ‘orientations’, ‘interfaces’ or ‘all’ for both.
  • verbosity (int) – Number of properties shown
Returns:

Data frame with the raw data

Return type:

pandas.core.frame.DataFrame

gempy_front.get_grid(geo_data)[source]
Parameters:geo_data (gempy.DataManagement.InputData object) –
Returns:Return series and formations relations
Return type:numpy.array
gempy_front.get_kriging_parameters(interp_data, verbose=0)[source]

Print the kringing parameters

Parameters:
  • interp_data (gempy.DataManagement.InterpolatorInput) –
  • verbose (int) – if > 0 print all the shape values as well.
Returns:

None

gempy_front.get_sequential_pile(geo_data)[source]

Visualize an interactive stratigraphic pile to move around the formations and the series. IMPORTANT NOTE: To have the interactive properties it is necessary the use of qt as interactive backend. (In notebook use: %matplotlib qt5)

Parameters:geo_data – gempy.DataManagement.InputData object
Returns:interactive Matplotlib figure
gempy_front.get_series(geo_gata)[source]
Parameters:geo_data (gempy.DataManagement.InputData object) –
Returns:Return series and formations relations
Return type:Pandas.DataFrame
gempy_front.get_surfaces(interp_data, potential_lith=None, potential_fault=None, n_formation=’all’, step_size=1, original_scale=True)[source]

compute vertices and simplices of the interfaces for its vtk visualization or further analysis

Parameters:
  • potential_lith (numpy.array) – 1D numpy array with the solution of the computation of the model containing the scalar field of potentials (second row of solution)
  • interp_data (gempy.DataManagement.InterpolatorInput) – Interpolator object.
  • n_formation (int or 'all') – Positive integer with the number of the formation of which the surface is returned. use method get_formation_number() to get a dictionary back with the values
  • step_size (int) – resolution of the method. This is every how many voxels the marching cube method is applied
  • original_scale (bool) – choosing if the coordinates of the vertices are given in the original or the rescaled coordinates
Returns:

vertices, simpleces

gempy_front.get_th_fn(interp_data, compute_all=True)[source]

Get theano function

Parameters:interp_data (gempy.DataManagement.InterpolatorInput) –

Rescaled data. compute_all (bool): If true the solution gives back the block model of lithologies, the potential field and the block model of faults. If False only return the block model of lithologies. This may be important to speed

up the computation. Default True
Returns:
Compiled function if C or CUDA which computes the interpolation given the input data
(XYZ of dips, dip, azimuth, polarity, XYZ ref interfaces, XYZ rest interfaces)
Return type:theano.function
gempy_front.plot_data(geo_data, direction=’y’, data_type=’all’, series=’all’, legend_font_size=6, **kwargs)[source]

Plot the projection of the raw data (interfaces and orientations) in 2D following a specific directions

Parameters:
  • direction (str) – xyz. Caartesian direction to be plotted
  • series (str) – series to plot
  • **kwargs – seaborn lmplot key arguments. (TODO: adding the link to them)
Returns:

None

gempy_front.plot_data_3D(geo_data)[source]

Plot in vtk all the input data of a model :param geo_data: Input data of the model :type geo_data: gempy.DataManagement.InputData

Returns:None
gempy_front.plot_scalar_field(geo_data, potential_field, cell_number, N=20, direction=’y’, plot_data=True, series=’all’, *args, **kwargs)[source]

Plot a potential field in a given direction.

Parameters:
  • cell_number (int) – position of the array to plot
  • potential_field (str) – name of the potential field (or series) to plot
  • n_pf (int) – number of the potential field (or series) to plot
  • direction (str) – xyz. Caartesian direction to be plotted
  • serieDeprecated
  • **kwargs – plt.contour kwargs
Returns:

None

gempy_front.plot_section(geo_data, block, cell_number, direction=’y’, **kwargs)[source]

Plot a section of the block model

Parameters:
  • cell_number (int) – position of the array to plot
  • direction (str) – xyz. Caartesian direction to be plotted
  • interpolation (str) – Type of interpolation of plt.imshow. Default ‘none’. Acceptable values are ‘none’
  • 'bilinear', 'bicubic', (,'nearest',) –
  • 'spline36', 'hanning', 'hamming', 'hermite', 'kaiser', ('spline16',) –
  • 'catrom', 'gaussian', 'bessel', 'mitchell', 'sinc', ('quadric',) –
  • 'lanczos'**kwargs: imshow keywargs
Returns:

None

gempy_front.plot_surfaces_3D(geo_data, vertices_l, simplices_l, alpha=1, plot_data=True, size=(1920, 1080), fullscreen=False, bg_color=None)[source]

Plot in vtk the surfaces. For getting vertices and simplices See gempy.get_surfaces

Parameters:
  • vertices_l (numpy.array) – 2D array (XYZ) with the coordinates of the points
  • simplices_l (numpy.array) – 2D array with the value of the vertices that form every single triangle
  • formations_names_l (list) – Name of the formation of the surfaces
  • formation_numbers_l (list) – Formation numbers (int)
  • alpha (float) – Opacity
  • plot_data (bool) – Default True
  • size (tuple) – Resolution of the window
  • fullscreen (bool) – Launch window in full screen or not
Returns:

None

gempy_front.plot_surfaces_3D_real_time(interp_data, vertices_l, simplices_l, alpha=1, plot_data=True, posterior=None, samples=None, size=(1920, 1080), fullscreen=False)[source]

Plot in vtk the surfaces in real time. Moving the input data will affect the surfaces. IMPORTANT NOTE it is highly recommended to have the flag fast_run in the theano optimization. Also note that the time needed to compute each model increases linearly with every potential field (i.e. fault or discontinuity). It may be better to just modify each potential field individually to increase the speed (See gempy.select_series).

Parameters:
  • vertices_l (numpy.array) – 2D array (XYZ) with the coordinates of the points
  • simplices_l (numpy.array) – 2D array with the value of the vertices that form every single triangle
  • formations_names_l (list) – Name of the formation of the surfaces
  • formation_numbers_l (list) – Formation numbers (int)
  • alpha (float) – Opacity
  • plot_data (bool) – Default True
  • size (tuple) – Resolution of the window
  • fullscreen (bool) – Launch window in full screen or not
Returns:

None

gempy_front.plot_topology(geo_data, G, centroids, direction=’y’)[source]

Plot the topology adjacency graph in 2-D.

Parameters:
  • geo_data (gempy.data_management.InputData) –
  • G (skimage.future.graph.rag.RAG) –
  • centroids (dict) – Centroid node coordinates as a dictionary with node id’s (int) as keys and (x,y,z) coordinates as values.
Keyword Args
direction (str): “x”, “y” or “z” specifying the slice direction for 2-D topology analysis. Default None.
Returns:Nothing, it just plots.
gempy_front.read_pickle(path)[source]

Read InputData object from python pickle.

Parameters:path (str) – path where save the pickle
Returns:Input Data object
Return type:gempy.DataManagement.InputData
gempy_front.rescale_data(geo_data, rescaling_factor=None)[source]

Rescale the data of a DataManagement object between 0 and 1 due to stability problem of the float32.

Parameters:
  • geo_data – Original gempy.DataManagement.InputData object
  • rescaling_factor (float) – factor of the rescaling. Default to maximum distance in one the axis
Returns:

Rescaled data

Return type:

gempy.data_management.InputData

gempy_front.rescale_factor_default(geo_data)[source]

Gives the default rescaling factor for a given geo_data

Parameters:geo_data – Original gempy.DataManagement.InputData object
Returns:rescaling factor
Return type:float
gempy_front.select_series(geo_data, series)[source]

Return the formations of a given serie in string

Parameters:series – list of int or list of str
Returns:formations of a given serie in string separeted by |
gempy_front.set_interfaces(geo_data, interf_Dataframe, append=False)[source]

Method to change or append a Dataframe to interfaces in place.

Parameters:
  • interf_Dataframe – pandas.core.frame.DataFrame with the data
  • append – Bool: if you want to append the new data frame or substitute it
gempy_front.set_interpolation_data(geo_data, **kwargs)[source]

InterpolatorInput is a class that contains all the preprocessing operations to prepare the data to compute the model. Also is the object that has to be manipulated to vary the data without recompile the modeling function.

Parameters:
  • geo_data (gempy.DataManagement.InputData) – All values of a DataManagement object
  • compile_theano (bool) – select if the theano function is compiled during the initialization. Default: True
  • compute_all (bool) –

    If true the solution gives back the block model of lithologies, the potential field and the block model of faults. If False only return the block model of lithologies. This may be important to speed

    up the computation. Default True
  • u_grade (list) – grade of the polynomial for the universal part of the Kriging interpolations. The value has to
  • either 0, 3 or 9 (be) –
  • on the number of points given as input to try to avoid singular matrix. NOTE (depends) – if during the computation
  • the model a singular matrix is returned try to reduce the u_grade of the series. (of) –
  • rescaling_factor (float) – rescaling factor of the input data to improve the stability when float32 is used. By
  • the rescaling factor is calculated to obtein values between 0 and 1. (defaut) –
Keyword Arguments:
 
  • dtype ('str') – Choosing if using float32 or float64. This is important if is intended to use the GPU
  • Also InterpolatorClass kwargs (See) –
geo_data

Original gempy.DataManagement.InputData object

geo_data_res

Rescaled data. It has the same structure has gempy.InputData

interpolator

Instance of the gempy.DataManagement.InterpolaorInput.InterpolatorClass. See Also gempy.DataManagement.InterpolaorInput.InterpolatorClass docs th_fn: Theano function which compute the interpolation

dtype

type of float

gempy_front.set_orientations(geo_data, orient_Dataframe, append=False)[source]

Method to change or append a Dataframe to orientations in place. A equivalent Pandas Dataframe with [‘X’, ‘Y’, ‘Z’, ‘dip’, ‘azimuth’, ‘polarity’, ‘formation’] has to be passed.

Parameters:
  • interf_Dataframe – pandas.core.frame.DataFrame with the data
  • append – Bool: if you want to append the new data frame or substitute it
gempy_front.set_series(geo_data, series_distribution=None, order_series=None, order_formations=None, verbose=0)[source]

Method to define the different series of the project.

Parameters:
  • series_distribution (dict) – with the name of the serie as key and the name of the formations as values.
  • order (Optional[list]) – order of the series by default takes the dictionary keys which until python 3.6 are random. This is important to set the erosion relations between the different series
Returns:

A pandas DataFrame with the series and formations relations self.interfaces: one extra column with the given series self.orientations: one extra column with the given series

Return type:

self.series

gempy_front.topology_compute(geo_data, lith_block, fault_block, cell_number=None, direction=None, compute_areas=False, return_label_block=False)[source]

Computes model topology and returns graph, centroids and look-up-tables.

Parameters:
  • geo_data (gempy.data_management.InputData) – GemPy’s data object for the model.
  • lith_block (np.ndarray) – Lithology block model.
  • fault_block (np.ndarray) – Fault block model.
Keyword Arguments:
 
  • cell_number (int) – Cell number for 2-D slice topology analysis. Default None.
  • direction (str) – “x”, “y” or “z” specifying the slice direction for 2-D topology analysis. Default None.
  • compute_areas (bool) – If True computes adjacency areas for connected nodes in voxel number. Default False.
  • return_label_block (bool) – If True additionally returns the uniquely labeled block model as np.ndarray. Default False.
Returns:

G: Region adjacency graph object (skimage.future.graph.rag.RAG) containing the adjacency topology graph

(G.adj).

centroids (dict): Centroid node coordinates as a dictionary with node id’s (int) as keys and (x,y,z) coordinates

as values. {node id (int): tuple(x,y,z)}

labels_unique (np.array): List of all labels used. lith_to_labels_lot (dict): Dictionary look-up-table to go from lithology id to node id. labels_to_lith_lot (dict): Dictionary look-up-table to go from node id to lithology id.

Return type:

tuple

This file is part of gempy.

gempy is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

gempy is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with gempy. If not, see <http://www.gnu.org/licenses/>.

class data_management.GridClass[source]

Class to generate grids to pass later on to a InputData class.

create_custom_grid(custom_grid)[source]

Give the coordinates of an external generated grid

Parameters:custom_grid (numpy.ndarray like) – XYZ (in columns) of the desired coordinates
Returns:Unraveled 3D numpy array where every row correspond to the xyz coordinates of a regular grid
Return type:numpy.ndarray
create_regular_grid_3d(extent, resolution)[source]

Method to create a 3D regular grid where is interpolated

Parameters:
  • extent (list) – [x_min, x_max, y_min, y_max, z_min, z_max]
  • resolution (list) – [nx, ny, nz].
Returns:

Unraveled 3D numpy array where every row correspond to the xyz coordinates of a regular grid

Return type:

numpy.ndarray

class data_management.InputData(extent, resolution=[50, 50, 50], path_i=None, path_o=None, path_f=None, **kwargs)[source]

Class to import the raw data of the model and set data classifications into formations and series. This objects will contain the main information of the model.

Parameters:
  • extent (list) – [x_min, x_max, y_min, y_max, z_min, z_max]
  • Resolution ((Optional[list])) – [nx, ny, nz]. Defaults to 50
  • path_i – Path to the data bases of interfaces. Default os.getcwd(),
  • path_o – Path to the data bases of orientations. Default os.getcwd()
extent

list – [x_min, x_max, y_min, y_max, z_min, z_max]

resolution

(Optional[list]) – [nx, ny, nz]

orientations

pandas.core.frame.DataFrame – Pandas data frame with the orientations data

Interfaces

pandas.core.frame.DataFrame – Pandas data frame with the interfaces data

series

pandas.core.frame.DataFrame – Pandas data frame which contains every formation within each series

add_interface(**kwargs)[source]

Adds interface to dataframe.

Parameters:**kwargs – X, Y, Z, formation, labels, order_series, series
Returns:None
add_orientation(**kwargs)[source]

Adds orientation to dataframe.

Parameters:**kwargs – G_x, G_y, G_z, X, Y, Z, azimuth, dip, formation, labels, order_series, polarity, series

Returns: Nothing

calculate_gradient()[source]

Calculate the gradient vector of module 1 given dip and azimuth to be able to plot the orientations

orientations

extra columns with components xyz of the unity vector.

calculate_orientations()[source]

Calculate and update the orientation data (azimuth and dip) from gradients in the data frame.

count_faults()[source]

Read the string names of the formations to detect automatically the number of faults.

data_to_pickle(path=False)[source]

Save InputData object to a python pickle (serialization of python). Be aware that if the dependencies versions used to export and import the pickle differ it may give problems

Parameters:path (str) – path where save the pickle
Returns:None
drop_interface(index)[source]

Drops interface from dataframe identified by index

Parameters:index – dataframe index
Returns:None
drop_orientations(index)[source]

Drops orientation from dataframe identified by index

Parameters:index – dataframe index
Returns:None
get_data(itype=’all’, numeric=False, verbosity=0)[source]

Method that returns the interfaces and orientations pandas Dataframes. Can return both at the same time or only one of the two

Parameters:
  • itype – input_data data type, either ‘orientations’, ‘interfaces’ or ‘all’ for both.
  • numeric (bool) – Return only the numerical values of the dataframe. This is much lighter database for storing traces
  • verbosity (int) – Number of properties shown
Returns:

Data frame with the raw data

Return type:

pandas.core.frame.DataFrame

get_formation_number()[source]

Get a dictionary with the key the name of the formation and the value their number

Returns:key the name of the formation and the value their number
Return type:dict
get_formations()[source]
Returns:Returns a list of formations
Return type:pandas.core.frame.DataFrame
import_data_csv(path_i, path_o, **kwargs)[source]

Method to import interfaces and orientations from csv. The format is the same as the export 3D model data of GeoModeller (check in the input_data data folder for an example).

Parameters:
  • path_i (str) – path to the csv table
  • path_o (str) – path to the csv table
  • **kwargs – kwargs of :func: ~pn.read_csv
orientations

pandas.core.frame.DataFrame – Pandas data frame with the orientations data

Interfaces

pandas.core.frame.DataFrame – Pandas data frame with the interfaces data

static load_data_csv(data_type, path=’/home/miguel/PycharmProjects/gempy/docs’, **kwargs)[source]

Method to load either interface or orientations data csv files. Normally this is in which GeoModeller exports it

Parameters:
  • data_type (str) – ‘interfaces’ or ‘orientations’
  • path (str) – path to the files. Default os.getcwd()
  • **kwargs – Arbitrary keyword arguments.
Returns:

Data frame with the raw data

Return type:

pandas.core.frame.DataFrame

modify_interface(index, **kwargs)[source]

Allows modification of the x,y and/or z-coordinates of an interface at specified dataframe index.

Parameters:
  • index – dataframe index of the orientation point
  • **kwargs – X, Y, Z (int or float)
Returns:

None

modify_orientation(index, recalculate_gradient=False, recalculate_orientations=False, **kwargs)[source]

Allows modification of orientation data at specified dataframe index.

Parameters:
  • index – dataframe index of the orientation point
  • **kwargs – G_x, G_y, G_z, X, Y, Z, azimuth, dip, formation, labels, order_series, polarity
Returns:

None

order_table()[source]

First we sort the dataframes by the series age. Then we set a unique number for every formation and resort the formations. All inplace

reset_indices()[source]

Resets dataframe indices for orientations and interfaces.

Returns:None
set_annotations()[source]

Add a column in the Dataframes with latex names for each input_data paramenter.

Returns:None
set_fault_relation_matrix(rel_matrix)[source]

Method to set the faults that offset a given sequence and therefore also another fault

Parameters:rel_matrix (numpy.array) – 2D Boolean array with the logic. Rows affect (offset) columns
set_faults(series_name)[source]

Set a flag to the series that are faults.

Parameters:series_name (list or array_like) – Name of the series which are faults
set_formation_number(formation_order=None)[source]

Set a unique number to each formation. NOTE: this method is getting deprecated since the user does not need to know it and also now the numbers must be set in the order of the series as well. Therefore this method has been moved to the interpolator class as preprocessing

Returns:Column in the interfaces and orientations dataframes
set_grid(custom_grid=None, extent=None, resolution=None, grid_type=None, **kwargs)[source]

Method to initialize the class GridClass. You can pass either a custom set of points or create a regular grid

Parameters:
  • grid_type (str) – regular_3D or None
  • custom_grid (array_like) – 2D array with XYZ columns. To exploit gempy functionality the indexing has to be ij (See Also numpy.meshgrid documentation)
  • **kwargs – Arbitrary keyword arguments.
Returns:

Object that contain different grids

Return type:

self.grid(gempy.GridClass)

set_interfaces(interf_Dataframe, append=False)[source]

Method to change or append a Dataframe to interfaces in place. A equivalent Pandas Dataframe with [‘X’, ‘Y’, ‘Z’, ‘formation’] has to be passed.

Parameters:
  • interf_Dataframe – pandas.core.frame.DataFrame with the data
  • append – Bool: if you want to append the new data frame or substitute it
set_orientations(foliat_Dataframe, append=False)[source]
Method to change or append a Dataframe to orientations in place. A equivalent Pandas Dataframe with

[‘X’, ‘Y’, ‘Z’, ‘dip’, ‘azimuth’, ‘polarity’, ‘formation’] has to be passed.

Args:
interf_Dataframe: pandas.core.frame.DataFrame with the data append: Bool: if you want to append the new data frame or substitute it
set_series(series_distribution=None, order=None)[source]

Method to define the different series of the project.

Parameters:
  • series_distribution (dict) – with the name of the serie as key and the name of the formations as values.
  • order (Optional[list]) – order of the series by default takes the dictionary keys which until python 3.6 are random. This is important to set the erosion relations between the different series
Returns:

A pandas DataFrame with the series and formations relations self.interfaces: one extra column with the given series self.orientations: one extra column with the given series

Return type:

self.series

class data_management.InterpolatorData(geo_data, geophysics=None, output=’geology’, compile_theano=False, theano_optimizer=’fast_compile’, u_grade=None, rescaling_factor=None, **kwargs)[source]

InterpolatorInput is a class that contains all the preprocessing operations to prepare the data to compute the model. Also is the object that has to be manipulated to vary the data without recompile the modeling function.

Parameters:
  • geo_data (gempy.data_management.InputData) – All values of a DataManagement object
  • geophysics (gempy.geophysics) – Object with the corresponding geophysical precomputations
  • compile_theano (bool) – select if the theano function is compiled during the initialization. Default: True
  • compute_all (bool) –

    If true the solution gives back the block model of lithologies, the potential field and the block model of faults. If False only return the block model of lithologies. This may be important to speed

    up the computation. Default True
  • u_grade (list) – grade of the polynomial for the universal part of the Kriging interpolations. The value has to
  • either 0, 3 or 9 (be) –
  • on the number of points given as input_data to try to avoid singular matrix. NOTE (depends) – if during the computation
  • the model a singular matrix is returned try to reduce the u_grade of the series. (of) –
  • rescaling_factor (float) – rescaling factor of the input_data data to improve the stability when float32 is used. By
  • the rescaling factor is calculated to obtein values between 0 and 1. (defaut) –
Keyword Arguments:
 
  • dtype ('str') – Choosing if using float32 or float64. This is important if is intended to use the GPU
  • Also InterpolatorClass kwargs (See) –
geo_data

Original gempy.DataManagement.InputData object

geo_data_res

Rescaled data. It has the same structure has gempy.InputData

interpolator

Instance of the gempy.DataManagement.InterpolaorInput.InterpolatorClass. See Also gempy.DataManagement.InterpolaorInput.InterpolatorClass docs th_fn: Theano function which compute the interpolation

dtype

type of float

class InterpolatorTheano(interp_data, **kwargs)[source]
Class which contain all needed methods to perform potential field implicit modelling in theano. Here there are methods to modify the shared parameters of the theano graph as well as the final preparation of the data from DataFrames to numpy arrays. This class is intended to be hidden from the user leaving the most useful calls into the InterpolatorData class
Parameters:
  • interp_data (gempy.data_management.InterpolatorData) – InterpolatorData: data rescaled plus geophysics and
  • additional data (other) –
Keyword Arguments:
 
  • range_var – Range of the variogram. Default None
  • c_o – Covariance at 0. Default None
  • nugget_effect – Nugget effect of the gradients. Default 0.01
  • u_grade – Grade of the polynomial used in the universal part of the Kriging. Default 2
  • rescaling_factor – Magic factor that multiplies the covariances). Default 2
  • verbose (list of str) – Level of verbosity during the execution of the functions. List of the strings with
  • parameters to be printed during the theano execution. TODO Make the list (the) –
data_prep(**kwargs)[source]

Ideally this method will only extract the data from the pandas dataframes to individual numpy arrays to be input_data of the theano function. However since some of the shared parameters are function of these arrays shape I also set them here

Returns:List of arrays which are the input_data for the theano function:
  • numpy.array: dips_position
  • numpy.array: dip_angles
  • numpy.array: azimuth
  • numpy.array: polarity
  • numpy.array: ref_layer_points
  • numpy.array: rest_layer_points
Return type:idl (list)
get_kriging_parameters(verbose=0)[source]

Print the kringing parameters

Parameters:verbose (int) – if > 0 print all the shape values as well.
Returns:None
order_table()[source]

First we sort the dataframes by the series age. Then we set a unique number for every formation and resort the formations. All inplace

set_densities(densities)[source]

WORKING IN PROGRESS – Set the weight of each voxel given a density :param densities:

Returns:

set_formation_number()[source]

Set a unique number to each formation. NOTE: this method is getting deprecated since the user does not need to know it and also now the numbers must be set in the order of the series as well. Therefore this method has been moved to the interpolator class as preprocessing

Returns:Column in the interfaces and orientations dataframes
set_theano_shared_parameteres(**kwargs)[source]

Here we create most of the kriging parameters. The user can pass them as kwargs otherwise we pick the default values from the DataManagement info. The share variables are set in place. All the parameters here are independent of the input_data data so this function only has to be called if you change the extent or grid or if you want to change one the kriging parameters.

Keyword Arguments:
 
  • u_grade (int) – Drift grade. Default to 2.
  • range_var (float) – Range of the variogram. Default 3D diagonal of the extent
  • c_o (float) – Covariance at lag 0. Default range_var ** 2 / 14 / 3. See my paper when I write it
  • nugget_effect (flaot) – Nugget effect of orientations. Default to 0.01
set_z_comp(tz, selected_cells)[source]

Set z component precomputation for the gravity. :param tz: :param selected_cells:

Returns:

compile_th_fn(output)[source]

Compile the theano function given the input_data data.

Parameters:compute_all (bool) –

If true the solution gives back the block model of lithologies, the potential field and the block model of faults. If False only return the block model of lithologies. This may be important to speed

up the computation. Default True
Returns:Compiled function if C or CUDA which computes the interpolation given the input_data data (XYZ of dips, dip, azimuth, polarity, XYZ ref interfaces, XYZ rest interfaces)
Return type:theano.function
get_input_data(u_grade=None)[source]

Get the theano variables that are input_data. This are necessary to compile the theano function or a theno op for pymc3 :param u_grade: grade of the polynomial for the universal part of the Kriging interpolations. The value has to :type u_grade: list

be either 0, 3 or 9 (number of equations) and the length has to be the number of series. By default the value depends on the number of points given as input_data to try to avoid singular matrix. NOTE: if during the computation of the model a singular matrix is returned try to reduce the u_grade of the series.
Returns:input_data nodes of the theano graph
Return type:theano.variables
rescale_data(geo_data, rescaling_factor=None)[source]

Rescale the data of a DataManagement object between 0 and 1 due to stability problem of the float32.

Parameters:
  • geo_data – Original gempy.DataManagement.InputData object
  • rescaling_factor (float) – factor of the rescaling. Default to maximum distance in one the axis
Returns:

Rescaled data

Return type:

gempy.data_management.InputData

set_geo_data_rescaled(geo_data, rescaling_factor=None)[source]

Set the rescale the data of a DataManagement object between 0 and 1 due to stability problem of the float32.

Parameters:
  • geo_data – Original gempy.DataManagement.InputData object
  • rescaling_factor (float) – factor of the rescaling. Default to maximum distance in one the axis
Returns:

Rescaled data

Return type:

gempy.data_management.InputData

set_gravity_precomputation(gravity_obj)[source]

Set a gravity object to the interpolator to pass the values to the theano graph

Parameters:gravity_obj (gempy.geophysics.GravityPreprocessing) – GravityPreprocessing (See Also gempy.geophysics.GravityPreprocessing documentation)
update_interpolator(geo_data=None, **kwargs)[source]

Method to update the constant parameters of the class interpolator (i.e. theano shared) WITHOUT recompiling. All the constant parameters for the interpolation can be passed as kwargs, otherwise they will take the default value (TODO: documentation of the dafault values)

Parameters:

geo_data – Rescaled gempy.DataManagement.InputData object. If None the stored geo_data_res will be used

Keyword Arguments:
 
  • range_var – Range of the variogram. Default None
  • c_o – Covariance at 0. Default None
  • nugget_effect – Nugget effect of the gradients. Default 0.01
  • u_grade – Grade of the polynomial used in the universal part of the Kriging. Default 2

DEP– I need to update this string Function that generates the symbolic code to perform the interpolation. Calling this function creates

both the theano functions for the potential field and the block.
returns:theano function for the potential field theano function for the block
class theano_graph.TheanoGraph(output=’geology’, optimizer=’fast_compile’, verbose=[0], dtype=’float32’)[source]

This class is used to help to divide the construction of the graph into sensical parts. All its methods build a part of the graph. Every method can be seen as a branch and collection of branches until the last method that will be the whole tree. Every part of the graph could be compiled separately but as we increase the complexity the input of each of these methods is more and more difficult to provide (if you are in a branch close to the trunk you need all the results of the branches above)

b_vector()[source]

Creation of the independent vector b to solve the kriging system

Parameters:verbose – -deprecated-
Returns:independent vector
Return type:theano.tensor.vector
block_series()[source]

Compute the part of the block model of a given series (dictated by the bool array yet to be computed)

Returns:Value of lithology at every interpolated point
Return type:theano.tensor.vector
compute_a_fault(len_i_0, len_i_1, len_f_0, len_f_1, n_form_per_serie_0, n_form_per_serie_1, u_grade_iter, final_block, fault_matrix)[source]

Function that loops each fault, generating a potential field for each on them with the respective block model

Parameters:
  • len_i_0 – Lenght of rest of previous series
  • len_i_1 – Lenght of rest for the computed series
  • len_f_0 – Lenght of dips of previous series
  • len_f_1 – Length of dips of the computed series
  • n_form_per_serie_0 – Number of formations of previous series
  • n_form_per_serie_1 – Number of formations of the computed series
Returns:

block model derived from the faults that afterwards is used as a drift for the “real” data

Return type:

theano.tensor.matrix

compute_a_series(len_i_0, len_i_1, len_f_0, len_f_1, n_form_per_serie_0, n_form_per_serie_1, u_grade_iter, final_block, fault_block)[source]

Function that loops each series, generating a potential field for each on them with the respective block model

Parameters:
  • len_i_0 – Lenght of rest of previous series
  • len_i_1 – Lenght of rest for the computed series
  • len_f_0 – Lenght of dips of previous series
  • len_f_1 – Length of dips of the computed series
  • n_form_per_serie_0 – Number of formations of previous series
  • n_form_per_serie_1 – Number of formations of the computed series
Returns:

final block model

Return type:

theano.tensor.matrix

cov_gradients(verbose=0)[source]

Create covariance function for the gradients

Returns:covariance of the gradients. Shape number of points in dip_pos x number of points in dip_pos
Return type:theano.tensor.matrix
cov_interface_gradients()[source]

Create covariance function for the gradiens :returns:

covariance of the gradients. Shape number of points in rest x number of
points in dip_pos
Return type:theano.tensor.matrix
cov_interfaces()[source]

Create covariance function for the interfaces

Returns:covariance of the interfaces. Shape number of points in rest x number of points in rest
Return type:theano.tensor.matrix
covariance_matrix()[source]

Set all the previous covariances together in the universal cokriging matrix

Returns:Multivariate covariance
Return type:theano.tensor.matrix
extend_dual_kriging()[source]

Tile the dual kriging vector to cover all the points to interpolate.So far I just make a matrix with the dimensions len(DK)x(grid) but in the future maybe I have to try to loop all this part so consume less memory

Returns:Matrix with the Dk parameters repeated for all the points to interpolate
Return type:theano.tensor.matrix
faults_contribution()[source]

Computation of the contribution of the faults drift at every point to interpolate. To get these we need to compute a whole block model with the faults data

Returns:Contribution of the faults drift (input) at every point to interpolate
Return type:theano.tensor.vector
faults_matrix()[source]

This function creates the part of the graph that generates the faults function creating a “block model” at the references and the rest of the points. Then this vector has to be appended to the covariance function

Returns:
  • theano.tensor.matrix: Drift matrix for the interfaces. Shape number of points in rest x n faults. This drif is a simple addition of an arbitrary number
  • theano.tensor.matrix: Drift matrix for the gradients. Shape number of points in dips x n faults. For discrete values this matrix will be null since the derivative of a constant is 0
Return type:list
gradient_contribution()[source]

Computation of the contribution of the foliations at every point to interpolate

Returns:Contribution of all foliations (input) at every point to interpolate
Return type:theano.tensor.vector
input_parameters_list()[source]

Create a list with the symbolic variables to use when we compile the theano function

Returns:
[self.dips_position_all, self.dip_angles_all, self.azimuth_all, self.polarity_all,
self.ref_layer_points_all, self.rest_layer_points_all]
Return type:list
interface_contribution()[source]

Computation of the contribution of the interfaces at every point to interpolate

Returns:Contribution of all interfaces (input) at every point to interpolate
Return type:theano.tensor.vector
matrices_shapes()[source]

Get all the lengths of the matrices that form the covariance matrix

Returns:length_of_CG, length_of_CGI, length_of_U_I, length_of_faults, length_of_C
scalar_field_at_all()[source]

Compute the potential field at all the interpolation points, i.e. grid plus rest plus ref :returns: Potential fields at all points :rtype: theano.tensor.vector

scalar_field_at_interfaces()[source]

Potential field at interfaces. To avoid errors I take all the points of rest that belong to one interface and make the average :returns: Potential field values at the interfaces of a given series :rtype: theano.tensor.vector

solve_kriging()[source]

Solve the kriging system. This has to get substituted by a more efficient and stable method QR decomposition in all likelihood

Returns:Dual kriging parameters
Return type:theano.tensor.vector
static squared_euclidean_distances(x_1, x_2)[source]

Compute the euclidian distances in 3D between all the points in x_1 and x_2

Parameters:
  • x_1 (theano.tensor.matrix) – shape n_points x number dimension
  • x_2 (theano.tensor.matrix) – shape n_points x number dimension
Returns:

Distancse matrix. shape n_points x n_points

Return type:

theano.tensor.matrix

universal_drift_contribution()[source]

Computation of the contribution of the universal drift at every point to interpolate

Returns:Contribution of the universal drift (input) at every point to interpolate
Return type:theano.tensor.vector
universal_matrix()[source]

Create the drift matrices for the potential field and its gradient

Returns:Drift matrix for the interfaces. Shape number of points in rest x 3**degree drift (except degree 0 that is 0)

theano.tensor.matrix: Drift matrix for the gradients. Shape number of points in dips x 3**degree drift (except degree 0 that is 0)

Return type:theano.tensor.matrix
x_to_interpolate(verbose=0)[source]

here I add to the grid points also the references points(to check the value of the potential field at the interfaces). Also here I will check what parts of the grid have been already computed in a previous series to avoid to recompute.

Returns:The 3D points of the given grid plus the reference and rest points
Return type:theano.tensor.matrix

This file is part of gempy.

gempy is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

gempy is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with gempy. If not, see <http://www.gnu.org/licenses/>.

@author: Alexander Schaaf

topology.check_adjacency(G, n1, n2)[source]

Check if n2 is adjacent/shares edge with n1.

topology.classify_edges(G, centroids, lith_block, fault_block)[source]

Classifies edges by type into stratigraphic or fault in “G.adj”. Accessible via G.adj[node1][node2][“edge_type”]

Parameters:
  • G (skimage.future.graph.rag.RAG) – Topology graph object.
  • centroids (dict) – Centroid dictionary {node id (int): tuple(x,y,z)}
  • lith_block (np.ndarray) – Shaped lithology block model.
  • fault_block (np.ndarray) – Shaped fault block model.

Returns:

topology.compare_graphs(G1, G2)[source]

Compares two graph objects using the Jaccard-Index.

Parameters:
  • G1 (skimage.future.graph.rag.RAG) – Topology graph object.
  • G2 (skimage.future.graph.rag.RAG) – Another topology graph object.
Returns:

Jaccard-Index

Return type:

(float)

topology.compute_adj_shape(n1, n2, labels_block)[source]

Compute shape of adjacency area: number of voxels in X, Y and Z (column height). (Fabian)

topology.compute_areas(G, labels_block)[source]

Computes adjacency areas and stores them in G.adj[n1][n2][“area”].

Parameters:
  • G (skimage.future.graph.rag.RAG) – Topology graph object.
  • labels_block (np.ndarray) – Uniquely labeled block model.
topology.get_centroids(label_block)[source]

Get node centroids in 2d and 3d as {node id (int): tuple(x,y,z)}.

topology.labels_lithology_lot(labels_unique, labels, block_original, verbose=0)[source]

Create LOT from label to lithology id.

topology.lithology_labels_lot(lithologies, labels, block_original, labels_unique, verbose=0)[source]

Create LOT from lithology id to label.

topology.topology_analyze(lith_block, fault_block, n_faults, areas_bool=False, return_block=False)[source]

Analyses the block models adjacency topology. Every lithological entity is described by a uniquely labeled node (centroid) and its connections to other entities by edges.

Parameters:
  • lith_block (np.ndarray) – Lithology block model
  • fault_block (np.ndarray) – Fault block model
  • n_faults (int) – Number of faults.
Keyword Arguments:
 
  • areas_bool (bool) – If True computes adjacency areas for connected nodes in voxel number. Default False.
  • return_block (bool) – If True additionally returns the uniquely labeled block model as np.ndarray.
Returns:

G: Region adjacency graph object (skimage.future.graph.rag.RAG) containing the adjacency topology graph

(G.adj).

centroids (dict): Centroid node coordinates as a dictionary with node id’s (int) as keys and (x,y,z) coordinates

as values.

labels_unique (np.array): List of all labels used. lith_to_labels_lot (dict): Dictionary look-up-table to go from lithology id to node id. labels_to_lith_lot (dict): Dictionary look-up-table to go from node id to lithology id.

Return type:

tuple

@author: Alexander Schaaf, Miguel de la Varga

posterior_analysis.calcualte_information_entropy(lith_prob)[source]

Calculates information entropy for the given probability array.

posterior_analysis.calculate_gradient(dip, az, pol)[source]

Calculates the gradient from dip, azimuth and polarity values.

posterior_analysis.calculate_information_entropy_total(ie, absolute=False)[source]

Calculate total information entropy (float) from an information entropy array.

posterior_analysis.change_input_data(db, interp_data, i, input_data_trace=’input_data’)[source]

Changes input data in interp_data to posterior input data at iteration i.

Parameters:
  • interp_data (gempy.data_management.InterpolationData) – An interp_data object with the structure we want to
  • compute.
  • i (int) – Iteration we want to recompute
Returns:

interp_data with the data of the given iteration

Return type:

gempy.data_management.InterpolationData

posterior_analysis.compute_posterior_models_all(db, interp_data, r=None, u_grade=None, get_potential_at_interfaces=False)[source]

Computes block models from stored input parameters for all iterations.

Parameters:
  • db
  • interp_data – GemPy interpolator object
  • r (tuple or list/array, optional) – The default value ‘None’ computes for models for the entire trace. If ‘tuple’ like (start, end) computes all models in the given range. If ‘np.array’ or ‘list’ with indices it computes the models for all given indices
  • u_grade (list, optional) –
  • get_potential_at_interfaces

Returns:

posterior_analysis.compute_probability_lithology(lith_blocks)[source]

Blocks must be just the lith blocks!

posterior_analysis.modify_plane_dip(dip, group_id, data_obj)[source]

Modify a dip angle of a plane identified by a group_id, recalculate the gradient and move the points vertically. Currently only supports the modification of dip angle - azimuth and polarity will stay the same.

Parameters:
  • dip (float) – Desired dip angle of the plane.
  • group_id (str) – Group id identifying the data points belonging to the plane.
  • ( (data_obj) – obj:): Data object to be modified (geo_data or interp_data.geo_data_res)
Returns:

Directly modifies the given data object.

posterior_analysis.move_plane_points(normal, centroid, data_obj, interf_f)[source]

Moves interface points to fit plane of given normal and centroid in data object.

Module with classes and methods to visualized structural geology data and potential fields of the regional modelling based on the potential field method. Tested on Ubuntu 14

Created on 23/09/2016

@author: Miguel de la Varga

class visualization.PlotData2D(geo_data, color_lot={0: [1.0, 0.6549019607843137, 0.14901960784313725], 1: [0.3607843137254902, 0.4196078431372549, 0.7529411764705882], 2: [0.9372549019607843, 0.3254901960784314, 0.3137254901960784], 3: [1.0, 0.9333333333333333, 0.34509803921568627], 4: [0.4, 0.7333333333333333, 0.41568627450980394], 5: [0.5529411764705883, 0.43137254901960786, 0.38823529411764707], 6: [0.25882352941176473, 0.6470588235294118, 0.9607843137254902], 7: [1.0, 0.792156862745098, 0.1568627450980392], 8: [0.9254901960784314, 0.25098039215686274, 0.47843137254901963], 9: [0.1607843137254902, 0.7137254901960784, 0.9647058823529412], 10: [0.8313725490196079, 0.8823529411764706, 0.3411764705882353], 11: [0.47058823529411764, 0.5647058823529412, 0.611764705882353], 12: [1.0, 0.4392156862745098, 0.2627450980392157], 13: [0.7411764705882353, 0.7411764705882353, 0.7411764705882353], 14: [0.14901960784313725, 0.7764705882352941, 0.8549019607843137], 15: [0.49411764705882355, 0.3411764705882353, 0.7607843137254902], 16: [0.6705882352941176, 0.2784313725490196, 0.7372549019607844], 17: [0.14901960784313725, 0.6509803921568628, 0.6039215686274509], 18: [0.611764705882353, 0.8, 0.396078431372549], 19: [0.9607843137254902, 0.48627450980392156, 0.0], 20: [0.18823529411764706, 0.24705882352941178, 0.6235294117647059], 21: [0.8274509803921568, 0.1843137254901961, 0.1843137254901961], 22: [0.984313725490196, 0.7529411764705882, 0.17647058823529413], 23: [0.2196078431372549, 0.5568627450980392, 0.23529411764705882], 24: [0.36470588235294116, 0.25098039215686274, 0.21568627450980393], 25: [0.09803921568627451, 0.4627450980392157, 0.8235294117647058], 26: [1.0, 0.6274509803921569, 0.0], 27: [0.7607843137254902, 0.09411764705882353, 0.3568627450980392], 28: [0.00784313725490196, 0.5333333333333333, 0.8196078431372549], 29: [0.6862745098039216, 0.7058823529411765, 0.16862745098039217], 30: [0.27058823529411763, 0.35294117647058826, 0.39215686274509803], 31: [0.9019607843137255, 0.2901960784313726, 0.09803921568627451], 32: [0.3803921568627451, 0.3803921568627451, 0.3803921568627451], 33: [0.0, 0.592156862745098, 0.6549019607843137], 34: [0.3176470588235294, 0.17647058823529413, 0.6588235294117647], 35: [0.4823529411764706, 0.12156862745098039, 0.6352941176470588], 36: [0.0, 0.4745098039215686, 0.4196078431372549], 37: [0.40784313725490196, 0.6235294117647059, 0.2196078431372549], 38: [1.0, 0.8784313725490196, 0.6980392156862745], 39: [0.7725490196078432, 0.792156862745098, 0.9137254901960784], 40: [1.0, 0.803921568627451, 0.8235294117647058], 41: [1.0, 0.9764705882352941, 0.7686274509803922], 42: [0.7843137254901961, 0.9019607843137255, 0.788235294117647], 43: [0.8431372549019608, 0.8, 0.7843137254901961], 44: [0.7333333333333333, 0.8705882352941177, 0.984313725490196], 45: [1.0, 0.9254901960784314, 0.7019607843137254], 46: [0.9725490196078431, 0.7333333333333333, 0.8156862745098039], 47: [0.7019607843137254, 0.8980392156862745, 0.9882352941176471], 48: [0.9411764705882353, 0.9568627450980393, 0.7647058823529411], 49: [0.8117647058823529, 0.8470588235294118, 0.8627450980392157], 50: [1.0, 0.8, 0.7372549019607844], 51: [0.9607843137254902, 0.9607843137254902, 0.9607843137254902], 52: [0.6980392156862745, 0.9215686274509803, 0.9490196078431372], 53: [0.8196078431372549, 0.7686274509803922, 0.9137254901960784], 54: [0.8823529411764706, 0.7450980392156863, 0.9058823529411765], 55: [0.6980392156862745, 0.8745098039215686, 0.8588235294117647], 56: [0.8627450980392157, 0.9294117647058824, 0.7843137254901961], 57: [1.0, 0.8, 0.5019607843137255], 58: [0.6235294117647059, 0.6588235294117647, 0.8549019607843137], 59: [0.9372549019607843, 0.6039215686274509, 0.6039215686274509], 60: [1.0, 0.9607843137254902, 0.615686274509804], 61: [0.6470588235294118, 0.8392156862745098, 0.6549019607843137], 62: [0.7372549019607844, 0.6666666666666666, 0.6431372549019608], 63: [0.5647058823529412, 0.792156862745098, 0.9764705882352941], 64: [1.0, 0.8784313725490196, 0.5098039215686274], 65: [0.9568627450980393, 0.5607843137254902, 0.6941176470588235], 66: [0.5058823529411764, 0.8313725490196079, 0.9803921568627451], 67: [0.9019607843137255, 0.9333333333333333, 0.611764705882353], 68: [0.6901960784313725, 0.7450980392156863, 0.7725490196078432], 69: [1.0, 0.6705882352941176, 0.5686274509803921], 70: [0.9333333333333333, 0.9333333333333333, 0.9333333333333333], 71: [0.5019607843137255, 0.8705882352941177, 0.9176470588235294], 72: [0.7019607843137254, 0.615686274509804, 0.8588235294117647], 73: [0.807843137254902, 0.5764705882352941, 0.8470588235294118], 74: [0.5019607843137255, 0.796078431372549, 0.7686274509803922], 75: [0.7725490196078432, 0.8823529411764706, 0.6470588235294118], 76: [1.0, 0.7176470588235294, 0.30196078431372547], 77: [0.4745098039215686, 0.5254901960784314, 0.796078431372549], 78: [0.8980392156862745, 0.45098039215686275, 0.45098039215686275], 79: [1.0, 0.9450980392156862, 0.4627450980392157], 80: [0.5058823529411764, 0.7803921568627451, 0.5176470588235295], 81: [0.6313725490196078, 0.5333333333333333, 0.4980392156862745], 82: [0.39215686274509803, 0.7098039215686275, 0.9647058823529412], 83: [1.0, 0.8352941176470589, 0.30980392156862746], 84: [0.9411764705882353, 0.3843137254901961, 0.5725490196078431], 85: [0.30980392156862746, 0.7647058823529411, 0.9686274509803922], 86: [0.8627450980392157, 0.9058823529411765, 0.4588235294117647], 87: [0.5647058823529412, 0.6431372549019608, 0.6823529411764706], 88: [1.0, 0.5411764705882353, 0.396078431372549], 89: [0.8784313725490196, 0.8784313725490196, 0.8784313725490196], 90: [0.30196078431372547, 0.8156862745098039, 0.8823529411764706], 91: [0.5843137254901961, 0.4588235294117647, 0.803921568627451], 92: [0.7294117647058823, 0.40784313725490196, 0.7843137254901961], 93: [0.30196078431372547, 0.7137254901960784, 0.6745098039215687], 94: [0.6823529411764706, 0.8352941176470589, 0.5058823529411764], 95: [1.0, 0.596078431372549, 0.0], 96: [0.24705882352941178, 0.3176470588235294, 0.7098039215686275], 97: [0.9568627450980393, 0.2627450980392157, 0.21176470588235294], 98: [1.0, 0.9215686274509803, 0.23137254901960785], 99: [0.2980392156862745, 0.6862745098039216, 0.3137254901960784], 100: [0.4745098039215686, 0.3333333333333333, 0.2823529411764706], 101: [0.12941176470588237, 0.5882352941176471, 0.9529411764705882], 102: [1.0, 0.7568627450980392, 0.027450980392156862], 103: [0.9137254901960784, 0.11764705882352941, 0.38823529411764707], 104: [0.011764705882352941, 0.6627450980392157, 0.9568627450980393], 105: [0.803921568627451, 0.8627450980392157, 0.2235294117647059], 106: [0.3764705882352941, 0.49019607843137253, 0.5450980392156862], 107: [1.0, 0.3411764705882353, 0.13333333333333333], 108: [0.6196078431372549, 0.6196078431372549, 0.6196078431372549], 109: [0.0, 0.7372549019607844, 0.8313725490196079], 110: [0.403921568627451, 0.22745098039215686, 0.7176470588235294], 111: [0.611764705882353, 0.15294117647058825, 0.6901960784313725], 112: [0.0, 0.5882352941176471, 0.5333333333333333], 113: [0.5450980392156862, 0.7647058823529411, 0.2901960784313726], 114: [0.984313725490196, 0.5490196078431373, 0.0], 115: [0.2235294117647059, 0.28627450980392155, 0.6705882352941176], 116: [0.8980392156862745, 0.2235294117647059, 0.20784313725490197], 117: [0.9921568627450981, 0.8470588235294118, 0.20784313725490197], 118: [0.2627450980392157, 0.6274509803921569, 0.2784313725490196], 119: [0.42745098039215684, 0.2980392156862745, 0.2549019607843137], 120: [0.11764705882352941, 0.5333333333333333, 0.8980392156862745], 121: [1.0, 0.7019607843137254, 0.0], 122: [0.8470588235294118, 0.10588235294117647, 0.3764705882352941], 123: [0.011764705882352941, 0.6078431372549019, 0.8980392156862745], 124: [0.7529411764705882, 0.792156862745098, 0.2], 125: [0.32941176470588235, 0.43137254901960786, 0.47843137254901963], 126: [0.9568627450980393, 0.3176470588235294, 0.11764705882352941], 127: [0.4588235294117647, 0.4588235294117647, 0.4588235294117647], 128: [0.0, 0.6745098039215687, 0.7568627450980392], 129: [0.3686274509803922, 0.20784313725490197, 0.6941176470588235], 130: [0.5568627450980392, 0.1411764705882353, 0.6666666666666666], 131: [0.0, 0.5372549019607843, 0.4823529411764706], 132: [0.48627450980392156, 0.7019607843137254, 0.25882352941176473], 133: [1.0, 0.9529411764705882, 0.8784313725490196], 134: [0.9098039215686274, 0.9176470588235294, 0.9647058823529412], 135: [1.0, 0.9215686274509803, 0.9333333333333333], 136: [1.0, 0.9921568627450981, 0.9058823529411765], 137: [0.9098039215686274, 0.9607843137254902, 0.9137254901960784], 138: [0.9372549019607843, 0.9215686274509803, 0.9137254901960784], 139: [0.8901960784313725, 0.9490196078431372, 0.9921568627450981], 140: [1.0, 0.9725490196078431, 0.8823529411764706], 141: [0.9882352941176471, 0.8941176470588236, 0.9254901960784314], 142: [0.8823529411764706, 0.9607843137254902, 0.996078431372549], 143: [0.9764705882352941, 0.984313725490196, 0.9058823529411765], 144: [0.9254901960784314, 0.9372549019607843, 0.9450980392156862], 145: [0.984313725490196, 0.9137254901960784, 0.9058823529411765], 146: [0.9803921568627451, 0.9803921568627451, 0.9803921568627451], 147: [0.8784313725490196, 0.9686274509803922, 0.9803921568627451], 148: [0.9294117647058824, 0.9058823529411765, 0.9647058823529412], 149: [0.9529411764705882, 0.8980392156862745, 0.9607843137254902], 150: [0.8784313725490196, 0.9490196078431372, 0.9450980392156862], 151: [0.9450980392156862, 0.9725490196078431, 0.9137254901960784], 152: [0.9372549019607843, 0.4235294117647059, 0.0], 153: [0.1568627450980392, 0.20784313725490197, 0.5764705882352941], 154: [0.7764705882352941, 0.1568627450980392, 0.1568627450980392], 155: [0.9764705882352941, 0.6588235294117647, 0.1450980392156863], 156: [0.1803921568627451, 0.49019607843137253, 0.19607843137254902], 157: [0.3058823529411765, 0.20392156862745098, 0.1803921568627451], 158: [0.08235294117647059, 0.396078431372549, 0.7529411764705882], 159: [1.0, 0.5607843137254902, 0.0], 160: [0.6784313725490196, 0.0784313725490196, 0.3411764705882353], 161: [0.00784313725490196, 0.4666666666666667, 0.7411764705882353], 162: [0.6196078431372549, 0.615686274509804, 0.1411764705882353], 163: [0.21568627450980393, 0.2784313725490196, 0.30980392156862746], 164: [0.8470588235294118, 0.2627450980392157, 0.08235294117647059], 165: [0.25882352941176473, 0.25882352941176473, 0.25882352941176473], 166: [0.0, 0.5137254901960784, 0.5607843137254902], 167: [0.27058823529411763, 0.15294117647058825, 0.6274509803921569], 168: [0.41568627450980394, 0.10588235294117647, 0.6039215686274509], 169: [0.0, 0.4117647058823529, 0.3607843137254902], 170: [0.3333333333333333, 0.5450980392156862, 0.1843137254901961], 171: [0.9019607843137255, 0.3176470588235294, 0.0], 172: [0.10196078431372549, 0.13725490196078433, 0.49411764705882355], 173: [0.7176470588235294, 0.10980392156862745, 0.10980392156862745], 174: [0.9607843137254902, 0.4980392156862745, 0.09019607843137255], 175: [0.10588235294117647, 0.3686274509803922, 0.12549019607843137], 176: [0.24313725490196078, 0.15294117647058825, 0.13725490196078433], 177: [0.050980392156862744, 0.2784313725490196, 0.6313725490196078], 178: [1.0, 0.43529411764705883, 0.0], 179: [0.5333333333333333, 0.054901960784313725, 0.30980392156862746], 180: [0.00392156862745098, 0.3411764705882353, 0.6078431372549019], 181: [0.5098039215686274, 0.4666666666666667, 0.09019607843137255], 182: [0.14901960784313725, 0.19607843137254902, 0.2196078431372549], 183: [0.7490196078431373, 0.21176470588235294, 0.047058823529411764], 184: [0.12941176470588237, 0.12941176470588237, 0.12941176470588237], 185: [0.0, 0.3764705882352941, 0.39215686274509803], 186: [0.19215686274509805, 0.10588235294117647, 0.5725490196078431], 187: [0.2901960784313726, 0.0784313725490196, 0.5490196078431373], 188: [0.0, 0.30196078431372547, 0.25098039215686274], 189: [0.2, 0.4117647058823529, 0.11764705882352941]}, cmap=<matplotlib.colors.ListedColormap object>, norm=<matplotlib.colors.BoundaryNorm object>, **kwargs)[source]

Class to make the different plot related with gempy

Parameters:
  • geo_data (gempy.InputData) – All values of a DataManagement object
  • block (numpy.array) – 3D array containing the lithology block
  • **kwargs – Arbitrary keyword arguments.
Keyword Arguments:
 
  • scalar_field (numpy.ndarray) – 3D array containing a individual potential field
  • verbose (int) – Level of verbosity during the execution of the functions (up to 5). Default 0
static annotate_plot(frame, label_col, x, y, **kwargs)[source]

Annotate the plot of a given DataFrame using one of its columns

Should be called right after a DataFrame or series plot method, before telling matplotlib to show the plot.

Parameters:
  • frame (pandas.DataFrame) –
  • plot_col (str) – The string identifying the column of frame that was plotted
  • label_col (str) – The string identifying the column of frame to be used as label
  • kwargs – Other key-word args that should be passed to plt.annotate
Returns:

Return type:

None

Notes

After calling this function you should call plt.show() to get the results. This function only adds the annotations, it doesn’t show them.

plot_block_section(cell_number=13, block=None, direction=’y’, interpolation=’none’, plot_data=False, **kwargs)[source]

Plot a section of the block model

Parameters:
  • cell_number (int) – position of the array to plot
  • direction (str) – xyz. Caartesian direction to be plotted
  • interpolation (str) – Type of interpolation of plt.imshow. Default ‘none’. Acceptable values are ‘none’
  • 'bilinear', 'bicubic', (,'nearest',) –
  • 'spline36', 'hanning', 'hamming', 'hermite', 'kaiser', ('spline16',) –
  • 'catrom', 'gaussian', 'bessel', 'mitchell', 'sinc', ('quadric',) –
  • 'lanczos'
  • **kwargs – imshow keywargs
Returns:

Block plot

plot_data(direction=’y’, data_type=’all’, series=’all’, legend_font_size=10, **kwargs)[source]

Plot the projecton of the raw data (interfaces and orientations) in 2D following a specific directions

Parameters:
  • direction (str) – xyz. Caartesian direction to be plotted
  • data_type (str) – type of data to plot. ‘all’, ‘interfaces’ or ‘orientations’
  • series (str) – series to plot
  • **kwargs – seaborn lmplot key arguments. (TODO: adding the link to them)
Returns:

Data plot

plot_scalar_field(scalar_field, cell_number, N=20, direction=’y’, plot_data=True, series=’all’, *args, **kwargs)[source]

Plot a scalar field in a given direction.

Parameters:
  • cell_number (int) – position of the array to plot
  • scalar_field (str) – name of the scalar field (or series) to plot
  • n_pf (int) – number of the scalar field (or series) to plot
  • direction (str) – xyz. Caartesian direction to be plotted
  • serieDeprecated
  • **kwargs – plt.contour kwargs
Returns:

scalar field plot

class visualization.vtkVisualization(geo_data, ren_name=’GemPy 3D-Editor’, verbose=0, color_lot={0: [1.0, 0.6549019607843137, 0.14901960784313725], 1: [0.3607843137254902, 0.4196078431372549, 0.7529411764705882], 2: [0.9372549019607843, 0.3254901960784314, 0.3137254901960784], 3: [1.0, 0.9333333333333333, 0.34509803921568627], 4: [0.4, 0.7333333333333333, 0.41568627450980394], 5: [0.5529411764705883, 0.43137254901960786, 0.38823529411764707], 6: [0.25882352941176473, 0.6470588235294118, 0.9607843137254902], 7: [1.0, 0.792156862745098, 0.1568627450980392], 8: [0.9254901960784314, 0.25098039215686274, 0.47843137254901963], 9: [0.1607843137254902, 0.7137254901960784, 0.9647058823529412], 10: [0.8313725490196079, 0.8823529411764706, 0.3411764705882353], 11: [0.47058823529411764, 0.5647058823529412, 0.611764705882353], 12: [1.0, 0.4392156862745098, 0.2627450980392157], 13: [0.7411764705882353, 0.7411764705882353, 0.7411764705882353], 14: [0.14901960784313725, 0.7764705882352941, 0.8549019607843137], 15: [0.49411764705882355, 0.3411764705882353, 0.7607843137254902], 16: [0.6705882352941176, 0.2784313725490196, 0.7372549019607844], 17: [0.14901960784313725, 0.6509803921568628, 0.6039215686274509], 18: [0.611764705882353, 0.8, 0.396078431372549], 19: [0.9607843137254902, 0.48627450980392156, 0.0], 20: [0.18823529411764706, 0.24705882352941178, 0.6235294117647059], 21: [0.8274509803921568, 0.1843137254901961, 0.1843137254901961], 22: [0.984313725490196, 0.7529411764705882, 0.17647058823529413], 23: [0.2196078431372549, 0.5568627450980392, 0.23529411764705882], 24: [0.36470588235294116, 0.25098039215686274, 0.21568627450980393], 25: [0.09803921568627451, 0.4627450980392157, 0.8235294117647058], 26: [1.0, 0.6274509803921569, 0.0], 27: [0.7607843137254902, 0.09411764705882353, 0.3568627450980392], 28: [0.00784313725490196, 0.5333333333333333, 0.8196078431372549], 29: [0.6862745098039216, 0.7058823529411765, 0.16862745098039217], 30: [0.27058823529411763, 0.35294117647058826, 0.39215686274509803], 31: [0.9019607843137255, 0.2901960784313726, 0.09803921568627451], 32: [0.3803921568627451, 0.3803921568627451, 0.3803921568627451], 33: [0.0, 0.592156862745098, 0.6549019607843137], 34: [0.3176470588235294, 0.17647058823529413, 0.6588235294117647], 35: [0.4823529411764706, 0.12156862745098039, 0.6352941176470588], 36: [0.0, 0.4745098039215686, 0.4196078431372549], 37: [0.40784313725490196, 0.6235294117647059, 0.2196078431372549], 38: [1.0, 0.8784313725490196, 0.6980392156862745], 39: [0.7725490196078432, 0.792156862745098, 0.9137254901960784], 40: [1.0, 0.803921568627451, 0.8235294117647058], 41: [1.0, 0.9764705882352941, 0.7686274509803922], 42: [0.7843137254901961, 0.9019607843137255, 0.788235294117647], 43: [0.8431372549019608, 0.8, 0.7843137254901961], 44: [0.7333333333333333, 0.8705882352941177, 0.984313725490196], 45: [1.0, 0.9254901960784314, 0.7019607843137254], 46: [0.9725490196078431, 0.7333333333333333, 0.8156862745098039], 47: [0.7019607843137254, 0.8980392156862745, 0.9882352941176471], 48: [0.9411764705882353, 0.9568627450980393, 0.7647058823529411], 49: [0.8117647058823529, 0.8470588235294118, 0.8627450980392157], 50: [1.0, 0.8, 0.7372549019607844], 51: [0.9607843137254902, 0.9607843137254902, 0.9607843137254902], 52: [0.6980392156862745, 0.9215686274509803, 0.9490196078431372], 53: [0.8196078431372549, 0.7686274509803922, 0.9137254901960784], 54: [0.8823529411764706, 0.7450980392156863, 0.9058823529411765], 55: [0.6980392156862745, 0.8745098039215686, 0.8588235294117647], 56: [0.8627450980392157, 0.9294117647058824, 0.7843137254901961], 57: [1.0, 0.8, 0.5019607843137255], 58: [0.6235294117647059, 0.6588235294117647, 0.8549019607843137], 59: [0.9372549019607843, 0.6039215686274509, 0.6039215686274509], 60: [1.0, 0.9607843137254902, 0.615686274509804], 61: [0.6470588235294118, 0.8392156862745098, 0.6549019607843137], 62: [0.7372549019607844, 0.6666666666666666, 0.6431372549019608], 63: [0.5647058823529412, 0.792156862745098, 0.9764705882352941], 64: [1.0, 0.8784313725490196, 0.5098039215686274], 65: [0.9568627450980393, 0.5607843137254902, 0.6941176470588235], 66: [0.5058823529411764, 0.8313725490196079, 0.9803921568627451], 67: [0.9019607843137255, 0.9333333333333333, 0.611764705882353], 68: [0.6901960784313725, 0.7450980392156863, 0.7725490196078432], 69: [1.0, 0.6705882352941176, 0.5686274509803921], 70: [0.9333333333333333, 0.9333333333333333, 0.9333333333333333], 71: [0.5019607843137255, 0.8705882352941177, 0.9176470588235294], 72: [0.7019607843137254, 0.615686274509804, 0.8588235294117647], 73: [0.807843137254902, 0.5764705882352941, 0.8470588235294118], 74: [0.5019607843137255, 0.796078431372549, 0.7686274509803922], 75: [0.7725490196078432, 0.8823529411764706, 0.6470588235294118], 76: [1.0, 0.7176470588235294, 0.30196078431372547], 77: [0.4745098039215686, 0.5254901960784314, 0.796078431372549], 78: [0.8980392156862745, 0.45098039215686275, 0.45098039215686275], 79: [1.0, 0.9450980392156862, 0.4627450980392157], 80: [0.5058823529411764, 0.7803921568627451, 0.5176470588235295], 81: [0.6313725490196078, 0.5333333333333333, 0.4980392156862745], 82: [0.39215686274509803, 0.7098039215686275, 0.9647058823529412], 83: [1.0, 0.8352941176470589, 0.30980392156862746], 84: [0.9411764705882353, 0.3843137254901961, 0.5725490196078431], 85: [0.30980392156862746, 0.7647058823529411, 0.9686274509803922], 86: [0.8627450980392157, 0.9058823529411765, 0.4588235294117647], 87: [0.5647058823529412, 0.6431372549019608, 0.6823529411764706], 88: [1.0, 0.5411764705882353, 0.396078431372549], 89: [0.8784313725490196, 0.8784313725490196, 0.8784313725490196], 90: [0.30196078431372547, 0.8156862745098039, 0.8823529411764706], 91: [0.5843137254901961, 0.4588235294117647, 0.803921568627451], 92: [0.7294117647058823, 0.40784313725490196, 0.7843137254901961], 93: [0.30196078431372547, 0.7137254901960784, 0.6745098039215687], 94: [0.6823529411764706, 0.8352941176470589, 0.5058823529411764], 95: [1.0, 0.596078431372549, 0.0], 96: [0.24705882352941178, 0.3176470588235294, 0.7098039215686275], 97: [0.9568627450980393, 0.2627450980392157, 0.21176470588235294], 98: [1.0, 0.9215686274509803, 0.23137254901960785], 99: [0.2980392156862745, 0.6862745098039216, 0.3137254901960784], 100: [0.4745098039215686, 0.3333333333333333, 0.2823529411764706], 101: [0.12941176470588237, 0.5882352941176471, 0.9529411764705882], 102: [1.0, 0.7568627450980392, 0.027450980392156862], 103: [0.9137254901960784, 0.11764705882352941, 0.38823529411764707], 104: [0.011764705882352941, 0.6627450980392157, 0.9568627450980393], 105: [0.803921568627451, 0.8627450980392157, 0.2235294117647059], 106: [0.3764705882352941, 0.49019607843137253, 0.5450980392156862], 107: [1.0, 0.3411764705882353, 0.13333333333333333], 108: [0.6196078431372549, 0.6196078431372549, 0.6196078431372549], 109: [0.0, 0.7372549019607844, 0.8313725490196079], 110: [0.403921568627451, 0.22745098039215686, 0.7176470588235294], 111: [0.611764705882353, 0.15294117647058825, 0.6901960784313725], 112: [0.0, 0.5882352941176471, 0.5333333333333333], 113: [0.5450980392156862, 0.7647058823529411, 0.2901960784313726], 114: [0.984313725490196, 0.5490196078431373, 0.0], 115: [0.2235294117647059, 0.28627450980392155, 0.6705882352941176], 116: [0.8980392156862745, 0.2235294117647059, 0.20784313725490197], 117: [0.9921568627450981, 0.8470588235294118, 0.20784313725490197], 118: [0.2627450980392157, 0.6274509803921569, 0.2784313725490196], 119: [0.42745098039215684, 0.2980392156862745, 0.2549019607843137], 120: [0.11764705882352941, 0.5333333333333333, 0.8980392156862745], 121: [1.0, 0.7019607843137254, 0.0], 122: [0.8470588235294118, 0.10588235294117647, 0.3764705882352941], 123: [0.011764705882352941, 0.6078431372549019, 0.8980392156862745], 124: [0.7529411764705882, 0.792156862745098, 0.2], 125: [0.32941176470588235, 0.43137254901960786, 0.47843137254901963], 126: [0.9568627450980393, 0.3176470588235294, 0.11764705882352941], 127: [0.4588235294117647, 0.4588235294117647, 0.4588235294117647], 128: [0.0, 0.6745098039215687, 0.7568627450980392], 129: [0.3686274509803922, 0.20784313725490197, 0.6941176470588235], 130: [0.5568627450980392, 0.1411764705882353, 0.6666666666666666], 131: [0.0, 0.5372549019607843, 0.4823529411764706], 132: [0.48627450980392156, 0.7019607843137254, 0.25882352941176473], 133: [1.0, 0.9529411764705882, 0.8784313725490196], 134: [0.9098039215686274, 0.9176470588235294, 0.9647058823529412], 135: [1.0, 0.9215686274509803, 0.9333333333333333], 136: [1.0, 0.9921568627450981, 0.9058823529411765], 137: [0.9098039215686274, 0.9607843137254902, 0.9137254901960784], 138: [0.9372549019607843, 0.9215686274509803, 0.9137254901960784], 139: [0.8901960784313725, 0.9490196078431372, 0.9921568627450981], 140: [1.0, 0.9725490196078431, 0.8823529411764706], 141: [0.9882352941176471, 0.8941176470588236, 0.9254901960784314], 142: [0.8823529411764706, 0.9607843137254902, 0.996078431372549], 143: [0.9764705882352941, 0.984313725490196, 0.9058823529411765], 144: [0.9254901960784314, 0.9372549019607843, 0.9450980392156862], 145: [0.984313725490196, 0.9137254901960784, 0.9058823529411765], 146: [0.9803921568627451, 0.9803921568627451, 0.9803921568627451], 147: [0.8784313725490196, 0.9686274509803922, 0.9803921568627451], 148: [0.9294117647058824, 0.9058823529411765, 0.9647058823529412], 149: [0.9529411764705882, 0.8980392156862745, 0.9607843137254902], 150: [0.8784313725490196, 0.9490196078431372, 0.9450980392156862], 151: [0.9450980392156862, 0.9725490196078431, 0.9137254901960784], 152: [0.9372549019607843, 0.4235294117647059, 0.0], 153: [0.1568627450980392, 0.20784313725490197, 0.5764705882352941], 154: [0.7764705882352941, 0.1568627450980392, 0.1568627450980392], 155: [0.9764705882352941, 0.6588235294117647, 0.1450980392156863], 156: [0.1803921568627451, 0.49019607843137253, 0.19607843137254902], 157: [0.3058823529411765, 0.20392156862745098, 0.1803921568627451], 158: [0.08235294117647059, 0.396078431372549, 0.7529411764705882], 159: [1.0, 0.5607843137254902, 0.0], 160: [0.6784313725490196, 0.0784313725490196, 0.3411764705882353], 161: [0.00784313725490196, 0.4666666666666667, 0.7411764705882353], 162: [0.6196078431372549, 0.615686274509804, 0.1411764705882353], 163: [0.21568627450980393, 0.2784313725490196, 0.30980392156862746], 164: [0.8470588235294118, 0.2627450980392157, 0.08235294117647059], 165: [0.25882352941176473, 0.25882352941176473, 0.25882352941176473], 166: [0.0, 0.5137254901960784, 0.5607843137254902], 167: [0.27058823529411763, 0.15294117647058825, 0.6274509803921569], 168: [0.41568627450980394, 0.10588235294117647, 0.6039215686274509], 169: [0.0, 0.4117647058823529, 0.3607843137254902], 170: [0.3333333333333333, 0.5450980392156862, 0.1843137254901961], 171: [0.9019607843137255, 0.3176470588235294, 0.0], 172: [0.10196078431372549, 0.13725490196078433, 0.49411764705882355], 173: [0.7176470588235294, 0.10980392156862745, 0.10980392156862745], 174: [0.9607843137254902, 0.4980392156862745, 0.09019607843137255], 175: [0.10588235294117647, 0.3686274509803922, 0.12549019607843137], 176: [0.24313725490196078, 0.15294117647058825, 0.13725490196078433], 177: [0.050980392156862744, 0.2784313725490196, 0.6313725490196078], 178: [1.0, 0.43529411764705883, 0.0], 179: [0.5333333333333333, 0.054901960784313725, 0.30980392156862746], 180: [0.00392156862745098, 0.3411764705882353, 0.6078431372549019], 181: [0.5098039215686274, 0.4666666666666667, 0.09019607843137255], 182: [0.14901960784313725, 0.19607843137254902, 0.2196078431372549], 183: [0.7490196078431373, 0.21176470588235294, 0.047058823529411764], 184: [0.12941176470588237, 0.12941176470588237, 0.12941176470588237], 185: [0.0, 0.3764705882352941, 0.39215686274509803], 186: [0.19215686274509805, 0.10588235294117647, 0.5725490196078431], 187: [0.2901960784313726, 0.0784313725490196, 0.5490196078431373], 188: [0.0, 0.30196078431372547, 0.25098039215686274], 189: [0.2, 0.4117647058823529, 0.11764705882352941]}, real_time=False, bg_color=None)[source]

Class to visualize data and results in 3D. Init will create all the render properties while the method render model will lunch the window. Using set_interfaces, set_orientations and set_surfaces in between can be chosen what will be displayed.

Parameters:
  • geo_data (gempy.InputData) – All values of a DataManagement object
  • ren_name (str) – Name of the renderer window
  • verbose (int) – Verbosity for certain functions
renWin

vtk.vtkRenderWindow()

camera_list

list – list of cameras for the distinct renderers

ren_list

list – list containing the vtk renderers

create_foliation(X, Y, Z, fn, Gx, Gy, Gz, n_plane=0, n_render=0, n_index=0, alpha=0.5)[source]

Method to create a plane given a foliation

Parameters:
  • X – X coord
  • Y – Y coord
  • Z – Z corrd
  • fn (int) – Formation number
  • Gx (str) – Component of the gradient x
  • Gy (str) – Component of the gradient y
  • Gz (str) – Component of the gradient z
  • n_plane (int) – Number of the plane
  • n_render (int) – Number of the render where the plane belongs
  • n_index (int) – index value in the PandasDataframe of InupData.interfaces
  • alpha – Opacity of the plane
Returns:

vtk.vtkPlaneWidget

create_ren_list()[source]

Create a list of the 4 renderers we use. One general view and 3 cartersian projections :returns: list of renderers :rtype: list

create_sphere(X, Y, Z, fn, n_sphere=0, n_render=0, n_index=0, r=0.03)[source]

Method to create the sphere that represent the interfaces points :param X: X coord :param Y: Y coord :param Z: Z corrd :param fn: Formation number :type fn: int :param n_sphere: Number of the sphere :type n_sphere: int :param n_render: Number of the render where the sphere belongs :type n_render: int :param n_index: index value in the PandasDataframe of InupData.interfaces :type n_index: int :param r: radio of the sphere :type r: float

Returns:vtk.vtkSphereWidget
create_surface(vertices, simplices, fn, alpha=1)[source]

Method to create the polydata that define the surfaces

Parameters:
  • vertices (numpy.array) – 2D array (XYZ) with the coordinates of the points
  • simplices (numpy.array) – 2D array with the value of the vertices that form every single triangle
  • fn (int) – Formation number
  • alpha (float) – Opacity
Returns:

vtk.vtkActor, vtk.vtkPolyDataMapper, vtk.vtkPolyData

static create_surface_points(vertices)[source]

Method to create the points that form the surfaces :param vertices: 2D array (XYZ) with the coordinates of the points :type vertices: numpy.array

Returns:with the coordinates of the points
Return type:vtk.vtkPoints
static create_surface_triangles(simplices)[source]

Method to create the Triangles that form the surfaces :param simplices: 2D array with the value of the vertices that form every single triangle :type simplices: numpy.array

Returns:vtk.vtkTriangle
static export_vtk_lith_block(geo_data, lith_block, path=None)[source]

Export data to a vtk file for posterior visualizations

Parameters:
  • geo_data (gempy.InputData) – All values of a DataManagement object
  • block (numpy.array) – 3D array containing the lithology block
  • path (str) – path to the location of the vtk
Returns:

None

static export_vtk_surfaces(vertices, simplices, path=None, name=’_surfaces’, alpha=1)[source]

Export data to a vtk file for posterior visualizations

Parameters:
  • geo_data (gempy.InputData) – All values of a DataManagement object
  • block (numpy.array) – 3D array containing the lithology block
  • path (str) – path to the location of the vtk
Returns:

None

planesCallback(obj, event)[source]

Function that rules what happens when we move a plane. At the moment we update the other 3 renderers and update the pandas data frame

render_model(size=(1920, 1080), fullscreen=False)[source]

Method to launch the window

Parameters:
  • size (tuple) – Resolution of the window
  • fullscreen (bool) – Launch window in full screen or not

Returns:

set_camera_backcolor(color=None)[source]

define background colors of the renderers

set_interfaces()[source]
Create all the interfaces points and set them to the corresponding renders for their posterior visualization
with render_model
Returns:None
set_orientations()[source]

Create all the orientations and set them to the corresponding renders for their posterior visualization with render_model :returns: None

set_surfaces(vertices, simplices, alpha=1)[source]

Create all the surfaces and set them to the corresponding renders for their posterior visualization with render_model

Parameters:
  • vertices (list) – list of 3D numpy arrays containing the points that form each plane
  • simplices (list) – list of 3D numpy arrays containing the verticies that form every triangle
  • formations (list) – ordered list of strings containing the name of the formations to represent
  • fns (list) – ordered list of formation numbers (int)
  • alpha – Opacity of the plane
Returns:

None

sphereCallback(obj, event)[source]

Function that rules what happens when we move a sphere. At the moment we update the other 3 renderers and update the pandas data frame