search

gempy-project / gempy

GemPy is an open-source, Python-based 3-D structural geological modeling software, which allows the implicit (i.e. automatic) creation of complex geological models from interface and orientation data. It also offers support for stochastic modeling to address parameter and model uncertainties.
https://gempy.org
European Union Public License 1.2
941 stars 232 forks source link

Legacy vs current GEMPY version #897

Closed Al3xdevito closed 1 month ago

Al3xdevito commented 2 months ago

Hi,

I have been trying to add a velocity model to the geological model, as the Transform2021 tutorial (https://github.com/devitocodes/transform2021/blob/main/creating_synthetics.ipynb). Looks like the current version of GEMPY does not support this functions, for example: AttributeError Traceback (most recent call last) Cell In[26], line 1 ----> 1 geo_model.add_surface_values([1.5,1.5, 1.75, 2.5,4.5], ['vp']) AttributeError: 'GeoModel' object has no attribute 'add_surface_values'

My question is: Should I try the Gempy legacy first or try to debug the current gempy version.

javoha commented 2 months ago

Hi, good question - I am not 100% sure but I was also not able to find the exact equivalent for this in gempy v3. Let me check again and come back to you. For now: You can access the results and set a numpy array with the shape of the geo_model.solutions.raw_arrays.lith_block and assign values based on this lith block. You can basically use the same logic used in the notebook you provided afterwards. I tried this quickly just using the lith block for another model:

shape=[50, 20, 20] # corresponding to resolution of gempy model

# Reshaping our data to the shape required by Devito
reshaped = np.reshape(geo_model.solutions.raw_arrays.lith_block, shape, order='C')
reshaped.shape

import matplotlib.pyplot as plt

# Take the center slice in the y direction in this case
# Remember that in Devito, indexing convention is [x, y, z] (need to flip for correct imshow display)
plt.imshow(reshaped[:,5,:].T, cmap='viridis', origin='lower')
plt.colorbar()
plt.show()

This gave me something like this (sorry for low res), you would still need to create the array with your physical properties based on the lith block array.

issue

Does this work for your case? Of course we would like you to work with gempy v3.

Cheers Jan

javoha commented 1 month ago

Hi @Al3xdevito, I will close this issue for now, please feel free to reopen if you need further assistance.