search

gimli-org / gimli

Geophysical Inversion and Modeling Library :earth_africa:
https://www.pygimli.org
Other
344 stars 131 forks source link

Challenges with ERT field data import - question #680

Closed NoteboomM closed 1 month ago

NoteboomM commented 1 month ago

Problem description

I've had some success with ERT import and inversion when I have Res2dInv formatted .dat files, but I'm hoping to be able to use pygimli as an alternative. The current source I have from the field is in .stg files from a Supersting instrument, but import is not working properly, and when I convert to what (to me) looks like the unified data format, import doesn't work at all.

Your environment

  Date: Thu Mar 21 12:13:45 2024 W. Europe Standard Time

                OS : Windows
            CPU(s) : 20
           Machine : AMD64
      Architecture : 64bit
               RAM : 15.6 GiB
       Environment : IPython

  Python 3.9.18 | packaged by conda-forge | (main, Aug 30 2023, 03:40:31) [MSC
  v.1929 64 bit (AMD64)]

           pygimli : 1.5.0
            pgcore : 1.5.0
             numpy : 1.26.4
        matplotlib : 3.8.0
             scipy : 1.11.2
           IPython : 8.18.1
           pyvista : 0.43.4
--------------------------------------------------------------------------------

Steps to reproduce

The code is probably not very helpful, as the issue is simply with:

data = ert.load("DD4384(1-3).ohm") # Load data to container

However, sample input may be helpful? My attempt at the UDF-ERT looks like this:

890
# a b m n r ip/ms err/%
0    20   40   60   -5.01792E-01  1.20190E-02   5.00000E-01
0    20   60   80   -9.18413E-02  1.64594E-02   2.00000E-01
0    20   80   100  -4.80586E-02  7.63345E-03   5.00000E-01
0    20   100  120  -2.48010E-02  -2.37149E-02  2.00000E-01
0    20   120  140  -1.30468E-02  4.77906E-02   4.00000E-01
0    20   140  160  -1.14125E-02  1.52930E-02   1.50000E+00
0    20   160  200  -1.77302E-02  -1.72300E-02  1.20000E+00
0    20   200  260  -2.24212E-02  -7.66751E-02  8.00000E-01
0    40   120  140  -3.44304E-02  2.21064E-02   5.00000E-01
0    40   140  180  -4.43341E-02  1.07752E-02   3.00000E-01
0    40   180  240  -5.13266E-02  -5.07978E-03  2.00000E-01
0    40   240  320  -2.62669E-02  -1.15708E-02  3.00000E-01
0    40   320  420  -2.35751E-02  1.33933E-02   9.00000E-01
0    40   420  540  -1.65428E-02  2.50861E-02   4.00000E-01
0    40   540  680  -9.74725E-03  -8.71220E-02  5.00000E-01
20   40   60   80   -3.15258E-01  1.02085E-02   2.00000E-01
20   40   80   100  -8.73706E-02  8.93586E-03   1.00000E-01
20   40   100  120  -3.66347E-02  -1.18161E-02  3.00000E-01
20   40   120  140  -1.88496E-02  2.97665E-02   4.00000E-01
20   40   140  160  -1.44397E-02  1.15169E-02   4.00000E-01
20   40   160  180  -1.08570E-02  -3.75377E-03  4.00000E-01
...

and the supersting file looks like this:

