R2 basic tutorial

In this tutorial you will learn how to use the Python API of R* codes (http://www.es.lancs.ac.uk/people/amb/Freeware/R2/R2.htm). Start by importing the Project master class from the API (Application Programming Interface).

1 Basics imports

Just import basic packages and the R2 API as a module (note : you will need to change the path for it, here we assume you launched the jupyter from inside the /examples/jupyter-notebook folder).

    import warnings
warnings.filterwarnings('ignore')
import os
import sys
sys.path.append((os.path.relpath('../src'))) # add here the relative path of the API folder
testdir = '../src/examples/dc-2d/'
from resipy import Project

API path =  /builds/hkex/resipy/src/resipy
ResIPy version =  3.4.1
cR2.exe found and up to date.
R3t.exe found and up to date.
cR3t.exe found and up to date.

it looks like wine32 is missing, you should install it.
multiarch needs to be enabled first.  as root, please
execute "dpkg --add-architecture i386 && apt-get update &&
apt-get install wine32"

2 Create an ‘Project’ object, import data and plot pseudo section

The Project class was referred to as R2 class in older version of ResIPy.

The first step is to create an object out of the Project class, let’s call it k . This is the main object we are going to interact with. Then the second step is to read the data from a survey file. Here we choose a csv file from the Syscal Pro that contains resistivity data only. Note then when importing the survey data, the object automatically search for reciprocal measurements and will compute a reciprocal error with the ones it finds.

k = Project(typ='R2') # create a Project object in a working directory (can also set using k.setwd())
k.createSurvey(testdir + 'syscal.csv', ftype='Syscal') # read the survey file

Working directory is: /builds/hkex/resipy/src/resipy
clearing dirname
308/344 reciprocal measurements found.
0 measurements error > 20 %

We can plot the pseudosection and display errors based on reciprocal measurements.

k.showPseudo()
k.showError() # plot the reciprocal errors
../_images/nb_01getting-started_0.png ../_images/nb_01getting-started_1.png

3 Data filtering

Below are a few examples of data filtering routines that can be used: - k.filterUnpaired() to remove unpaired measurements (so measurements with no reciprocal) -> those could be dummy measurements in a dipole-dipole configuration - k.filterElec([5]) to remove a specific electrode (e.g. here all quadrupoles with electrode 5 are deleted) - k.filterRecip(20) to remove measurements based on their relative reciprocal error (e.g. all quadrupoles with a reciprocal error > 20% are discarded). More advanced data filtering can be achieved using the filterData() method from the Survey class. This method allows to filter out specific quadrupoles. An interactive version of it can be access using the filterManual() method which produces an interactive pseudo-section in the UI.

k.filterUnpaired()
k.showPseudo() # this actually removed the dummy measurements in this dipole-dipole survey (added for speed optimization)

removeUnpaired:filterData: 36 / 344 quadrupoles removed.
../_images/nb_01getting-started_2.png 
k.filterElec([5]) # remove all quadrupoles associated with electrode 5
k.showPseudo()

filterData: 48 / 308 quadrupoles removed.
48 measurements removed!
../_images/nb_01getting-started_3.png 
k.filterRecip(percent=20) # in this case this only removes one quadrupoles with reciprocal error bigger than 20 percent
k.showPseudo()

0 measurements with greater than 20.0% reciprocal error removed!
../_images/nb_01getting-started_4.png

4 Fitting an error model

Different errors models are available to be fitted for DC data: - a simple linear model: k.fitErrorLin() - a power law model: k.fitErrorPwl() - a linear mixed effect model: k.fitErrorLME() (on Linux only with an R kernel installed) Each of those will create a new error column in the Survey object that will be used in the inversion if k.err = True.

k.fitErrorLin()

Error model is R_err = 0.006*R_avg + 0 (R^2 = 0.996)
../_images/nb_01getting-started_5.png 
k.fitErrorPwl()

Error model is R_err = 0.005 R_avg^1.127 (R^2 = 0.992)
../_images/nb_01getting-started_6.png

5 Mesh

Two types of mesh can be created in 2D: - a quadrilateral mesh (k.createMesh('quad')) - a triangular mesh (k.createmesh('trian')) For 3D, only tetrahedral mesh can be created using k.createMesh('tetra').

k.createMesh(typ='quad') # generate quadrilateral mesh (default for 2D survey)
k.showMesh()

Creating quadrilateral mesh...done (3306 elements)
../_images/nb_01getting-started_7.png 
k.createMesh('trian', show_output=False) # this actually call gmsh.exe to create the mesh
k.showMesh()

Creating triangular mesh...done (1770 elements)
../_images/nb_01getting-started_8.png

7 Inversion

The inversion takes place in the specify working directory of the R2 object specified the first time the k = R2(<workingDirectory>) is called. It can be changed after by using k.setwd(<newWorkingDirectory>). The parameters of the inversion are defined in a dictionnary in k.param and ca be changed manually by the user (e.g. k.param['a_wgt'] = 0.01. All parameters have a default values and their names follow the R2 manual. The .in file is written automatically when the k.invert() method is called.

k.param['data_type'] = 1 # using log of resistitivy
k.err = True # if we want to use the error from the error models fitted before
k.invert() # this will do the inversion

Writing .in file and protocol.dat... done

--------------------- MAIN INVERSION ------------------

it looks like wine32 is missing, you should install it.
multiarch needs to be enabled first.  as root, please
execute "dpkg --add-architecture i386 && apt-get update &&
apt-get install wine32"

 >> R  2    R e s i s t i v i t y   I n v e r s i o n   v4.02 <<

 >> D a t e : 25 - 10 - 2022
 >> My beautiful survey
 >> I n v e r s e   S o l u t i o n   S e l e c t e d <<
 >> Determining storage needed for finite element conductance matrix
 >> Generating index array for finite element conductance matrix
 >> Reading start resistivity from res0.dat
 >> R e g u l a r i s e d   T y p e <<
 >>   L i n e a r    F i l t e r    <<
 >> L o g - D a t a   I n v e r s i o n <<
 >> N o r m a l   R e g u l a r i s a t i o n <<
 >> D a t a   w e i g h t s   w i l l   b e  m o d i f i e d <<
 >> D a t a   w e i g h t   t o   b e   r e a d   f r o m   d a t a   f i l e <<


 Processing dataset   1


 Measurements read:   130     Measurements rejected:     0
   Geometric mean of apparent resistivities:  0.52480E+02

 >> Total Memory required is:          0.017 Gb

   Iteration   1
     Initial RMS Misfit:       125.18       Number of data ignored:     0
     Alpha:        4318.670   RMS Misfit:        7.62  Roughness:        1.663
     Alpha:        2004.549   RMS Misfit:        6.52  Roughness:        2.360
     Alpha:         930.429   RMS Misfit:        5.61  Roughness:        3.457
     Alpha:         431.867   RMS Misfit:        4.93  Roughness:        5.095
     Alpha:         200.455   RMS Misfit:        4.59  Roughness:        7.158
     Alpha:          93.043   RMS Misfit:        4.64  Roughness:        9.464
     Step length set to   1.00000
     Final RMS Misfit:        4.59
     Updated data weights

   Iteration   2
     Initial RMS Misfit:         3.45       Number of data ignored:     0
     Alpha:         114.845   RMS Misfit:        1.61  Roughness:        8.473
     Alpha:          53.306   RMS Misfit:        1.10  Roughness:       10.746
     Alpha:          24.743   RMS Misfit:        0.79  Roughness:       12.992
     Step length set to   1.00000
     Final RMS Misfit:        0.79
     Final RMS Misfit:        1.02

 Solution converged - Outputing results to file

 Calculating sensitivity map


 Processing dataset   2


 End of data:  Terminating
1/1 results parsed (1 ok; 0 failed)

All ok

8 Results visualisation and post-processing

Results can be show with k.showResults(). Multiple arguments can be passed to the method in order rescale the colorbar bar, view the sensitivity or not, change the attribute or plot contour. The errors from the inversion can also be plotted using either k.pseudoError() or k.showInvError().

k.showResults(attr='Resistivity(ohm.m)', sens=False, contour=True, vmin=30, vmax=100)
../_images/nb_01getting-started_9.png 
k.showPseudoInvError() # allow to see if some electrodes get higher error
../_images/nb_01getting-started_10.png 
k.showInvError() # all errors should be between -3 and 3
../_images/nb_01getting-started_11.png

In a nutshell

Below is a minimal example which imports the data, plots a pseudo section and inverts using all default parameters.

k = Project(typ='R2') # create an Project object in a working directory (can also set using k.setwd())
k.createSurvey(testdir + 'syscal.csv', ftype='Syscal') # read the survey file
k.showPseudo() # plot pseudo section
k.invert(iplot=True) # does the inversion (generate quand mesh and use default R2.in settings)

Working directory is: /builds/hkex/resipy/src/resipy
clearing dirname
308/344 reciprocal measurements found.
0 measurements error > 20 %
Creating triangular mesh...done (1928 elements)
Writing .in file and protocol.dat... done

--------------------- MAIN INVERSION ------------------

it looks like wine32 is missing, you should install it.
multiarch needs to be enabled first.  as root, please
execute "dpkg --add-architecture i386 && apt-get update &&
apt-get install wine32"

 >> R  2    R e s i s t i v i t y   I n v e r s i o n   v4.02 <<

 >> D a t e : 25 - 10 - 2022
 >> My beautiful survey
 >> I n v e r s e   S o l u t i o n   S e l e c t e d <<
 >> Determining storage needed for finite element conductance matrix
 >> Generating index array for finite element conductance matrix
 >> Reading start resistivity from res0.dat
 >> R e g u l a r i s e d   T y p e <<
 >>   L i n e a r    F i l t e r    <<
 >> L o g - D a t a   I n v e r s i o n <<
 >> N o r m a l   R e g u l a r i s a t i o n <<
 >> D a t a   w e i g h t s   w i l l   b e  m o d i f i e d <<


 Processing dataset   1


 Measurements read:   190     Measurements rejected:     0
   Geometric mean of apparent resistivities:  0.52987E+02

 >> Total Memory required is:          0.020 Gb

   Iteration   1
     Initial RMS Misfit:        28.72       Number of data ignored:     0
     Alpha:         408.075   RMS Misfit:        1.94  Roughness:        1.992
     Alpha:         189.412   RMS Misfit:        1.52  Roughness:        2.951
     Alpha:          87.917   RMS Misfit:        1.20  Roughness:        4.152
     Alpha:          40.807   RMS Misfit:        0.99  Roughness:        5.616
     Step length set to   1.00000
     Final RMS Misfit:        0.99

 Cannot fit quadratic through step lengths
     Final RMS Misfit:        0.99

 Solution converged - Outputing results to file

 Calculating sensitivity map


 Processing dataset   2


 End of data:  Terminating
1/1 results parsed (1 ok; 0 failed)

All ok
../_images/nb_01getting-started_12.png ../_images/nb_01getting-started_13.png

Download python script: nb_01getting-started.py

Download Jupyter notebook: nb_01getting-started.ipynb

View the notebook in the Jupyter nbviewer

Run this example interactively: binder