Variograms

The typical variogram modeling workflow is an iterative process:

  1. Examine the data to determine a reasonable set of directions and data spacing values for lag distances
  2. Calculate experimental variograms in the chosen directions (major, minor, vertical)
    • Refine these direction choices iteratively, considering omnidirectional and alternative variogram types (PR, Correlogram etc)
  3. Evaluate and interpret the calculated variograms
  4. Model the calculated variograms with a consistent 3D model
    • Tune the variogram model with the assistance of automatic fitting
  5. Use the variogram, and possible come back around to compare with simulation results

Introduction

Variogram modeling in pygeostat centers around the gs.Variogram. A brief overview on python class objects can be found in the the python introduction section.

As a quick summary, the gs.Variogram class contains:

  • Wrappers for variogram calculation, modeling, and plotting functions
  • Container for experimental variogram calculation parameters for both point and gridded data
  • Container for variogram modeling parameters
  • Container for experimental (point and gridded data) and modeled variogram data
  • Functions to aid in setting various variogram parameters

As a very quick example, lets set-up a variogram object:

>>> vario = gs.Variogram(datafl=datafl, var=var)

On the user-end, to call variables stored within the an instance of the variogram class, you would replace self with vario. For example, lets say the major directions azimuth, dip, and tilt were already saved in the object; to see what it is, you would call the variable containing that information by:

>>> vario.major
    [40, 0, 0]

There are various functions stored within the Variogram objects as well. For example, once the variogram parameters and directions are set, the experimental variogram calculation function can be used by calling:

>>> vario.varcalc()

Potentially considering many directions, and possibly different tolerances and types for a single direction makes establishing a generic class which does everything challenging. This class is setup for 2 common use cases:

  1. Defined major, minor and vertical directions - used for calculation and modeling

    • In a 3D omnidirectional case, directional data is stored in the major
  2. Generic set of specified parameters - flexible but requires more manual setup

There are two iPython Notebooks within the pygeostat distribution. Please review them for example workflows. They are located in:

>>> ./examples/variograms

Variogram Object Attributes

Variables and Data

All of the variables and data stored in the gs.Variogram class are as follows:

Variables Description
self.datafl The gs.DataFile class that contains the required data
self.ndim The number of dimensions of the system being investigated
self.ndir The number of principle directions being calculated
self.omnihz Indicator if a horizontal omni-directional variogram is being used
self.omni3d Indicator if a 3D omni-directional variogram is being used
self.calcpars Dictionary containing variogram calculation parameters
self.modelpars Dictionary containing variogram modeling parameters
self.simpars Dictionary containing gridded data variogram calculation parameters
self.directions List containing boolean indicators if the [Major, Minor, and Vertical] directions are being used
self.major List containing the azimuth, dip, and tilt of major direction
self.minor List containing the azimuth, dip, and tilt of minor direction
self.vert List containing the azimuth, dip, and tilt of vertical direction
self.majorexp Pandas dataframe containing experimental variogram information
self.minorexp Pandas dataframe containing experimental variogram information
self.vertexp Pandas dataframe containing experimental variogram information
self.model gs.VariogramModel class containing the variogram model
self.simdatfl The gs.DataFile class that contains the simulated data or a string pointing to the location of the GSLIB output simulation file (if it is too large to load to python)
self.majorsim The rotation matrix equivalent of the grid steps used by gs.varsim() to calculate the Major direction
self.minorsim The rotation matrix equivalent of the grid steps used by gs.varsim() to calculate the Minor direction
self.vertsim The rotation matrix equivalent of the grid steps used by gs.varsim() to calculate the Vertical direction
self.simexp Pandas dataframe containing the experimental variograms of all simulated realizations
self.simavgexp Pandas dataframe containing the average experimental variograms of all simulated realizations

Calcpars

As indicated in the above table, the the dictionaries self.calcpars, self.modelpars, and self.simpars contain all of the parameters used by gs.varcalc(), gs.varmodel(), and gs.varsim(). Please refer to the their documentation on specific details. The bellow tables provide information as to what each parameter is, how to update it, and if it is automatically populated, where it’s information is sourced. They are only automatically filled if they are not explicitly set by the user. All parameters mimic the GSLIB programs they are derived from.

varcalc Description
'variotypes' A single, permiciple variogram type enumerator. Default value of 1. Update using self.update_calcpars().
'variotails' A single string indicating which columns in self.datafl contains the tail variable being used. Update using self.update_calcpars().
'varioheads' A single string indicating which columns in self.datafl contains the head variable being used. Update using self.update_calcpars().
'azm' A list of as many azimuth values as their are directions being calculated. Update using either self.inferdirections, self.setdirections,or self.update_calcpars().
'dip' A list of as many dip values as their are directions being calculated. Update using either self.inferdirections, self.setdirections,or self.update_calcpars().
'tilt' A list of as many tilt values as their are directions being calculated. Update using either self.inferdirections, self.setdirections,or self.update_calcpars().
'azmtol' A list of as many azimuth tolerance values as their are directions being calculated. Update using either self.settols or self.update_calcpars().
'diptol' A list of as many dip tolerance values as their are directions being calculated. Update using either self.settols or self.update_calcpars().
'nlags' A list of as many number of lags to calculate as there are directions being calculated. Update using either self.settols or self.update_calcpars().
'lagdist' A list of as many lag distances to use as there are directions being calculated. Update using either self.settols or self.update_calcpars().
'lagtol' A list of as many lag distance tolerance values as their are directions being calculated. Update using either self.settols or self.update_calcpars().
'bandhorz' A list of as many horizontal bandwidth values as there are directions being calculated. Deafult value of 1e+21. Update using either self.settols or self.update_calcpars().
'bandvert' A list of as many vertical bandwidth values as there are directions being calculated. Deafult value of 1e+21. Update using either self.settols or self.update_calcpars().
'tmin' A single minimum trimming limit. Default value of -999.0. Update using self.update_calcpars().
'tmax' A single maximum trimming limit. Default value of 1.0e21. Update using self.update_calcpars().
'variostd' A single boolean value indicating if the variograms should be standardized. Default value of False. Update using self.update_calcpars().
'variosills' A single sill value to standardize to if 'variostd' is True. Default value of None indicates that varcalc calculate an appropriate sill. Update using self.update_calcpars().
'strict' Logical indicating whether to run some checks for common errors. Setting to true will slightly slow execution but could catch common user errors. Default value of True. Update using self.update_calcpars().
'variocutcat' Cutoff/categories for 'variotypes' 9-12 otherwise these values will be ignored. Default value of None will be transformed to a list of zeroes. Update using self.update_calcpars().

Modelpars

varmodel Description
'vargstr' List indicating an acceptable range that varmodel can use for the sill. Default value is [0.01, 1.0e21] if self.calcpars['variostd'] is False, otherwise it defaults to [1.0, 1.0]. Update using self.update_modelpars().
'nst' List indicating an acceptable number of variogram structure(s) that varmodel can use. Default value is 1. Update using self.update_modelpars().
'c0' List indicating an acceptable range that varmodel can use for the nugget effect. Default value is [1.0e-6, sill]. Update using self.update_modelpars().
'it' List indicating an acceptable range of variogram structure types that varmodel can use. Default value is [1, 3]. Update using self.update_modelpars().
'cc' List indicating an acceptable range of values that varmodel can use for the covariance contribution of each structure. Default value is [0.1, sill]. Update using self.update_modelpars().
'ang1' List indicating an acceptable range of values that varmodel can use for the azimuth rotation angle. Default value is sourced from self.calcpars['azm']. Update using self.update_modelpars().
'ang2' List indicating an acceptable range of values that varmodel can use for the dip rotation angle. Default value is sourced from self.calcpars['dip']. Update using self.update_modelpars().
'ang3' List indicating an acceptable range of values that varmodel can use for the tilt rotation angle. Default value is sourced from self.calcpars['tilt']. Update using self.update_modelpars().
'ahmax' List indicating an acceptable range of values that varmodel can use for the major direction range. Default range is [1, 1.0e21]. Update using self.update_modelpars().
'ahmin' List indicating an acceptable range of values that varmodel can use for the minor direction range. Default range is [1, 1.0e21]. Update using self.update_modelpars().
'avert' List indicating an acceptable range of values that varmodel can use for the vertical direction range. Default range is [1, 1.0e21]. Update using self.update_modelpars().