Advanced Geosciences, Inc. SuperSting R8-IP Resistivity meter. S/N: SS1404154 Type: 3D
Firmware version: 01.33.78E Survey period: 20230508 Records: 890
Unit: meter
   1,USER   ,20230508,10:06:31,-5.01792E-01,   0,321, 1.89171E+02,DD4384   , 0.00000E+00, 0.00000E+00, 0.00000E+00, 2.00000E+01, 0.00000E+00, 0.00000E+00, 4.00000E+01, 0.00000E+00, 0.00000E+00, 6.00000E+01, 0.00000E+00, 0.00000E+00,IP:, 0, 2000, 3.75790E-03, 2.38829E-03, 1.86652E-03, 1.55905E-03, 1.29818E-03, 1.14909E-03, 1.20190E-02,Cmd=1,HV=400,Cyk=2,MTime=8.0,Gain=10,Ch=1
   2,USER   ,20230508,10:06:31,-9.18413E-02,   2,321, 1.38493E+02,DD4384   , 0.00000E+00, 0.00000E+00, 0.00000E+00, 2.00000E+01, 0.00000E+00, 0.00000E+00, 6.00000E+01, 0.00000E+00, 0.00000E+00, 8.00000E+01, 0.00000E+00, 0.00000E+00,IP:, 0, 2000, 8.52673E-03, 2.31620E-03, 1.84950E-03, 1.57814E-03, 1.11984E-03, 1.06902E-03, 1.64594E-02,Cmd=1,HV=400,Cyk=2,MTime=8.0,Gain=20,Ch=2
   3,USER   ,20230508,10:06:31,-4.80586E-02,   0,321, 1.81177E+02,DD4384   , 0.00000E+00, 0.00000E+00, 0.00000E+00, 2.00000E+01, 0.00000E+00, 0.00000E+00, 8.00000E+01, 0.00000E+00, 0.00000E+00, 1.00000E+02, 0.00000E+00, 0.00000E+00,IP:, 0, 2000, 3.45029E-03, 1.90674E-03, 1.07660E-03, 6.03139E-04, 3.43743E-04, 2.52934E-04, 7.63345E-03,Cmd=1,HV=400,Cyk=2,MTime=8.0,Gain=50,Ch=3
   4,USER   ,20230508,10:06:31,-2.48010E-02,   2,321, 1.86995E+02,DD4384   , 0.00000E+00, 0.00000E+00, 0.00000E+00, 2.00000E+01, 0.00000E+00, 0.00000E+00, 1.00000E+02, 0.00000E+00, 0.00000E+00, 1.20000E+02, 0.00000E+00, 0.00000E+00,IP:, 0, 2000,-5.37887E-03,-3.98394E-03,-4.07198E-03,-3.74496E-03,-3.51904E-03,-3.01613E-03,-2.37149E-02,Cmd=1,HV=400,Cyk=2,MTime=8.0,Gain=50,Ch=4
   5,USER   ,20230508,10:06:31,-1.30468E-02,   4,321, 1.72149E+02,DD4384   , 0.00000E+00, 0.00000E+00, 0.00000E+00, 2.00000E+01, 0.00000E+00, 0.00000E+00, 1.20000E+02, 0.00000E+00, 0.00000E+00, 1.40000E+02, 0.00000E+00, 0.00000E+00,IP:, 0, 2000, 1.02227E-02, 9.16704E-03, 8.85518E-03, 7.67883E-03, 6.34765E-03, 5.51924E-03, 4.77906E-02,Cmd=1,HV=400,Cyk=2,MTime=8.0,Gain=100,Ch=5
...

Expected behavior

I'm trying to produce a data container that I can use for inversion, but what I'm getting has issues

Actual behavior

From the file that I think is UDF-ERT, I get the following response:

21/03/24 - 12:20:03 - pyGIMLi - INFO - could not read unified data format for ERT ... try res2dinv
21/03/24 - 12:20:03 - pyGIMLi - INFO - could not read res2dinv ... try Ascii columns
21/03/24 - 12:20:03 - pyGIMLi - INFO - Failed importing Ascii column file. Consider using pybert.
21/03/24 - 12:20:03 - pyGIMLi - INFO - No electrode positions found!
21/03/24 - 12:20:04 - pyGIMLi - INFO - imported:  Data: Electrodes: 0 data: 0
760  780  800  820  -4.03606E+01  -2.48530E-02  2.00000E-01

followed by a wave of:

./core/src/datacontainer.cpp:281        virtual int GIMLI::DataContainer::load(const std::string&, bool, bool)  Warning! format description unknown: format[0] = a column ignored.
./core/src/datacontainer.cpp:281        virtual int GIMLI::DataContainer::load(const std::string&, bool, bool)  Warning! format description unknown: format[1] = b column ignored.
...

for every line/column combination, and then a wave of hundreds of:

Warning: Duplicated sensor position found at: 0 0   0

and at the end, there is no data container. Even if this is not proper UDF-ERT, it is clearly ASCII in columns...but doesn't read as that either.

If I try to import the supersting file, I get this response:

