Small Bug with non-constant density contrast

[LL-Geo]

Describe the bug
I had some fun playing the package. Well done and nice design @mdtanker 👍!

Just find a small bug that when the input density is non-constant, it will still force the density contrast to be constant. And of course the topography is very off in these wrong density region.

Full code that generate the error

# the density contrast is between rock (~2670 kg/m3) and air (~1 kg/m3)
density_contrast1 = 2670 - 1
density_contrast2 = 1322

# prisms are created between the mean topography value and the height of the topography
zref = true_topography.values.mean()

# prisms above zref have positive density contrast and prisms below zref have negative
# density contrast
density_grid = xr.where(true_topography > zref, density_contrast1, -density_contrast2)

That should be easy to fix In

def update_prisms_ds(
    prisms_ds: xr.Dataset,
    correction_grid: xr.DataArray,
    zref: float,
) -> xr.Dataset:

Change

    density_contrast = ds.density.values.max()

To

   density_contrast = np.fabs(ds.density)

This should solve the issue.😁

[mdtanker]

Hey @LL-Geo thanks for checking out the package! Yes that's a good catch, I've only been using constant densities so far, but yes we should allow variable densities. The only other place it may need some work to implement would be the density cross validation. I think we would just need to clarify in the documentation that the density CV only works for constant density values. Maybe there could be some scheme for a CV for finding the optimal density values for different regions, i.e. do a CV to find optimal density for 1 region based on a mask, then do another CV for the remaining areas.

I'll try and add some documentation for now.

[LL-Geo]

Hey @LL-Geo thanks for checking out the package! Yes that's a good catch, I've only been using constant densities so far, but yes we should allow variable densities. The only other place it may need some work to implement would be the density cross validation. I think we would just need to clarify in the documentation that the density CV only works for constant density values. Maybe there could be some scheme for a CV for finding the optimal density values for different regions, i.e. do a CV to find optimal density for 1 region based on a mask, then do another CV for the remaining areas.

I'll try and add some documentation for now.

Yeh~ We can try something like this.. 😁
https://academic.oup.com/gji/article/221/3/1896/5808816?login=true
Basically, divide the region into several blocks and get the density contrast. That might be a bit tricky for bathymetry inversion as we don't have much control points, haha. (but we can do a clustering to separate domains based on gravity and mag).

[mdtanker]

Yeah that sounds great, happy to implement that at some point! Thanks for sharing that paper, don't think I had seen that one yet.

[LL-Geo]

Yeah that sounds great, happy to implement that at some point! Thanks for sharing that paper, don't think I had seen that one yet.

haha.. looking forward to it 😉