Simpars

varsim Description
'lithfl' File with lithology information. Update using self.update_simpars().
'lithcol' Column within 'lithfl' that the lithology information is located. Update using self.update_simpars().
'lithcode' Code indicating if lithology information is used (0=not used). Update using self.update_simpars().
'datafl' A gs.DataFile class, a valid path to a data file containing simulated realizations values, or a list of file locations each with only 1 realization. Update using self.addsimdatafl(). If a gs.DataFile is passed, a gs.GridDef must be attached to it.
'nvar' Number of variables used to calculate the variogram values. Taken from self.calcpars unless explicitly set using self.update_simpars().
'varcols' Column(s) within 'datafl' that the variable(s) being used to calculate the variogram are located in. Only needed if 'datafl' is not a gs.DataFile, otherwise it is automatically set to '1 1' or '1 2' depending on what 'nvar' is. Update using self.addsimdatafl().
'tmin' Lower trimming limit. Default value of -999.0. Update using self.update_simpars().
'tmax' Upper trimming limit. Default value of 1.0e21. Update using self.update_simpars().
'outvargs' Output file for variograms of realizations. Defaults to a temporary file that is read into python then deleted. Update using self.update_simpars().
'outavg' Output file for average variogram. Defaults to a temporary file that is read into python then deleted. Update using self.update_simpars().
'griddef' gs.GridDef class. No default value. Will grab the gs.GridDef attached to a gs.DataFile passed into self.addsimdatafl(). Otherwise, use self.addgriddef().
'nreal' Number of realizations within the simulated realization file. Default value is calculated if 'datafl' is a gs.DataFile class and 'griddef' has been set. Update using self.update_simpars().
'ndir' The number of principle directions being calculated. Default value is set to the value of self.ndir. Update using self.update_simpars().
'nlags' Number of lags to calculate. Default value is calculated from two values. The first being the maximum lag distance used for experimental variogram calculation, as determined from values in self.calcpars. The second values used is based on the maximum lag distance calculated from the grid steps for each direction being calculated. Update using self.update_simpars().
'dirs' A string containing lines of grid cell shift values for each direction. The default values are calculated using gs.get_varsimpars() which use what directional data exists within self.major, self.minor, and self.vert. Update using self.update_simpars().
'stdsill' Boolean indicator if the variograms should be standardized or not. Default value is False. Update using self.update_simpars().
'nvargs' Number of required variograms. The gs.Variogram class is only set up to handle a value of 1, which is what the default value is. Update using self.update_simpars().
'vargpars' The used variable enumerators and variogram type. The default value uses 'variotypes' from self.calcpars['variotypes']. In the case where 'nvar' is 1, the value is set to '1 1 variotypes'. In the case where 'nvar' is 2, the value is to '1 2 variotypes'. Update using self.update_simpars().
'tailcol' If 'datafl' is a gs.DataFile class, 'tailcol' needs to be the column name which contains the tail variable data. If 'datafl' is a file location, 'tailcol' needs to be a column number, referencing the tail variable data. Update using self.addsimdatafl().
'headcol' If 'datafl' is a gs.DataFile class, 'headcol' needs to be the column name which contains the head variable data. If 'datafl' is a file location, 'headcol' needs to be a column number, referencing the head variable data. Update using self.addsimdatafl().

Experimental Variogram Types

Value Experimental variogram type (from VarCalc)
1 traditional semivariogram
2 traditional cross semivariogram
3 covariance (-3 calculates variance (provided sill) -covariance)
4 correlogram (-4 calculates 1-correlation)
5 general relative semivariogram
6 pairwise relative semivariogram
7 semivariogram of logarithms
8 semimadogram
9 indicator semivariogram - continuous - requires a cutoff
10 indicator semivariogram - categorical - requires a category
11 traditional cross indicator semivariogram - continuous
12 traditional cross indicator semivariogram - categorical

Model Variogram Types

Value Model variogram type (from VarModel) Equation
1 Spherical \(\gamma(h) = c \left[ 1.5\frac{h}{a} - 0.5\left(\frac{h}{a}\right)^3 \right]\)
2 Exponential \(\gamma(h) = c \left[ 1 - \exp \frac{-3h}{a} \right]\)
3 Gaussian \(\gamma(h) = c \left[ 1 - \exp \left(\frac{-3h}{a}\right)^2 \right]\)
4 Power Model \(\gamma(h) = c h^\omega\)
5 Hole Effect \(\gamma(h) = c \left[ 1 - \cos\left( \frac{h}{a}\pi \right) \right]\)

Which of the Major, Minor, and Vertical Direction Variables are Used?

Using the arguments passed during initialization, the required principle directions are determined. The following parameters are what are used: self.ndir, self.ndim, and self.omnihz.

Using those three parameters, the variable self.directions is set up as follows:

if self.ndir == 3:
    self.directions = [True, True, True]
elif self.ndir == 2 and self.omnihz and self.ndim == 3:
    self.directions = [True, False, True]
elif self.ndir == 2 and self.ndim == 2:
    self.directions = [True, True, False]
elif self.ndir == 1 and self.omnihz and self.ndim == 2:
    self.directions = [True, False, False]
elif self.ndir == 1:
    self.directions = [True, False, False]

The outcome self.directions also indicates which of the variables self.major, self.minor, and self.vert will contain information and in turn… self.majorexp, self.majorsim, ect.

As an example of how the code would work above, in the case where a horizontal omni-directional variogram and a vertical variogram is being calculated, the directions self.major and self.vert are being used while self.minor is left to a value of None.

Variogram Class

class pygeostat.variograms.variogram.Variogram(datafl=None, var=None, ndim=3, ndir=3, omnihz=False, omni3d=False, variostd=False, catpars=None, calcpars=None, modelpars=None, simpars=None, model=None, mute=False)

Variogram class that enables a simplified variogram work flow due to the storage of parameters within the class and it’s integration of various variogram utilities.

Parameters:
  • datafl (gs.DataFile) – A gs.DataFile class that must contain the appropriate coordinate information (i.e., x, y, and z attributes).
  • var (list) – 1 or 2 column name(s) indicating which columns in datafl contain the variables being used. The tail variable MUST be listed first then the head variable if the variogram class is for a cross-variogram
  • ndim (int) – The number of dimensions of the system being investigated. If your working in a 3D system but are only using two principle directions, you’ll still need to put 3. The variogram class will mishandle variograms otherwise
  • ndir (int) – The number of principle directions being calculated
  • omnihz (bool) – Indicate if a horizontal omni-directional variogram is being used
  • omni3d (bool) – Indicate if a 3D omni-directional variogram is being used
  • variostd (bool) – Indicate if the variograms should be standardized
  • catpars (tuple) – Indicate if the current variogram should consider a category column and a category defined in that column, e.g. catpars = ("categories", 0) would indicate to consider the "categories" column, and category 0 in that column.
  • calcpars (dict) – Optional dictionary of experimental variogram calculation parameters Does not have to include all parameters, will update the default values
  • modelpars (dict) – Optional dictionary of variogram modeling parameters. Does not have to include all parameters, will update the default values
  • simpars (dict) – Optional dictionary of experimental realization variogram calculation parameters. Does not have to include all parameters, will update the default values
  • model (gs.VariogramModel) – Optional gs.VariogramModel class of a already defined variogram model
  • mute (bool) – Indicate if printouts from the variogram class functions should be muted. Only recommended if looping variogram classes. Printouts were deliberatively placed so the user can verify automation at all steps of variography.

