Run Gravity Inversion by using real field data

[Carlos Osmin Pocasangre Jimenez]

Hello,

I'm working on Gravity data from my study are so I want to run an inversion modeling. I have checked the Cookbook and all examples use synthetic data for running the inversion modeling. How can I use the fatiando package and my field data?

best regards

Carlos

[Leonardo Uieda]

Hi Carlos,

Yes you can. Everything works with field data. The examples all use synthetic because it's easier to show what the method does when we know the ground truth.

You'll have to load your data somehow.

I recommend using an ASCII xyz file and using the numpy.loadtxt function (https://docs.scipy.org/doc/numpy/reference/generated/numpy.loadtxt.html).

Leo

--
Leonardo Uieda

Visiting Researcher | University of Hawaii at Manoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

--
You received this message because you are subscribed to the Google Groups "Fatiando a Terra" group.
To unsubscribe from this group and stop receiving emails from it, send an email to fatiando+unsubscribe@googlegroups.com.
To post to this group, send email to fati...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Carlos Jimenez]

Hello Leo, 

Thanks for the rapid response. I wrote a code using based on your examples. However, I got some issues. First, the minimum RMS is 2.014, second, the MayaVI 3d plot is transposed. I copied a piece of code next. I imported the data as you said xyz.

if __name__ == "__main__": print '\n---- 3D GRAVITY INVERSION ----' # upload data xp, yp, gz = np.loadtxt('local.csv', comments='#', skiprows=1, delimiter=',', unpack=True) zp = np.empty(len(xp)) zp.fill(-1.0) xp = xp - min(xp) yp = yp - min(yp) bounds = [0, 11000, 0, 11000, 0, 4000] prism_shape = (20, 40, 40) mesh = mesher.PrismMesh(bounds, prism_shape) seeds = harvester.sow( [(4000, 8689, 2000, {'density': 400}), (4181.91, 7825, 2000, {'density': -400}), (2349.18, 6725.36, 2000, {'density': 400}), (935.351, 4473.71, 2000, {'density': -400}), (3972.46, 1488.96, 2000, {'density': 400}), (6590.65, 494.051, 2000, {'density': -400}), (10177.6, 7458.45, 2000, {'density': -400}), (6067.01, 3190.79, 2000, {'density': 400}), (7009.56, 1881.69, 2000, {'density': 400}), (7978.3, 3374.07, 2000, {'density': -400}), (4836.46, 5966.08, 2000, {'density': -400}), (6695.38, 5547.17, 2000, {'density': 400}), (8449.57, 6227.9, 2000, {'density': -400}), (9627.76, 8165.36, 2000, {'density': 400}), (8632.85, 9212.64, 2000, {'density': -400})], mesh) data = [harvester.Gz(x, y, z, g_obs)]

estimate, predicted = harvester.harvest(data, seeds, mesh, compactness=0.08, threshold=0.00001) mesh.addprop('density', estimate['density']) bodies = mesher.vremove(0, 'density', cube_mesh)

[Leonardo Uieda]

Hi Carlos,

Your code got a bit mangled in the email.

Please send the script and some screenshots of the problem.

It's kind of difficult to diagnose without seeing the output and the full code.

--
Leonardo Uieda

Visiting Researcher | University of Hawaii at Manoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

[Filippo Spadavechia]

leonardo do you have a scrip to visualize sgy files using python 3.5 ?

[Leonardo Uieda]

Hi Filippo,

We don't but maybe this can do the job: https://github.com/Statoil/segyviewer

--
Leonardo Uieda

Visiting Research Scholar | University of Hawaii at Manoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

[Filippo Spadavechia]

ok let me try it! thanks! 

[accep handyarso]

Hello Leo,,

I am very interested with your Planting Seeds algorithm. I have installed fatiando and tried to run the Gravity Inversion scripts but I get stuck with importing "vremove", could you please to help??

I don't know what missing, perhaps some scripts but I am not sure

thanks for your kindly help

regards

[Leonardo Uieda]

Hi,

Could you please tell me which version you installed, what OS, which Python (Anaconda, etc) and which script you ran?

--
Leonardo Uieda

Visiting Research Scholar

Department of Earth Sciences - SOEST

University of Hawai'i at Mānoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

--
You received this message because you are subscribed to the Google Groups "Fatiando a Terra" group.

