Inversion with topography
In this notebook, we will learn how to add topography and perform 2D
inversion with rectangular and triangular meshes. The files needed can
be found in examples/dc-2d-topo/.
import warnings
warnings.filterwarnings('ignore') # just to make it cleaner in the notebook
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-topo/'
import numpy as np # this will be used to read the topography file
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"
First, let’s create an Project object and import the survey as
usual.
k = Project(typ='R2') # create new Project object and use default working directory
k.createSurvey(os.path.join(testdir, 'syscal.csv'))
Working directory is: /builds/hkex/resipy/src/resipy
clearing dirname
600/636 reciprocal measurements found.
10 measurements error > 20 %
We can also plot the pseudo-section. Note that this one remains flat event when there is topogragraphy involved.
k.showPseudo()
If we take a look at the electrodes position by plotting R2.elec we
can see that there is no topography yet (all the values in the third
column (=z) are 0.0).
k.elec
| x | y | z | remote | buried | label | |
|---|---|---|---|---|---|---|
| 0 | 0.00 | 0.0 | 0.0 | False | False | 1 |
| 1 | 0.25 | 0.0 | 0.0 | False | False | 2 |
| 2 | 0.50 | 0.0 | 0.0 | False | False | 3 |
| 3 | 0.75 | 0.0 | 0.0 | False | False | 4 |
| 4 | 1.00 | 0.0 | 0.0 | False | False | 5 |
| 5 | 1.25 | 0.0 | 0.0 | False | False | 6 |
| 6 | 1.50 | 0.0 | 0.0 | False | False | 7 |
| 7 | 1.75 | 0.0 | 0.0 | False | False | 8 |
| 8 | 2.00 | 0.0 | 0.0 | False | False | 9 |
| 9 | 2.25 | 0.0 | 0.0 | False | False | 10 |
| 10 | 2.50 | 0.0 | 0.0 | False | False | 11 |
| 11 | 2.75 | 0.0 | 0.0 | False | False | 12 |
| 12 | 3.00 | 0.0 | 0.0 | False | False | 13 |
| 13 | 3.25 | 0.0 | 0.0 | False | False | 14 |
| 14 | 3.50 | 0.0 | 0.0 | False | False | 15 |
| 15 | 3.75 | 0.0 | 0.0 | False | False | 16 |
| 16 | 4.00 | 0.0 | 0.0 | False | False | 17 |
| 17 | 4.25 | 0.0 | 0.0 | False | False | 18 |
| 18 | 4.50 | 0.0 | 0.0 | False | False | 19 |
| 19 | 4.75 | 0.0 | 0.0 | False | False | 20 |
| 20 | 5.00 | 0.0 | 0.0 | False | False | 21 |
| 21 | 5.25 | 0.0 | 0.0 | False | False | 22 |
| 22 | 5.50 | 0.0 | 0.0 | False | False | 23 |
| 23 | 5.75 | 0.0 | 0.0 | False | False | 24 |
Then we can load a csv file that will add topography for each electrode. The csv file should be of the form X,Y,Z and no headers are needed.
k.importElec(os.path.join(testdir + 'elec.csv'))
print(k.elec)
x y z remote buried label
0 0.00 0.0 29.499 False False 1
1 0.25 0.0 29.504 False False 2
2 0.50 0.0 29.509 False False 3
3 0.75 0.0 29.516 False False 4
4 1.00 0.0 29.478 False False 5
5 1.25 0.0 29.461 False False 6
6 1.50 0.0 29.454 False False 7
7 1.75 0.0 29.428 False False 8
8 2.00 0.0 29.416 False False 9
9 2.25 0.0 29.411 False False 10
10 2.50 0.0 29.398 False False 11
11 2.75 0.0 29.362 False False 12
12 3.00 0.0 29.329 False False 13
13 3.25 0.0 29.245 False False 14
14 3.50 0.0 29.159 False False 15
15 3.75 0.0 29.083 False False 16
16 4.00 0.0 29.010 False False 17
17 4.25 0.0 28.929 False False 18
18 4.50 0.0 28.872 False False 19
19 4.75 0.0 28.761 False False 20
20 5.00 0.0 28.672 False False 21
21 5.25 0.0 28.593 False False 22
22 5.50 0.0 28.509 False False 23
23 5.75 0.0 28.412 False False 24
We can now create a mesh either triangular or rectangular.
k.createMesh(typ='trian') # or trian for triangular mesh
k.showMesh()
Creating triangular mesh...done (2356 elements)
We can finally invert the data on this mesh. Note that it might take a while so be patient.
k.invert()
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: 336 Measurements rejected: 0
Geometric mean of apparent resistivities: 0.35349E+03
>> Total Memory required is: 0.032 Gb
Iteration 1
Initial RMS Misfit: 33.53 Number of data ignored: 0
Alpha: 1224.408 RMS Misfit: 4.57 Roughness: 7.952
Alpha: 568.320 RMS Misfit: 4.77 Roughness: 10.106
Step length set to 1.00000
Final RMS Misfit: 4.57
Updated data weights
Iteration 2
Initial RMS Misfit: 3.32 Number of data ignored: 0
Alpha: 734.014 RMS Misfit: 2.47 Roughness: 7.520
Alpha: 340.699 RMS Misfit: 2.03 Roughness: 8.829
Alpha: 158.139 RMS Misfit: 1.77 Roughness: 10.252
Alpha: 73.401 RMS Misfit: 1.62 Roughness: 11.911
Alpha: 34.070 RMS Misfit: 1.55 Roughness: 13.778
Alpha: 15.814 RMS Misfit: 1.52 Roughness: 15.765
Alpha: 7.340 RMS Misfit: 1.52 Roughness: 18.023
Alpha: 3.407 RMS Misfit: 1.54 Roughness: 21.352
Step length set to 1.00000
Final RMS Misfit: 1.52
Attempted to update data weights and caused overshoot
treating as converged
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
k.showResults(contour=True, sens=True, zlim=[28,29.6]) # with contour
In a nutshell
k = Project(typ='R2')
k.createSurvey(testdir + 'syscal.csv')
k.importElec(testdir + 'elec.csv')
k.invert()
k.showResults(contour=True, sens=True, zlim=[28,29.6])
Working directory is: /builds/hkex/resipy/src/resipy
clearing dirname
600/636 reciprocal measurements found.
10 measurements error > 20 %
Creating triangular mesh...done (2356 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: 336 Measurements rejected: 0
Geometric mean of apparent resistivities: 0.35349E+03
>> Total Memory required is: 0.032 Gb
Iteration 1
Initial RMS Misfit: 33.53 Number of data ignored: 0
Alpha: 1224.408 RMS Misfit: 4.57 Roughness: 7.952
Alpha: 568.320 RMS Misfit: 4.77 Roughness: 10.106
Step length set to 1.00000
Final RMS Misfit: 4.57
Updated data weights
Iteration 2
Initial RMS Misfit: 3.32 Number of data ignored: 0
Alpha: 734.014 RMS Misfit: 2.47 Roughness: 7.520
Alpha: 340.699 RMS Misfit: 2.03 Roughness: 8.829
Alpha: 158.139 RMS Misfit: 1.77 Roughness: 10.252
Alpha: 73.401 RMS Misfit: 1.62 Roughness: 11.911
Alpha: 34.070 RMS Misfit: 1.55 Roughness: 13.778
Alpha: 15.814 RMS Misfit: 1.52 Roughness: 15.765
Alpha: 7.340 RMS Misfit: 1.52 Roughness: 18.023
Alpha: 3.407 RMS Misfit: 1.54 Roughness: 21.352
Step length set to 1.00000
Final RMS Misfit: 1.52
Attempted to update data weights and caused overshoot
treating as converged
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
Download python script: nb_02dc-topo.py
Download Jupyter notebook: nb_02dc-topo.ipynb
View the notebook in the Jupyter nbviewer
Run this example interactively: 