21/03/24 - 12:24:25 - pyGIMLi - INFO - could not read unified data format for ERT ... try res2dinv
21/03/24 - 12:24:25 - pyGIMLi - INFO - could not read res2dinv ... try Ascii columns
21/03/24 - 12:24:25 - pyGIMLi - INFO - Failed importing Ascii column file. Consider using pybert.
21/03/24 - 12:24:25 - pyGIMLi - INFO - No electrode positions found!
21/03/24 - 12:24:25 - pyGIMLi - INFO - imported:  Data: Electrodes: 42 data: 889

Data validity check: found 1 invalid data. 
Data validity check: remove invalid data.
Data validity check: found 1 invalid data. 
Data validity check: remove invalid data.

and I then have a data container that appears to have the electrode positions and resistivity, but if I do a print(data), I don't get the information about which valid columns are there, and I can't do ert.show(data['??']) for any parameter except rhoa.

halbmy commented 1 month ago

I would always prefer original instrument files as they bear more information (e.g. errors, current and voltage and not just apparent resistivity like in res2dinv). About the SuperSting import, there is an import routine in the pyBERT project (https://gitlab.com/resistivity-net/bert): https://gitlab.com/resistivity-net/bert/-/blob/master/python/pybert/importer/importData.py?ref_type=heads#L1845 Please try pb.importData('myfile.stg) and report any error by raising an issue in the BERT project, attaching your file, so that we can fix it.

halbmy commented 1 month ago

The problem in your file is the following: the unified BERT/GIMLi format consists of a list of electrodes and then the ABMN columns are just indices into the electrode list, so instead of 

890
# a b m n r ip/ms err/%
0    20   40   60   -5.01792E-01  1.20190E-02   5.00000E-01
0    20   60   80   -9.18413E-02  1.64594E-02   2.00000E-01
0    20   80   100  -4.80586E-02  7.63345E-03   5.00000E-01

you would have

6
#x y 
0 0
20 0
40 0
60 0
80 0
100 0
890
# a b m n r ip/ms err/%
1    2   3   4   -5.01792E-01  1.20190E-02   5.00000E-01
1    2   4   5   -9.18413E-02  1.64594E-02   2.00000E-01
1    3   5   6  -4.80586E-02  7.63345E-03   5.00000E-01
halbmy commented 1 month ago

I copied the specified 8 lines into a file test.stg and read it successfully with

import pybert as pb
data = pb.importData("test.stg")
print(data)

with the result

Data: Sensors: 8 data: 5, nonzero entries: ['a', 'b', 'err', 'i', 'ip', 'ip1', 'ip2', 'ip3', 'ip4', 'ip5', 'ip6', 'm', 'n', 'rhoa', 'u', 'valid']
NoteboomM commented 1 month ago

Thanks @halbmy. After a break for lunch, I re-read the BERT site and worked out what I had wrong there. And I agree, it's always good to go to the original where possible, but in my script to make my own BERT file, I carried over a few of those extra fields.

After fixing the BERT format, and with using pybert for the stg file, both import and seem to work with anything I try, thanks!

But maybe it's a separate issue, for both import approaches, I still only get this from print(data):

Data: Electrodes: 42 data: 889

I don't know that it is critical for other functions, but it's useful to have that list of nonzero fields.

halbmy commented 1 month ago

Can you attach your stg file?

NoteboomM commented 1 month ago

Here it is (zipped as github rejected the stg file). I've been trying a few more things, and different approaches/inversion settings, and everything seems fine with it practically speaking.

DD4384(1-3).zip

halbmy commented 1 month ago

I did

import pybert as pb
data = pb.importData("DD4384(1-3).stg")
print(data, data.tokenList()) 
data.show()

and it works well

Data: Electrodes: 42 data: 889 SensorIdx: a b m n  Data: err i ip ip1 ip2 ip3 ip4 ip5 ip6 iperr k r rhoa u valid

image

halbmy commented 1 month ago

I realized that pyBERT was overwriting the DataContainerERT.__repr__ function so that the included tokens were not shown. This is now corrected in the pyBERT code.

NoteboomM commented 1 month ago

Thanks Thomas - can I pull that update from git or somewhere? I think I installed pybert from an anaconda prompt.

halbmy commented 1 month ago

You don't really need this update if you print(data.tokenList()) which was just not included in print(data) when importing from pyBERT. Just a detail. But of course it is generally better to update pyBERT by git instead of conda (https://gitlab.com/resistivity-net/bert)