To unsubscribe from this group and stop receiving emails from it, send an email to fatiando+u...@googlegroups.com.

[accep handyarso]

Hi Leo

I would like to apologize for my late reply.

the code (fatiando) is already works Leo, thank you for your kindly attention and help

at that time, I wiped out my Ubuntu 18.04 and reinstall everythings.

and now everything works well,

my OS is Ubuntu 18.04

python 2.7 (i use anaconda 2)

the scripts I run is gravity harvester

and the result is as below

Leo, may I know your pseudo code for the seeds inversion??

I hava confussions with the determination of the direction "where the seeds should grow?"

do you applied some constrains? or boundary for your model? (suppose, the seeds should not exceed certain boundaries or any other rule may be..)

or it is purely calculation such as from the gradient method (suppose gauss newton etc)

I mean, which one first, you move the next seeds and calculates the best movement from several movements before or

you calculate to determine which the best position or direction for the new seeds and then after that you put the next seeds there?

Thank you for your kindly help Leo

regards

acceph

[Leonardo Uieda]

Hi Acceph,

> Leo, may I know your pseudo code for the seeds inversion??

The source code for the seed inversion is at:

And the algorithm is described in these two papers:

* http://www.leouieda.com/papers/paper-planting-anomalous-densities-2012.html

* http://www.leouieda.com/talks/seg2012.html

> do you applied some constrains? or boundary for your model?

The boundaries are determined by the boundaries of the input mesh.

You can specify if you don't want seeds to grow in a given direction.

This is implemented by removing the neighbors in that direction from the search algorithm.

> or it is purely calculation such as from the gradient method

There is no gradient calculation in this method.

> I mean, which one first, you move the next seeds and calculates the best movement from several movements before or

> you calculate to determine which the best position or direction for the new seeds and then after that you put the next seeds there?

I'm not sure I understand what you mean by that.

Could you please clarify?

--
Leonardo Uieda

Visiting Research Scholar

Department of Earth Sciences - SOEST

University of Hawai'i at Mānoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

[acceph]

Hello Leo

thank you for your reply and those several paper links..

I think I have found the answer inside your code. it is in harvester.harvest isn't it?

I still confused with the compactness which is incorporate inside this algorithm...

looks like more difficult for me rather than in the formal form of inversion

regards

acceph

[Leonardo Uieda]

It's definitely a different kind of algorithm.

I have to admit that the code in harvester is not the cleanest implementation.

I'm planning an updated version for Harmonica (https://github.com/fatiando/harmonica) which will be a lot better.

--
Leonardo Uieda

Visiting Research Scholar

Department of Earth Sciences - SOEST

University of Hawai'i at Mānoa

ORCID: 0000-0001-6123-9515

www.leouieda.com

[acceph]

thank you for the information. nice work leo, keep going...

[Fabiano Della Justina]

Hi Leo, 

How have you been?

I have had the same issue pointed out by Acceph - Something related to vremove.

However, I am using windows 10 and like him python 2.7.

Actually, the more I have tried to resolve the problem the worse it has become. The first time that I installed it, I ran your below code:

https://www.fatiando.org/v0.5/cookbook/gravmag_harvester_grav.html

and I could see the data and the predicted data, but I couldn't get the model solution, then I tried to solve that, and now I can not see any plot, just errors.

I really like your method, it is really intuitive and I am going to present it in a final project to an inverse problem course. Although I was supposed to apply to real data, at this time I would be happy only apply in synthetic ones. Please see below the error that has popped out:

C:\Users\fjustina\Anaconda2\lib\site-packages\fatiando\vis\mpl.py:70: UserWarning: This module will be removed in v0.6. We recommend the use of matplotlib.pyplot module directly. Some of the fatiando specific functions will remain.
  "specific functions will remain.")

ImportErrorTraceback (most recent call last)
<ipython-input-1-981edd511661> in <module>()
      5 from fatiando import gridder, utils
      6 from fatiando.gravmag import polyprism, harvester
----> 7 from fatiando.mesher import PolygonalPrism, PrismMesh, vremove
      8 from fatiando.vis import mpl, myv
      9 

ImportError: cannot import name vremove

Please help me out,

Fabiano

To unsubscribe from this group and stop receiving emails from it, send an email to fati...@googlegroups.com.