Example

Setup the datafile containing the point data:

>>> datafl = gs.DataFile(flname='../data/3DDecor.dat', readfl=True, dftype='point',
...                      fltype='gslib')

Generate the variogram class from the imported data

>>> vario = gs.Variogram(datafl, 'Var1', variostd=True)

Code author: Warren E. Black, Jared Deutsch, and Tyler Acorn - 2015-11-09

Add a gs.GridDef

Variogram.addgriddef(griddef)

In the case where a simulated data file or files are loaded into the Variogram class using self.addsimdatfl() is a string indicating the file location or a list of file locations, you will need to add a pointer to a gs.GridDef class manually. This function will store a pointer in the self.simpars['griddef'] parameter.

Parameters:griddef (GridDef) – A properly loaded gs.GridDef class

Code author: Warren E. Black - 2016-04-07

Add Simulated Data File

Variogram.addsimdatfl(simdatfl, varcols=None)

Add and validate the gridded data file and variable column information and save the information to the appropriate self.simpars parameters. simdatfl is used to update self.simpars['datafl'] and varcols is used to update self.simpars['varcols'].

The function accepts three types of inputs passed to the simdatfl argument:

  • A gs.DataFile class
  • A string containing the location of a single GSLIB formated file. All realizations for each variable needed must be in a single column. In the case where two variables are required for cross-variograms, the file must contain a column for each simulated variable
  • A string containing a wildcard ‘*’, where it is assumed that that nreal data files are named through replacing the wild card with i = 1,…,nreal integer values
  • A list containing numerous GSLIB formated files. Each file contains the required variables in separate columns, but must only contain a single realization.
  • A python generator object that yields a 1-D numpy array. Note: python generators cannot be pickled.
Parameters:simdatfl – Input simulation data in a form specified above.
Keyword Arguments:
 varcols (list) – The column(s) which contain the variable(s) being studied. Required when simdatfl is not a python generator

Example

Add the simulated file to the variogram class object ‘vario’ created previously:

>>> vario.addsimdatfl(simdatfl='../data/3DDecor_PVar1_sim.out', varcols=1)

All realizations were writen to a single HDF5 file. In this example, all the realizations are in the root directory of the HDF5 file and each dataset contains the string “Real”. A pygeostat H5Store class is created and the iteritems method is utilized:

>>> # Create the H5Store object
>>> h5dat = gs.H5Store('./sgsim.h5')
>>> simdat = h5dat.iteritems(wildcard="Real")
>>>
>>> # Pass the generator to addsimdatfl
>>> variogram.addsimdatfl(simdatfl=simdat)

In this example, multivariate simulation was completed where each realization was written to its own file. A python generator is created and passed:

>>> # Create the generator
>>> def iter_simreals():
>>>     for ireal in range(100):
>>>         simfl = './real_%s.h5' % ireal
>>>         data = gs.H5Store(simfl)
>>>         yeild data['NS_AU']
>>>
>>> # Pass the generator to addsimdatfl
>>> variogram.addsimdatfl(iter_simreals())

Code author: Warren E. Black - 2015-11-09

Fit Variogram Model

Variogram.fitmodel(expdata=None, nst=None, sill=None, maxiter=250000, rseed=132421, minpairs=10, npairswt=True, invdistwt=True, fixhmaxvertanis=False, hmaxvertanis=1.0, fixhminhmaxanis=False, hminhmaxanis=1.0)

Fit a variogram model to the experimental variogram data using the fortran modules

Parameters:
  • expdata (list?) – Variogram data to be fit
  • nst (int) – Number of nested structures to fit
  • sill (float) – Variogram sill
  • maxiter (int) – Max autofit iterations
  • rseed (int) – Random number seed
  • minpairs (int) – Mininmum number of pairs to fit a point to
  • npairswt (bool) – Weight exp points by the number of pairs found for each point
  • invdistwt (bool) – Weight exp points by IDW to origin
  • fixhmaxvertanis (bool) – Fix the max/vert anisotropy
  • hmaxvertanis (float) – Fixed max/vert anisotropy ratio (this may be vert/max.. not sure!)
  • fixhminhmaxanis (bool) – Fix the max/min anisotropy
  • hminhmaxanis (float) – Fixed max/min anisotropy ratio (this may be min/max.. not sure!)

Code author: Jared Deutsch - 2015-11-09

Interpret VarModel Parameters

Variogram.fitparsfromdict(modelpars=None, sill=None, nst=None)

Interpret a dictionary of model parameters for fitting parameters

Parameters:
  • modelpars (dict) – Dictionary of model parameters to be updated?
  • sill (float) – Sill value to update
  • nst (int) – Number of nested structures
Returns:

nst, c0, it, cc, ang1, ang2, ang3, ahmax, ahmin, avert, optsill

Return type:

All parameters parsed from the self.modelpars dictionary

Infer Directions

Variogram.inferdirections(azm=None, dip=None, tilt=None)

Sets the variogram calculation ‘azm’, ‘dip’, and ‘tilt’ parameters for all directions simultaneously. Meant for use when a single set of final azimuth, dip, and tilt parameters can be used to extrapolate all of the direction(s) direction parameters. Only single values can be passed for each parameter.

Parameters:
  • azm (float) – The azimuth of the major direction
  • dip (float) – The dip of the major direction
  • tilt (float) – The tilt of the minor or vertical

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

A simple call with it’s output:

>>> vario.inferdirections(azm=40, dip=0, tilt=0)
The major direction has been set to a [azm, dip, tilt] value of: [40.0, -90.0, 0]
The minor direction has been set to a [azm, dip, tilt] value of: [130.0, -90.0, 10.0]
The vert direction has been set to a [azm, dip, tilt] value of: [40.0, 0.0, 10.0]

Code author: Warren E. Black - 2015-11-09

Initialize from File

Variogram.init_from_file(file, major=None, minor=None, vert=None)

Initialize the experimental data from an externally generated VarCalc formatted datafile. The experimental data is added to the major, minor or vertical variogram depending on what variogram numbers are passed for each.

Parameters:
  • file (str) – the path to the experimental output
  • major (int) – the varg number in the file for major data
  • minor (int) – the varg number in the file for minor data
  • vert (int) – the varg number in the file for vert data

Code author: Ryan Martin - 2017-02-24

Legacy Output

Variogram.legacyout(outfl, vari=None, append=False)

Wrapper for gs.legacyout()

Plotting

Variogram.plot(experimental=True, model=False, sim=False, sill=None, legend=None, titles=True, suptitle=None, trimylim=True, pltstyle=None, cust_style=None, outfl=None, xlim=None, ylim=None, colors=None, printmodel=False, figsize=None, separate_plts=True, axes=None, varplt_kws=None, **kwargs)

Quick plotting function to view either experimental variograms, modeled variograms, or both.Several style related kwargs in varplt, including unit, grid, axis_xy are not provided as kwargs here, but may be accessed indirectly through setting their associated gsParams: gsParams['plotting.grid'], gsParams['plotting.axis_xy'], gsParams['plotting.unit']

Note

This function returns the figure as a mpl.fig object. If you are using iPython and plotting figures inline, two figures may appear. You can display only one image at a time be doing the following:

