Skip to content

Is it possible to fit a variogram on 2D image and generate an equivalent 3D field using gstools? #286

@rkenko

Description

@rkenko

Hi,

I am new in geostats and i wonder if it is possible to fit a variogram on a 2D field then generate a 3D equivalent field.

This is the code i am using:

#importing the image from which we seek the variogram 
image       = plt.imread("field.png")[:,:,0]
# spreading values of imported image between ]-inf;+inf[
image       = np.tan(((10/np.pi)*(image))+(np.pi/2)) 


# estimating empirical variogram
gamma = gs.vario_estimate_structured(image)
bin_center = np.arange(len(gamma))


# fiiting the empirical variogram at small range, 128 px
bins        = np.arange(128)
fit_model   = gs.Gaussian(dim=2)
fit_model.fit_variogram(bin_center[:128], gamma[:128], nugget=False)
# printing all the parameters of the covariance model fitted
print(fit_model)


# generating 3D field from fitted variogram with size 256
x  = y  = z = range(256)
model       = gs.Gaussian(dim=3, var=fit_model.var, len_scale=fit_model.len_scale)
# generating the spatial random field
srf         = gs.SRF(model)
srf((x, y, z), mesh_type = 'structured')
srf         = srf.field


# passing the field through atan in order to spread distribution between [0;1]
srf         = (np.pi/10)*(np.arctan(srf)-(np.pi/2))+1
print(np.amin(srf))
print(np.amax(srf))
plt.imsave("x_middle_slice.png",srf[:,:,128],cmap="gray",vmin=0,vmax=1)
# saving to a VTK file for visualization
srf         = pv.wrap(srf)
srf.save('gs.vtk')

My question is: am i doing right ? the purpose of the tanh is to ensure values stays in the range [0,1] since they represents pixels values in a picture. without the tanh, the field isn't realistic in terms of colors.

I observed the slices of the generated 3D volume dont have the same variogram than the 2D picture of the begining. I expected them to be the same... Is there a way to correct my code please ?

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type
    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions