GeostatsPy

repository·master·Indexed 20 days ago

https://github.com/geostatsguy/geostatspy

An open-source Python package for spatial data analytics and geostatistics that provides Python-native reimplementations and wrappers of GSLIB (Geostatistical Library) functions. It includes modules for numerical geostatistical methods such as variograms, distribution transformations (nscore, affine), and spatial estimation/simulation (kb2d, sgsim, sisim), as well as utilities for data conversion between pandas DataFrames and GSLIB Geo-EAS files and matplotlib-based spatial visualizations.

Tokens
26.5K
Snippets
82
Records
103
Agent score
68%

What's inside geostatspy

  1. Understand the GeostatsPy package structure

    master

    The package is divided into two main modules designed to bring GSLIB (Geostatistical Library) functionality to Python:

    1. geostatspy.geostats: Contains GSLIB functions rewritten directly in Python. This includes variograms, distribution transformations, and spatial estimation/simulation methods (e.g., kb2d, sgsim, nscore).
    2. geostatspy.GSLIB: Contains reimplementations of GSLIB visualizations and "low-tech" wrappers. These wrappers write parameter files, execute the original GSLIB FORTRAN binaries, and read the results back into Python. Note: Using these wrappers requires the GSLIB executables to be present in your working directory.
  2. Configure the GeostatsPy environment and dependencies

    master

    GeostatsPy requires Python 3.10 - 3.13 and numpy < 2.2.x (due to numba compatibility).

    While many dependencies like numpy, pandas, scipy, and matplotlib are included in standard Anaconda distributions, you may need to install tqdm and statsmodels manually.

    If you encounter import errors, you can install missing packages using: python -m pip install [package-name]

    pip install tqdm
    pip install statsmodels
  3. Set up GSLIB executables for geostatspy

    master

    To use certain functions in geostatspy that interface with GSLIB, you must ensure the following GSLIB executables are present in your current working directory or in a directory listed in your system's PATH environment variable:

    • nscore.exe
    • declus.exe
    • gam.exe
    • gamv.exe
    • vmodel.exe
    • kb2d.exe
    • sgsim.exe

    Note: Some functions like nscore, gam, gamv, and kb2d have been rewritten in Python and are available directly within the geostatspy package, reducing the dependency on external executables for those specific tasks.

    # Example: Setting the working directory to ensure executables are found
    import os
    os.chdir("c:/PGE337")
  4. Install GeostatsPy via pip

    master

    GeostatsPy is available on PyPI. You can install it using pip. To ensure you have the most recent version or to fix potential installation issues, you can use the --force-reinstall flag.

    # Standard installation
    pip install geostatspy
    
    # Update to the most recent version
    pip install --upgrade --force-reinstall geostatspy
  5. Setup GeostatsPy and GSLIB requirements

    master

    GeostatsPy provides Python wrappers for the GSLIB (Geostatistical Library) executables. To use the core geostatistical functions, you must ensure the following GSLIB executables are present in your current working directory:

    • nscore.exe
    • declus.exe
    • gam.exe
    • gamv.exe
    • vmodel.exe
    • kb2d.exe
    • sgsim.exe

    Note on Operating Systems:

    Required Python Libraries:

    • numpy
    • pandas
    • matplotlib
    # Example of typical imports used in GeostatsPy workflows
    import pandas as pd
    import os
    import numpy as np
    import matplotlib.pyplot as plt
    import random as rand
  6. Configure Kriging types in sgsim()

    master

    The ktype parameter in sgsim determines the kriging method used to calculate the conditional mean and standard deviation during simulation. Supported types include:

    • 0: Simple Kriging
    • 1: Ordinary Kriging
    • 2: Locally Varying Mean (uses sec_map as a secondary variable)
    • 3: External Drift (uses secondary data as a drift component)
    • 4: Collocated Cosimulation (uses a secondary variable with a correlation factor colocorr)
  7. Setup GSLIB executables for geostatspy

    master

    To use certain functions in geostatspy that act as wrappers for GSLIB (like sgsim), you must have the GSLIB executables in your working directory.

    Required executables:

    • nscore.exe, declus.exe, gam.exe, gamv.exe, vmodel.exe, kb2d.exe, and sgsim.exe.

    Note: Many of these (like nscore, gam, gamv, and kb2d) have been rewritten in Python and are available directly in the geostatspy package, reducing the need for external executables. However, sgsim.exe is still required for the GSLIB_sgsim_2d_uncond wrapper.

    Downloads:

  8. Access GSLIB functionality in Python

    master

    The geostatspy package provides Python implementations of several GSLIB (Geostatistical Library) functions. Some functions are direct Python rewrites, while others may require GSLIB executables to be present in the working directory.

    Available Python-native functions/utilities:

    • nscore: Normal score transformation.
    • gam, gamv: Variogram calculation.
    • vmodel: Variogram modeling.
    • kb2d: 2D Kriging.
    • sgsim: Sequential Gaussian Simulation.
    • hist: Histogram calculation (reimplemented with MatPlotLib).
    • locmap: Location map (reimplemented with MatPlotLib).

    Note: For some legacy workflows, GSLIB executables (nscore.exe, declus.exe, gam.exe, gamv.exe, vmodel.exe, kb2d.exe, and sgsim.exe) must be located in the working directory.

  9. Configure GSLIB executables for GeostatsPy

    master

    Several GeostatsPy functions act as wrappers for GSLIB Fortran executables. To use functions like nscore, gam, gamv, vmodel, kb2d, or sgsim, the corresponding .exe files must be present in either your current working directory or in a directory included in your system's PATH environment variable.

    Required executables:

    • nscore.exe
    • declus.exe
    • gam.exe
    • gamv.exe
    • vmodel.exe
    • kb2d.exe
    • sgsim.exe

    Executables can be obtained from GSLIB.com. For Mac OS X users, compiled versions are available at https://github.com/GeostatsGuy/GSLIB_MacOS.

  10. Example: Declustering workflow with GeostatsPy

    master

    This example demonstrates creating a synthetic 2D simulation, introducing sampling bias, and using the declus method to calculate weights and debias the data.

    import geostatspy.GSLIB as GSLIB                          # GSLIB utilities, viz and wrapped functions
    import geostatspy.geostats as geostats                    # GSLIB converted to Python
    import matplotlib.pyplot as plt                           # plotting
    import scipy.stats                                        # summary stats of ndarrays
    
    # Make a 2d simulation
    nx = 100; ny = 100; cell_size = 10                        # grid number of cells and cell size
    xmin = 0.0; ymin = 0.0;                                   # grid origin
    xmax = xmin + nx * cell_size; ymax = ymin + ny * cell_size# calculate the extent of model
    seed = 74073                                              # random number seed  for stochastic simulation    
    range_max = 1800; range_min = 500; azimuth = 65           # Porosity variogram ranges and azimuth
    vario = GSLIB.make_variogram(0.0,nst=1,it1=1,cc1=1.0,azi1=65,hmaj1=1800,hmin1=500) # assume variogram model
    mean = 10.0; stdev = 2.0                                  # Porosity mean and standard deviation
    vmin = 4; vmax = 16; cmap = plt.cm.plasma                 # color min and max and using the plasma color map
    
    # calculate a stochastic realization with standard normal distribution
    sim = GSLIB.sgsim_uncond(1,nx,ny,cell_size,seed,vario,"simulation") # 2d unconditional simulation
    sim = GSLIB.affine(sim,mean,stdev)                        # correct the distribution to a target mean and standard deviation
    
    # extract samples from the 2D realization 
    sampling_ncell = 10  # sample every 10th node from the model
    samples = GSLIB.regular_sample(sim,xmin,xmax,ymin,ymax,sampling_ncell,10,10,nx,ny,'Realization')
    
    # remove samples to create a sample bias (preferentially removed low values to bias high)
    samples_cluster = samples.drop([80,79,78,73,72,71,70,65,64,63,61,57,56,54,53,47,45,42]) # this removes specific rows (samples)
    samples_cluster = samples_cluster.reset_index(drop=True)  # we reset and remove the index (it is not sequential anymore)
    GSLIB.locpix(sim,xmin,xmax,ymin,ymax,cell_size,vmin,vmax,samples_cluster,'X','Y','Realization','Porosity Realization and Regular Samples','X(m)','Y(m)','Porosity (%)',cmap,"Por_Samples")
    
    # apply the declus program convert to Python
    wts,cell_sizes,averages = geostats.declus(samples_cluster,'X','Y','Realization',iminmax=1,noff=5,ncell=100,cmin=1,cmax=2000)
    samples_cluster['wts'] = wts            # add the weights to the sample data
    samples_cluster.head()
    
    # plot the results and diagnostics for the declustering
    plt.subplot(321)
    GSLIB.locmap_st(samples_cluster,'X','Y','wts',xmin,xmax,ymin,ymax,0.0,2.0,'Declustering Weights','X (m)','Y (m)','Weights',cmap)
    
    plt.subplot(322)
    GSLIB.hist_st(samples_cluster['wts'],0.0,2.0,log=False,cumul=False,bins=20,weights=None,xlabel="Weights",title="Declustering Weights")
    plt.ylim(0.0,20)
    
    plt.subplot(323)
    GSLIB.hist_st(samples_cluster['Realization'],0.0,20.0,log=False,cumul=False,bins=20,weights=None,xlabel="Porosity",title="Naive Porosity")
    plt.ylim(0.0,20)
    
    plt.subplot(324)
    GSLIB.hist_st(samples_cluster['Realization'],0.0,20.0,log=False,cumul=False,bins=20,weights=samples_cluster['wts'],xlabel="Porosity",title="Naive Porosity")
    plt.ylim(0.0,20)
    
    # Plot the declustered mean vs. cell size to check the cell size selection
    plt.subplot(325)
    plt.scatter(cell_sizes,averages, c = "black", marker='o', alpha = 0.2, edgecolors = "none")
    plt.xlabel('Cell Size (m)')
    plt.ylabel('Porosity Average (%)')
    plt.title('Porosity Average vs. Cell Size')
    plt.ylim(8,12)
    plt.xlim(0,2000)
    
    print(scipy.stats.describe(wts))
    
    plt.subplots_adjust(left=0.0, bottom=0.0, right=2.0, top=3.5, wspace=0.2, hspace=0.2)
    plt.show()
  11. Visualize Kriging results with `locpix_st` and `pixelplt_st`

    master

    To visualize kriging estimates and variance maps, use the following plotting utilities:

    • locpix_st: Plots the kriging estimates (kmap) overlaid with the original data points.
    • pixelplt_st: Plots the kriging variance map (vmap) as a pixelated grid.

    Common arguments for both:

    • map_array: The 2D array to plot (kmap or vmap).
    • xmin, xmax, ymin, ymax: Spatial extent of the grid.
    • xsiz: Cell size.
    • vmin, vmax: Color scale limits.
    • df (for locpix_st only): The DataFrame containing the data points.
    • x, y, vr (for locpix_st only): Column names for data coordinates and values.
    • cmap: Matplotlib colormap (e.g., plt.cm.plasma).
    import matplotlib.pyplot as plt
    
    # Assuming kmap, vmap, xmin, xmax, ymin, ymax, xsiz, df, etc. are defined
    plt.subplot(121)
    locpix_st(kmap, xmin, xmax, ymin, ymax, xsiz, 0.0, 3.0, df, 'x', 'y', 'vr', 'Kriging Estimate', 'X(m)', 'Y(m)', 'Porosity (%)', cmap)
    
    plt.subplot(122)
    pixelplt_st(vmap, xmin, xmax, ymin, ymax, xsiz, 0.0, 1.0, 'Kriging Variance', 'X(m)', 'Y(m)', 'Porosity (%^2)', cmap)
    
    plt.show()