>>> fig = Variogram.plot()
Keyword Arguments:
 
  • experimental (bool) – Used to specify experimental style of the plot
  • model (bool) – Used to specify model style of the plot
  • sim (bool) – Used to specify sim style of the plot
  • sill (float) – value for plotting the sill. Set to False if you don’t want to plot the sill. Default value pulled from self.calcpars['variosills'] if it is not also None.
  • legend (int) – A permissible matplotlib legend() location code. Set to False if no legend is wanted. Default value is 4 (lower right) and 0 (best) if the sill is negative.
  • titles (bool, str, or list) – If set to True, titles will be automatically generated, for each direction and variogram type. A single string value can be passed if a single plot is going to be generated. Or, a list can be passed containing unique titles for each of the axis plots; values of the list can be set to a value of None if a title isn’t desired for a particular axis.
  • suptitle (str) – Add a centered title to the figure. Will only work if seperate_plts is set to True, otherwise
  • trimylim (bool) – Indicate if realization plots should cease plotting once they pass the limits specified by ylim.
  • pltstyle (str) – Name of class to use. default 'ccgpaper'
  • cust_style (str OR dict) – set the plot style. 'ccgpaper' or 'thesis' or custom dictionary are your options. See gs.set_style() for more information.
  • outfl (str) – Name of file to save figure too.
  • xlim (tuple) – default is None. Can pass a list of tuples [(min1, max1), (min2, max2)]
  • ylim (tuple) – default is None. Can only pass a single tuple that is applied to all plots
  • colors (str or list) – Can take either a string (name of color map to draw from or matplotlib color) or a list of colors to use in this order [major, minor, vert]. if a list of length 1 is passed it will copy that color to all other directions.
  • printmodel (bool) – print the variogram model on one of the plotting axes, usually the last one which extends the figure size
  • figsize (tuple) – Tuple of float values specifiying the figure size in inches
  • separate_plts (bool) – Boolean to control separate plots for each variogram
  • axes (mpl.axis or 2-axes) – Pass a single axis and a 2-long axis tuple or list. The single axis is only be used if separate_plts is False and the set of 2 axes is used for a case where ndir == 2,
  • varplt_kws

    Nested dictionary of permissible keyword arguments to pass to gs.varplt() or gs.varpltsim(). As there are 3 types of variograms plotting 'experimental', 'model', 'sim', this dictionary is nested. On the root level of the dictionary only these three strings can exist. They then point to the required keyword arguments to then pass on to the respective gs.varplt() call for 'expierimental' and 'model' or gs.varpltsim() for 'sim'. For example:

    >>> varplt_kws = {'experimental': {'ms': 5},
    ...               'model': {'lw': 2}},
    ...               'sim': {'alpha': 0.25}}
    
  • **kwargs – Optional permissible keyword arguments to pass to. gs.exportimg()
Returns:

Figure object containing a single axis or multiple axes. Use fig.get_axes() if you’d like to modify one or more of the axes

Return type:

fig (mpl.fig)

Example

Use the previously created variogram object and its data to generate a set of plots for the major, minor and vertical directions

>>> vario.plot(experimental=False, model=True, sim=True, xlim=(0,10))

Code author: Warren E. Black and Tyler Acorn - 2015-11-09

Set Directions

Variogram.setdirections(major=None, minor=None, vert=None)

Sets (or updates) the specified principal directions (i.e., Major, Minor, and/or Vert) stored within the variogram class. An alternative to inferdirections() if it does not allow enough flexibility or if directions want to be set explicitly.

A valid direction is a list with up to 3 components (i.e., [ Azimuth, Dip, Tilt ]) or a dictionary with up to 3 directions (i.e., { ‘azm’ ‘dip’ ‘tilt’ })

For the 3D definition, see: http://geostatisticslessons.com/lessons/anglespecification

Parameters:
  • major (list or dict) – The components of the major principal direction
  • minor (list or dict) – The components of the minor principal direction
  • vert (list or dict) – The components of the vert principal direction

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

Simple usage:

>>> vario.setdirections(major=[300.0, 90.0, 0], minor=[30.0, 90.0, 0.0],
                        vert=[300.0, 0.0, 0.0])

Check to make sure the directions were saved properly to the vario gs.Variogram class:

>>> vario.major
[300.0, 90.0, 0]
>>> vario.minor
[30.0, 90.0, 0.0]
>>> vario.vert
[300.0, 0.0, 0.0]

Code author: Jared Deutsch - 2015-11-09

Set Tolerances

Variogram.settols(nlags=None, lagdist=None, lagtol=None, azmtol=None, diptol=None, bandhorz=None, bandvert=None)

Sets (or updates) the variogram calculation tolerance parameters. Not all parameters have to be specified at the same time. If one is already set, it is not replaced unless its corresponding keyword argument is passed.

Parameters passed should have either one value that can be used for all directions, or a list of parameters for each direction being calculated.

This function is an alternative to gs.update_calcpars() if it is too unstructured, and allows users to help break up the parameter set-up workflow.

Parameters:
  • nlags (list) – Number of lags to calculate
  • lagdist (list) – Lag distance
  • lagtol (list) – Lag distance tolerance
  • azmtol (list) – Azimuth tolerance
  • diptol (list) – Dip tolerance
  • bandhorz (list) – Horizontal bandwidth
  • bandvert (list) – Vertical bandwidth

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

A simple call:

>>> vario.settols(nlags=[30,30,10], lagdist=[75,75,3], lagtol=[37.5,37.5,1.5],
...               azmtol=[22.5,22.5,5.0], diptol=5.0, bandvert=5.0)

Update the lag distance, lag tolerance, and azimuth tolerance. Only the specified parameters are changed, the previously set parameters remain unchanged:

>>> vario.settols(lagdist=[50,50,3], lagtol=[25,25,1.5], azmtol=[15.5,15.5,5.0])

Code author: Jared Deutsch - 2015-11-09

Update Calculation Parameters

Variogram.update_calcpars(variotypes=None, azm=None, dip=None, tilt=None, azmtol=None, diptol=None, nlags=None, lagdist=None, lagtol=None, bandhorz=None, bandvert=None, variosills=None, strict=True, variocutcat=None, variostd=None)

Allow an unstructured amount of keywords arguments to be used to set (or update) the self.calcpars parameter dictionary. Please refer to the list of variogram calculation parameters to see which parameters can be passed.

As many or few parameters can be set at once.

Parameters:**kwargs (keyword arguments) – Any permissible variogram calculation parameters passed as keyword arguments

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

A simple call updating only the azimuth, dip, lag distance, and tilit parameters:

>>> vario.update_calcpars(azm=[40.0, 130.0, 40.0], dip=[0.0, 0.0, -90.0],
...                      lagdist=[75, 75, 3], tilt=0)

Permissible Calcpars

variotypes=1, variotails=None, varioheads=None, azm=None, dip=None, tilt=None,
azmtol=None, diptol=None, nlags=None, lagdist=None, lagtol=None, bandhorz=1e+21,
bandvert=1e+21, variostd=variostd, variosills=None, strict=True,
variocutcat=None, variostd=True

Code author: Warren E. Black - 2015-11-09

Update Model Parameters

Variogram.update_modelpars(vargstr=None, nst=None, c0=None, it=None, cc=None, ang1=None, ang2=None, ang3=None, ahmax=None, ahmin=None, avert=None, **kwargs)

Allow an unstructured amount of keywords arguments to be used to set (or update) the self.modelpars parameter dictionary. Please refer to the list of variogram modeling parameters to see which parameters can be passed.

As many or few parameters can be set at once.

Parameters:**kwargs (keyword arguments) – Any permissible variogram modeling parameters passed as keyword arguments

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

A simple call:

>>> vario.update_modelpars()

Permissible ModelPars

vargstr=None, nst=None, c0=None, it=None, cc=None,
ang1=None, ang2=None, ang3=None, ahmax=None, ahmin=None,
avert=None

Code author: Warren E. Black - 2015-11-09

Update VarSim Parameters

Variogram.update_simpars(lithfl=None, lithcol=None, lithcode=None, datafl=None, nvar=None, varcols=None, tmin=None, tmax=1e+21, outvargs=None, outavg=None, griddef=None, nreal=None, ndir=None, nlags=None, dirs=None, stdsill=None, nvargs=None, vargpars=None, tailcol=None, headcol=None)

Allow an unstructured amount of keywords arguments to be used to set (or update) the simpars parameter dictionary. Please refer to the list of simulation variogram parameters to see which parameters can be passed.

As many or few parameters can be set at once.

Parameters:**kwargs (keyword arguments) – Any permissible variogram modeling parameters passed as keyword arguments

Example

Set-up a gs.Variogram class:

>>> vario = gs.Variogram(datafl=datafl, var=var)

A simple call:

>>> vario.update_simpars()

Permissible SimPars

lithfl=None, lithcol=None, lithcode=None, datafl=None,
nvar=None, varcols=None, tmin=None, tmax=1e+21,
outvargs=None, outavg=None, griddef=None, nreal=None,
ndir=None, nlags=None, dirs=None, stdsill=None, nvargs=None,
vargpars=None, tailcol=None, headcol=None

Code author: Warren E. Black - 2015-11-09

Validate VarCalc Parameters

Variogram.validcalcpars()

Prior to running varcalc, a series of checks are completed to ensure the variogram calculation parameters are valid. These checks are as follows:

  1. Reject parameter values that are a dictionary and convert tuples to lists
  2. The parameters listed below do not have a value of None
>>> 'azm', 'dip', 'tilt', 'nlags', 'lagdist', 'lagtol', 'azmtol', 'diptol', 'bandhorz',
... 'bandvert'
  1. If one of the parameters listed below only has one value set but there are more than one direction being calculated, repeat that value for every direction.
>>> 'nlags', 'lagdist', 'lagtol', 'azmtol', 'diptol', 'bandhorz', 'bandvert'
  1. Make sure the directional parameters (i.e., 'azm', 'dip', and 'tilt') have as many values as there are directions
  2. Make sure there are only as many values for each parameter as there are directions
  3. Make sure the directions within the class are set if there is directional data in the calculation parameters dictionary
  4. Make sure the the directional information stored in the class (i.e., self.major ect.) matches what is within self.calcpars. If not, update them with what is found in self.calcpars
  5. Ensure a variogram type has been set
  6. If the head and tail variables are different, set 'variotype' to 2 (i.e., traditional cross semivariogram)
  7. If calculating a transitional cross semivariogram, set the sill to the covariance between the two variables.

Code author: Warren E. Black - 2015-11-09

Validate Directions

Variogram.validdirection(direction)

Validate a direction - return the corrected values.

Parameters:direction (dict or list) – 3 directions poteitally ..
Returns:input directions converted to float
Return type:validdir (list)

Code author: Jared Deutsch - 2015-11-09

Variogram Calculation Wrapper

Variogram.varcalc()

Calculates an experimental variogram. No parameters are passed to this function.

See also

The function self.validcalcpars() is called which completes a variety of checks. Please refer to it’s documentation for information on what is done.

Example

Calc the variograms for the created ‘vario’ class

>>> vario.varcalc()

Code author: Jared Deutsch - 2015-11-09

Variogram Grid Calculation Wrapper

Variogram.varsim(minstep=None, maxpct=None, maxstep=None, maxlags=500, minlags=40, degerr=0.5, vargrng=None, optlags=None, rngmult=1.25, dbg_fl=False, tmin=None, tmax=1e+21, nprocess=None)

Wrapper for varsim that extrapolates nearly all of the needed parameters. This is done by using the parameters previously set for experimental variogram calculation (i.e., the parameters set in self.calcpars). The only parameters that cannot be drawn automatically or do not have default values are:

>>> datafl, tailcol, headcol

If any of the above parameters need to be set, please use the function self.update_simpars as this function calls upon the self.simpars dictionary for its parameters. The function self.addsimdatfl() updates self.simpars['datafl'], self.simpars['tailcol'], and self.simpars['headcol'] when it is called.

Warning

The number of grid nodes that must be shifted in each direction, if not already set in the parameter self.simpars['dirs'] is calculated using the function gs.get_varsimpars() or gs.get_varsimpars_2d(). If the system is 3-D, make sure you verify the idx, idy, idx parameters found as it is an experimental function.

Note

In the case where the simulated realizations are too large to load into python, a string containing a single file location can be saved to self.simpars['datafl'] which will be used to draw the realization data from. If this method is used the following parameters (in addition to the previously stated parameters) will need to be set manually if required:

>>>  ``nreal``, ``lithfl``, ``lithcol``, ``lithcode``

These parameters can be set using self.update_simpars.

Parameters:
  • minstep (float) – Used if self.simpars['dirs'] is not set. Minimum steps allowed in each direction. This parameter is used by gs.get_varsimpars().
  • maxpct (float) – Used if self.simpars['dirs'] is not set. Maximum steps allowed in each direction as a percentage of the number of cells in the corresponding direction. This parameter is used by gs.get_varsimpars().
  • maxstep (float) – Used if self.simpars['dirs'] is not set. Maximum steps allowed in each direction. This parameter is used by gs.get_varsimpars().
  • maxlags (int) – Maximum number of lags to calculate before throwing an error. Helps when trying to tune the above three parameters. This parameter is used by gs.get_varsimpars() or gs.get_varsimpars_2d()
  • minlags (int) – Minimum number of lags to calcualte before throwing an error. Helps when trying to tune the top three parameters. This parameter is used by gs.get_varsimpars() or gs.get_varsimpars_2d()
  • degerr (float) – Allowable absolute deviation from the major direction. This parameter is used by gs.get_varsimpars_2d()
  • vargrng (float) – The desired final calculated range. Will not be multiples by the parameter rngmult and replaces the ahmax value pulled from the variogram model. Intended use is for when a variogram is going to be displayed on a plot that extends significantly past its maximum range. This parameter is used by gs.get_varsimpars_2d()
  • optlags (int) – Optimal number of lags that will cover the range of the major variogram range multiplied by the parameter rngmult. Only considered if there are values that fall within the range defined by degerr. This parameter is used by gs.get_varsimpars_2d()
  • rngmult (float) – Major direction range multiplier which is used to scale the major variogram range, which is the approximate maximum lag calculated. This parameter is used by gs.get_varsimpars_2d()
  • dbg_fl (bool) – If a string is passed then the varsim parameter file will be written out
  • tmin (float) – Trimming limits, which are required presently because input data is a GSLIB file for this wrapper
  • tmax (float) – Trimming limits, which are required presently because input data is a GSLIB file for this wrapper
  • verbose (float) – Report progress through the realizations if True
  • nprocess (int) – number of processes to use for parallelizing the variogram calculations across realizations. Drawn from gsParams[‘config.nprocess’].
  • None. (if) –
Returns:

  • simexp (pd.DataFrame) – Dataframe containing output realization data saved to self.simexp
  • simavgexp (pd.DataFrame) – Dataframe containing output averaged realization data saved to self.simavgexp

Code author: Warren E. Black - 2015-11-09 and Ryan M. Barnett - 2018-04-06

VariogramModel Class

class pygeostat.variograms.variogrammodel.VariogramModel(vargstr=None, nst=None, c0=None, it=None, cc=None, ang1=None, ang2=None, ang3=None, ahmax=None, ahmin=None, avert=None, parsestring=True)

Class containing variogram model and associated routines. This class is generated automatically by gs.Variogram class when calling varg.fitmodel(). This class is also extremely useful in scripting workflows since the string representation of the object is the model variogram, e.g. str(vargmodel) prints the GSLIB style variogram model that can be written to a parameter file.

