Loop3D/LoopStructural

adding unconformity and DTM - clipping the stratigraphy

Closed this issue · 1 comments

gpirot commented

I am not sure I have been able to add DTM data using the create_and_add_dtm function as it did not clip the stratigraphy nor in the surface rendering, nor in the volume rendering.

With LoopStructural version 1.5.3, the create_and_add_dtm function returns in my case:

DTM:
dtmdata dtm 
	0 regions
	0 faults.
	Fault enabled True

With LoopStructural version 1.4.9, , the create_and_add_dtm function returns in my case:

DTM: 
dtmdata dtm 
	1 regions
	 	<method-wrapper '__str__' of RegionEverywhere object at 0x000001692089FE50>
	0 faults.
	Fault enabled True

To get there, I adapted the following documentation example (see attache self-contained notebook ls_test_dtm.zip ) from: https://loop3d.org/LoopStructural/_auto_examples/1_basic/plot_3_model_visualisation.html#sphx-glr-auto-examples-1-basic-plot-3-model-visualisation-py

I noted that for both versions 1.5.3 or 1.4.9, the stratigraphic surfaces were also not clipped by the unconformity.

So, what shall I do to add a dtm properly to a model, and to get the stratigraphy clipped both in surface rendering and volume rendering?

gpirot commented

Integrating the dtm as a property to the model, such as the example below worked.

clip_on_dtm=True
if(clip_on_dtm):
    dtm = gdal.Open(os.path.join(model_name, "StupidGDALLocalFile.tif"))
    dtm_val = dtm.GetRasterBand(1).ReadAsArray()
    # Convert bounds to gdal raster bounds
    x=np.linspace(minx,maxx,dtm_val.shape[0])
    y=np.linspace(miny,maxy,dtm_val.shape[0])
    dtm_interpolator = RegularGridInterpolator((x,y),dtm_val)
    model.dtm = lambda xyz : dtm_interpolator(xyz[:,:2])

The updated notebook related to the initial case is here.
ls_test_dtm_new_way.zip

Stratigraphic surfaces are clipped by the unconformity in version 1.5.3 (but not 1.4.9).