Space domain upward and downward continuation

[leouieda]

Description of the desired feature

Just came across the paper Potential field continuation in spatial domain: A new kernel function and its numerical scheme by Zhikui Guo (@zguoch) and Chunhui Tao. The method looks really promising and would be a cool thing to have on Harmonica.

The paper has an accompanying C++ implementation that is on GitHub conti2d but there is no open-source license so we can't use it directly (which probably wouldn't be a good idea anyway).

Are you willing to help implement and maintain this feature? No

@zguoch would you have any interest in contributing a Python implementation of your algorithm to Harmonica? We can help you get started and work through the development process (making pull requests, testing, documentation, etc).

The algorithm can be made fast in Python using numba.

[drandykass]

I may be able to help out with this if we get a small project funded. Feel free to add me and keep me in the loop. I would like to become involved again.

[leouieda]

Awesome, @drandykass! You're already a member of the organization so I assigned you this issue. Please keep us updated and we'd love to see you at our Developer calls on Fridays 13:00 UTC 🙂

[drandykass]

Sure thing. I’d love to join and get back into it.

Sent from mobile device. Please excuse typos.

Den 12. okt. 2021 kl. 09.44 skrev Leonardo Uieda @.***>:

 Awesome, @drandykass! You're already a member of the organization so I assigned you this issue. Please keep us updated and we'd love to see you at our Developer calls on Fridays 13:00 UTC 🙂

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or unsubscribe. Triage notifications on the go with GitHub Mobile for iOS or Android.

[leouieda]

Awesome! More info on that here: https://www.fatiando.org/community/ Hope to see you there

[RichardScottOZ]

Speaking of downward continuation - something like?

import numpy as np
import xarray as xr

# Define the gravity anomaly function and contour
def g(z):
    return 1 / (1 + z ** 2)

# Define the Cauchy integral function for downward continuation
def cauchy_integral(z, z_depth):
    def integrand(t, z, z_depth):
        w = np.exp(1j * t)
        return g(w) * (w / (z - w) ** 2) * np.exp(-z_depth * np.abs(w))

    integrand_vec = np.vectorize(integrand)
    result = -integrand_vec(np.linspace(0, 2 * np.pi, 100), z, z_depth).sum()
    return result / (2 * np.pi) * np.exp(-z_depth * np.abs(z))

# Create a sample gravity grid using xarray
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

coords = {'x': x, 'y': y}
dims = ('y', 'x')

C = xr.DataArray(X + 1j * Y, coords=coords, dims=dims)

# Define the depth at which to continue the gravity field
z_depth = 1000  # meters

G = xr.apply_ufunc(cauchy_integral, C, z_depth,
                   input_core_dims=[['x', 'y'], []],
                   output_core_dims=[['x', 'y']],
                   output_dtypes=[np.complex128])

# Plot the original and downward continued gravity fields
G.real.plot()
``` @

[RichardScottOZ]

https://www.miningnews.net/partners/partner-content/1441144/sensore-intrepid-release-new-geophysics-products