Parameters:
  • vargstr (str) – a string containing the standard gslib-style variogram string
  • nst (int) – number of nested structures
  • c0 (float) – the nugget effect
  • it (int[nst]) – nst-long array of int providing the type variogram for each structure
  • cc (float[nst]) – nst-long array of float providing the variance contribution of each structure
  • ang1 (float[nst]) – nst-long array of float providing the ang1 rotation angles for each structure
  • ang2 (float[nst]) – nst-long array of float providing the ang2 rotation angles for each structure
  • ang3 (float[nst]) – nst-long array of float providing the ang3 rotation angles for each structure
  • ahmax (float[nst]) – nst-long array of float providing the range of the variogram in the horizontal major direction
  • ahmin (float[nst]) – nst-long array of float providing the range of the variogram in the horizontal minor direction
  • avert (float[nst]) – nst-long array of float providing the range of the variogram in the vertical direction
  • parsestring (bool) – default True to parse the vargstr

Examples

Generate a simple variogram model:

>>> variostring = """1 0.1
... 1 0.99 0 0 0
...      100 100 100"""
>>> vargmodel = gs.VariogramModel(variostring)
>>> print(vargmodel)
... 1    0.010000000              -  nst, nugget effect
... 1    0.990000000   0.0   0.0   0.0   - it,cc,ang1,ang2,ang3
...     100.0    100.0    100.0       - a_hmax, a_hmin, a_vert

To generate the set of points Distance and Variogram for plotting:

>>> vargmodel["major"]
...       Distance  Variogram
... 0      0.00000   0.000000
... 1      0.00001   0.000000
... 2      1.00000   0.024850
... 3      2.00000   0.039696
... 4      3.00000   0.054537

These points can also be generated by considering either calcpoints() or calcpoints_principle() functions:

>>> vargmodel.calcpoints(0, 0)
...       Distance  Variogram
... 0      0.00000   0.000000
... 1      0.00001   0.000000
... 2      1.00000   0.024850
... 3      2.00000   0.039696
... 4      3.00000   0.054537
... ....
>>> vargmodel.calcpoints_principle("major")
...       Distance  Variogram
... 0      0.00000   0.000000
... 1      0.00001   0.000000
... 2      1.00000   0.024850
... 3      2.00000   0.039696
... 4      3.00000   0.054537
... ....

Code author: Jared Deutsch - 2015-02-04

Calculate cmax

VariogramModel.cmax()

Returns cmax = c0 + sum(cc)

Code author: Jared Deutsch - 2015-02-04

Calculate Variogram Model Points

VariogramModel.calcpoints(azm, dip, lags=None, lagdist=None, nlags=None, returnstr=False, epslon=1e-05)

Calculates variogram model points along the vector specified by the azimuth (azm) and dip provided at either points given by either:

  1. lags - list of lags. ie: [0.0, 1.0e-5, 10.0, 20.0] note that 1.0e-5 is for the nugget effect - an arbitrarily small distance
  2. or lagdist + nlags - lag distance and number of lags to calc.
Parameters:
  • azm (float) – the azimuth to calculate points along
  • dip (float) – the dip to calculate points along
Keyword Arguments:
 
  • lags (list) – a list of lag distances (as above)
  • lagdist (float) – the lag distance
  • nlags (int) – the number of lags to calculate
  • returnstr (bool) – write the output points to a formatted string
  • epslon (float) – a small value
Returns:

  • vmodelstr (pd.DataFrame) – with ‘Distance’ and ‘Variogram’ if returnstr==False
  • vmodelstr (str) – in vmodel.var output format if returnstr==True

Note

this function returns the covariance along the vector which may not coincide with the principle directions implied by the GSLIB angles: ang1, ang2, ang3

Code author: Jared Deutsch - 2015-02-04

Calculate Variogram Model Points along Principle Directions

VariogramModel.calcpoints_principle(which, lags=None, lagdist=None, nlags=None, returnstr=False, epslon=1e-05)

Calculates variogram model points along the principle directions implied by the modeled parameters in the variogram object.

  1. lags - list of lags. ie: [0.0, 1.0e-5, 10.0, 20.0] note that 1.0e-5 is for the nugget effect - an arbitrarily small distance
  2. or lagdist + nlags - lag distance and number of lags to calc.
Keyword Arguments:
 
  • lags (list) – a list of lag distances (as above)
  • lagdist (float) – the lag distance
  • nlags (int) – the number of lags to calculate
  • returnstr (bool) – write the output points to a formatted string
  • epslon (float) – a small value
Returns:

  • vmodelstr (pd.DataFrame) – with ‘Distance’ and ‘Variogram’ if returnstr==False
  • vmodelstr (str) – in vmodel.var output format if returnstr==True

Note

this function returns the variogram values along the principle directions implied by the GSLIB directions: ang1, ang2, ang3

Code author: Jared Deutsch - 2015-02-04 - Modified by Ryan Martin - 2016-07-29

Copy

VariogramModel.copy()

Return a copy of this object

Get Covariance Between two Points

VariogramModel.getcova(x1, y1, z1, x2, y2, z2, variogram=False)

Uses GSLIB cova3 routine to calculate covariance using the variogram model. If variogram==True, returns variogram value instead. Requires (x,y,z)_1 and (x,y,z)_2.

Code author: Jared Deutsch - 2015-02-04

Get Pairwise Covariance

VariogramModel.pairwisecova(points)

Gets the pairwise covariance between points

Parameters:points (np.ndarray) – Should be a tidy array generated using datfl.locations.values or similar

Code author: Ryan Martin - 14-01-2018

Get Pointwise Covariance

VariogramModel.pointwisecova(point, points)

Gets the pointwise covariance between point and all points

Parameters:
  • point (np.ndarray) – A single point where point.shape[0] == points.shape[1]
  • points (np.ndarray) – Should be a tidy array generated using datfl[[datfl.x, datfl.y, datfl.z]].values or similar

Code author: Ryan Martin - 14-01-2018

Parse GSLIB-Style Model

VariogramModel.parsestr(vargstr=None)

Parses a GSLIB-style variogram model string

Code author: Jared Deutsch - 2015-02-04

Read Varfit Output

VariogramModel.readvarfit(varfitflname)

Reads a varfit output file and parses the variogram

Code author: Jared Deutsch - 2015-02-04

Read Varfit Model

VariogramModel.readvarmodel(varmodelflname)

Reads a varfit output file and parses the variogram

Code author: Jared Deutsch - 2015-02-04

Scale Variogram Model

VariogramModel.scalevar(variance)

Scales the variogram to the provided variance

Code author: Jared Deutsch - 2015-02-04

Set Rotation Matrix

VariogramModel.setcova()

Sets up the rotation matrix for cova3 (GSLIB setrot). This must be called before calling getcova

Code author: Jared Deutsch - 2015-02-04

Write Variogram Model

VariogramModel.write(outfile)

Write out the variogram in standard GSLIB format

Code author: Ryan Martin - 14-01-2018

Variogram Utility Functions

Discrete Gaussian Model

pygeostat.variograms.variogram_utils.dgm(evenvar, mean_x=None, var_x=None, numb_p=100, fcoef=None, accerr=1e-05, returnH=False)

Reads in a variogram file and returns a pygeostat data file or list of pygeostat DataFiles with the variograms

The only required parameters are evenvar and vargnum. All other parameters are optional.

The parameter vargnum can be one the following:

  • a single variogram number to return one pygeostat DataFile
  • a list of variogram numbers to return a list of DataFiles
  • the string ‘all’ to return a list with all variograms
Parameters:
  • evenvar[nd] – evenly spaced values of the variable to fit a DGM to: this means that if using declustering weights these must be evenly spaced samples of the declustered CDF!
  • vargnum – one of the the permissible vargnum parameters as specified above
  • mean_x – mean of variable
  • var_x – variance of variable
  • numb_p – number of Hermite polynomial coefficients to fit
  • accerr – acceptable error value for scaling
  • fcoef – equal to scale term if scaling the histogram
Returns:

fi coefficients used with Hermite polynomial terms

Return type:

fi[numb_p]

Returns:

values fit to evenly spaced quantiles

Return type:

fitvals[nd]

Returns:

Returned if required. Array of scaled values if fcoef not None

Return type:

scaledvals[nd]

Code author: Jared Deutsch - 2015-02-04

Get VarSim Parameters to Match Directions

pygeostat.variograms.variogram_utils.get_varsimpars(griddef, azm, dip=None, minstep=None, maxpct=None, maxstep=None, mute=False)

The program varsim, which is used to calculate experimental variograms of gridded data, requires the number of grid notes that must be shifted in each direction, equating to a lag distance. Finding the required shift values to mimic the variogram directions used to calculate the gridded data can be difficult. This function aims to determine some cell shift values (i.e., ixd, iyd, and izd) settings for varsim based on the desired azimuth and if required, dip.

Note

There is a newer 2-D version of this function that is more intuitive to use: gs.get_varsimpars_2d().

A maximum number of steps in each direction is calculated based on the passed gs.GridDef and the parameters maxpct and maxstep. The the number of cells determined by maxpct and the value of maxstep, which ever is smaller is used to cap the maximum number of steps in each direction. In other words, say there are 1000 cells in the x-direction, if maxpct is 0.05, then a maximum movement of 50 cells in that direction are allowed. If maxstep is set to 25, then that parameter will only allow moment in any directions by that many cells, as it is lower than the 50 previously determined. The default value of 10 for maxstep appears to give a reasonable range of possible directions. The default value of 5% for maxpct is to ensure that the final lag distance isn’t too large given the grid size, ensuring the experimental variogram calculated from the gridded data has an adequate resolution.

If you find you are getting enough resolution and the direction found is adequate (i.e., you want a larger lag distance), increase maxpct, maxstep, and minstep. If you require better precision in regards to getting the right azimuth, only increase maxpct and maxstep to allow as many different movements as possible.

For dip values greater than a value of |90| you may notice the returned izd value is indicating a dip opposite inclination of the direction specified. This is so the calculated ixd and iyd values do not have to be rotated. The calculations by varsim will still be completed on the input direction. For example, if a azimuth of 45° and a dip of -135° is passed, that is equivalent to an azimuth of 45° and a dip of 45° (for the purposes of varsim) which is what would be returned by the function.

Parameters:
  • griddef (GridDef) – A gs.GridDef containing the number of cells in each direction
  • azm (float) – Target azimuth value
  • dip (float) – Target dip value if required
  • minstep (int) – Minimum steps allowed in each direction
  • maxpct (float) – Maximum steps allowed in each direction as a percentage of the number of cells in the corresponding direction. Default value of 0.05 (i.e., 5%)
  • maxstep (int) – Maximum steps allowed in each direction. Default value of 25
  • mute (bool) – Indicate if printouts should be muted. Only recommended if looping the function. Printouts were deliberatively placed so the user can verify automation of variography.
Returns:

Azimuth value of the determined horizontal shift

Return type:

azm (float)

Returns:

Dip value of the determined horizontal and vertical shift

Return type:

dip (float or None)

Returns:

The determined shift in the x direction

Return type:

ixd (int)

Returns:

The determined shift in the y direction

Return type:

iyd (int)

Returns:

The determined shift in the z direction

Return type:

izd (int or None)

Code author: Warren E. Black - 2015-11-06

pygeostat.variograms.variogram_utils.get_varsimpars_2d(azm, ahmax, griddef, minlags=40, maxlags=200, degerr=0.5, vargrng=None, optlags=None, rngmult=1.25)

Warning

This is a 2-D Version of gs.get_varsimpars(). That function is wired and doesn’t work well. So far, have only recoded for 2D as I can’t think of a nice way to optimize 3D.

The program varsim, which is used to calculate experimental variograms of gridded data, requires the number of grid notes that must be shifted in each direction, equating to a lag distance. Finding the required shift values to mimic the variogram directions used to calculate the gridded data can be difficult. This function aims to determine some cell shift values (i.e., ixd, iyd) settings for varsim based on the desired azimuth.

The parameter degerr is the allowable absolute deviation from the major direction. For example, if the major direction is 20 degrees, a idx/idy configuration that results in a azimuth in the range of [19.5, 20.5] degrees is accepted.

Within this function, finding an optimal direction takes precedence over an optimal number of lags. The parameter optlags can be passed to specify an ideal number of lags to calculate. If a match is found within the allowable direction deviation range - as described above - it will be used.

If a final variogram range is desired the parameter vargrng can be passed. This will override the use of the variogram models maximum range for all directions. The parameter rngmult will also not be utilized.

Note

Assumes minor is perpendicular to major. Optimizes based on the major range only.

Parameters:
  • varmodel (gs.VariogramModel) – A gs.VariogramModel() class
  • griddef (GridDef) – A gs.GridDef class
  • minlags (int) – Minimum allowed lags that will cover the range of the major variogram range multiplied by the parameter rngmult
  • maxlags (int) – Maximum allowed lags that will cover the range of the major variogram range multiplied by the parameter rngmult
  • degerr (float) – Allowable absolute deviation from the major direction. See above for more detail
  • vargrng (float) – The desired final calculated range. Will not be multiples by the parameter rngmult and replaces the ahmax value pulled from the variogram model. Intended use is for when a variogram is going to be displayed on a plot that extends significantly past its maximum range
  • optlags (int) – Optimal number of lags that will cover the range of the major variogram range multiplied by the parameter rngmult. Only considered if there are values that fall within the range defined by degerr
  • rngmult (float) – Major direction range multiplier which is used to scale the major variogram range, which is the approximate maximum lag calculated. Only used if optdist is not used.
Returns:

  • major (tuple) (Tuple containing (ixd, idy) for the major direction)
  • minor (tuple) (Tuple containing (ixd, idy) for the minor direction)
  • majazm (float) (The real azimuth calculated in the major direction)
  • minazm (float) (The real azimuth calculated in the minor direction)
  • lag (int) (The lag distance covered by (ixd, idy))
  • nlag (int) (The number of lags to covered the desired distance)

Code author: Warren E. Black - 2016-04-01

Legacy Output

pygeostat.variograms.variogram_utils.legacyout(data, tail, head, direction, hdir, vdir, lags, vari, outfl, append=False)

Write a varcalc output DataFrame as a legacy gamv2004 and export it to a file.

Note

Will only except one variogram at a time. Designed for use within the pygeostat variogram class therefore it may not work outside.

Parameters:
  • tail (str) – Column header of the tail variable
  • head (str) – Column header of the head variable
  • direction (int) – Integer indicating the direction number
  • hdir (list) – A list containing the following parameters: [azm, azmtol, bandhorz]. Must be in that order.
  • vdir (list) – A list containing the following parameters: [dip, diptol, bandvert]. Must be in that order.
  • lags (list) – A list containing the following parameters: [nlags, lagdist, lagtol]. Must be in that order.
  • vari (list) – A list containing the following parameters: [variotypes, variotails, varioheads]. Must be in that order.
  • outfl (str) – Output file

Code author: Warren E. Black - 2016-02-22

Legacy Input

pygeostat.variograms.variogram_utils.legacyvarginput(file)

Read in legacy variogram output into a dataframe suitable for experimental data

Fully assumes the structure of the file, and assumes that ‘variogram’ is found in lines separating the different variogram models in the legacy file.

Parameters:file (str) – the variogram file to read in
Returns:
a dictionary of variogram models suitable for exp fitting with the
variogram class
Return type:variograms (dict)

Code author: Ryan Martin - 2017-02-24

Parse varfit_lmc Output File

pygeostat.variograms.variogram_utils.parseslmcvar(varfitfl)

Parse the output varfit_lmc.var file generate by the GSLIB varfit_lmc program.

Returns a list containing nested lists which contain the two variables the variogram model represents and the variogram string.

Example

A simple call:

>>> varfitfl = '../varfit_lmc.var'
>>> varstrs = gs.parselmcvar(varfitfl)

Print out the first variogram string parsed which contains the variable ID(s) which relate to the respective variogram model. In this example, the direct variogram model for the first variable is shown.

>>> print(varstrs[0][0])
    [1, 1]
>>> print(varstrs[0][1])
    2    0.020                               -  nst, nugget effect
    1    0.446   0.0   0.0   0.0             -  it,cc,ang1,ang2,ang3
         20000.00    20000.00    20000.00    -  a_hmax,a_hmin,a_vert
    2    0.533   0.0   0.0   0.0             -  it,cc,ang1,ang2,ang3
         45000.00    45000.00    45000.00    -  a_hmax,a_hmin,a_vert
Parameters:varfitfl – Location of the file containing the variogram model(s) fit by varfit_lmc
Returns:List containing nested lists with the variogram
Return type:varstrs (list)

Code author: Warren E. Black - 2016-02-26

Read Variogram Output

pygeostat.variograms.variogram_utils.readvarg(vargflname, vargnum)

Reads in a variogram file and returns a pygeostat data file or list of pygeostat DataFiles with the variograms

The parameter vargnum can be one of the following:

  • a single variogram number to return one pygeostat DataFile
  • a list of variogram numbers to return a list of DataFiles
  • the string ‘all’ to return a list with all variograms

The returned gs.DataFrame objects contain the following parameters:

  • data[‘Number’] - lag number
  • data[‘Distance’] - lag distance
  • data[‘Value’] - variogram value
  • data[‘Points’] - number of points
  • data[‘Head’] - head value of variogram
  • data[‘Tail’] - tail value of variogram
  • tailname - tail variable name
  • headname - head variable name
  • direction - direction number

or for gamv2004 variograms:

  • azm - azimuth
  • azmtol - azimuth tolerance
  • azmbandw - azimuth bandwidth
  • dip - dip
  • diptol - dip tolerance
  • dipbandw - dip bandwidth
  • lags - number of lags
  • lagdist - lag distance
  • lagtol - lag tolerance
  • vartype - variogram type
  • tailvar - tail variable number
  • headvar - head variable number
  • indcatcut - indicator category or cutoff (None if missing)
Parameters:
  • vargflname – name of variogram file in either gamv or gamv2004 style
  • vargnum – permissible vargnum parameter as specified above
Returns:

A gs.DataFile (or list of DataFiles) with data set to the variogram data as

as specified above

Return type:

vargs

Code author: Jared Deutsch - 2015-02-04

VarCalc Wrapper

pygeostat.variograms.variogram_utils.varcalc(datafl, variotypes, variotails, varioheads, azm, dip, tilt, azmtol, diptol, nlags, lagdist, lagtol, bandhorz, bandvert, variostd=False, variosills=None, strict=True, variocutcat=None)

Wrapper for the varcalc subroutine for calculating an experimental variogram

Variogram types are the same as in GSLIB:

  1. traditional semivariogram
  2. traditional cross semivariogram
  3. covariance (-3 calculates variance (provided sill) -covariance)
  4. correlogram (-4 calculates 1-correlation)
  5. general relative semivariogram
  6. pairwise relative semivariogram
  7. semivariogram of logarithms
  8. semimadogram
  9. indicator semivariogram - continuous - requires a cutoff
  10. indicator semivariogram - categorical - requires a category

Code author: Jared Deutsch - 2015-02-04

VarSim Wrapper

pygeostat.variograms.variogram_utils.varsim(lithfl=None, lithcol=None, lithcode=None, datafl=None, nvar=None, varcols=None, tmin=None, tmax=1e+21, outvargs=None, outavg=None, griddef=None, nreal=None, ndir=None, nlags=None, dirs=None, stdsill=None, nvargs=None, vargpars=None, dbg_fl=False)

Wrapper function that creates a text parameter file and passes it to the GSLIB varsim program that has been wrapped itself for use within python.

A temporary parameter file is written which is used by varsim; after which, it is deleted. No information is passed between python and varsim. Varsim is programed to only open ‘varsim.par’, which is generated by this function.

If the parameters outvargs and/or outavg are left at their default values of None, data will be returned as a pandas dataframe instead of writing to a file.

If a file name is passed to the dbg_fl parameter then a copy of the parameter file will written out using the file name

Warning

If any errors are encountered by varsim, python will crash. You can check the console python is running in to see any errors reported by varsim. A separate instance of python is opened so that if a crash occurs, the instance of python the user is currently using is protected.

The parameter file used by GSLIB varsim:

START OF PARAMETERS:
../data/lithology.dat        -file with lithology information
1   7                        -   lithology column (0=not used), code
../data/true.dat             -file with data
2   1   2                    -   number of variables, column numbers
-1.0e21     1.0e21           -   trimming limits
varsim_reals.out             -output file for variograms of realizations
varsim_avg.out               -output file for average variogram
50   0.5   1.0               -nx, xmn, xsiz
50   0.5   1.0               -ny, ymn, ysiz
 1   0.5   1.0               -nz, zmn, zsiz
 100                         -number of realizations
2  10                        -number of directions, number of lags
 1  0  0                     -ixd(1),iyd(1),izd(1)
 0  1  0                     -ixd(2),iyd(2),izd(2)
1                            -standardize sill? (0=no, 1=yes)
5                            -number of variograms
1   1   1                    -tail variable, head variable, variogram type
1   1   3                    -tail variable, head variable, variogram type
2   2   1                    -tail variable, head variable, variogram type
2   2   3                    -tail variable, head variable, variogram type
1   1   9  2.5               -tail variable, head variable, variogram type


type 1 = traditional semivariogram
     2 = traditional cross semivariogram
     3 = covariance
     4 = correlogram
     5 = general relative semivariogram
     6 = pairwise relative semivariogram
     7 = semivariogram of logarithms
     8 = semimadogram
     9 = indicator semivariogram - continuous
     10= indicator semivariogram - categorical

How the above parameter file relates to the variables used by this function:

START OF PARAMETERS:
{lithfl}                -file with lithology information
{lithcol}  {lithcode}   -   lithology column (0=not used), code
{datafl}                -file with data
{nvar}  {varcols}       -   number of variables, column numbers
{tmin}  {tmax}          -   trimming limits
{outvargs}              -output file for variograms of realizations
{outavg}                -output file for average variogram
{griddef}               -Full three line griddef string
{nreal}                 -number of realizations
{ndir}  {nlags}         -number of directions, number of lags
{dirs}                  -Lines of ixd(1),iyd(1),izd(1) parameters
{stdsill}               -standardize sill? (0=no, 1=yes)
{nvargs}                -number of variograms
{vargpars}              -Lines of tail variable, head variable, variogram type parameters

Code author: Warren E. Black - 2015-11-26 and Ryan M. Barnett - 2018-04-06

Calculate VarSim Lag Distance

pygeostat.variograms.variogram_utils.varsimlagdist(griddef, ixd, iyd=None, izd=None)

Calculate the lag distance of a set of cells shifts

Code author: Warren E. Black - 2015-11-06

Varvisual Wrapper

pygeostat.variograms.variogram_utils.varvisual(vtkfl, xyzorigin, azm, dip, tilt, azmtol, diptol, nlags, lagdist, lagtol, bandhorz, bandvert)

Wrapper for the varvisual subroutine for visualizing an experimental variogram definition

Code author: Jared Deutsch - 2015-02-04