Skip to content

Commit

Permalink
Change default value for depth in EquivalentSourcesGB (#515)
Browse files Browse the repository at this point in the history
Change the default value for the `depth` argument in
`EquivalentSourcesGB` to `"default"`. Set the depth argument in the same
way as `EquivalentSources` does, relying on the inherited methods.
Extend docstrings to reflect these changes.
  • Loading branch information
santisoler authored Jun 25, 2024
1 parent 5375826 commit 5091fcb
Show file tree
Hide file tree
Showing 2 changed files with 41 additions and 7 deletions.
20 changes: 15 additions & 5 deletions harmonica/_equivalent_sources/gradient_boosted.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,8 @@
"""
Gradient-boosted equivalent sources in Cartesian coordinates
"""
from __future__ import annotations

import warnings

import numpy as np
Expand Down Expand Up @@ -47,13 +49,16 @@ class EquivalentSourcesGB(EquivalentSources):
If None, will place one point source below each observation point at
a fixed relative depth below the observation point [Cooper2000]_.
Defaults to None.
depth : float
depth : float or "default"
Parameter used to control the depth at which the point sources will be
located.
Each source is located beneath each data point (or block-averaged
location) at a depth equal to its elevation minus the ``depth`` value.
If a value is provided, each source is located beneath each data point
(or block-averaged location) at a depth equal to its elevation minus
the ``depth`` value.
If set to ``"default"``, the depth of the sources will be estimated as
4.5 times the mean distance between first neighboring sources.
This parameter is ignored if *points* is specified.
Defaults to 500.
Defaults to ``"default"``.
block_size: float, tuple = (s_north, s_east) or None
Size of the blocks used on block-averaged equivalent sources.
If a single value is passed, the blocks will have a square shape.
Expand Down Expand Up @@ -87,6 +92,10 @@ class EquivalentSourcesGB(EquivalentSources):
The boundaries (``[W, E, S, N]``) of the data used to fit the
interpolator. Used as the default region for the
:meth:`~harmonica.EquivalentSources.grid` method.
depth_ : float or None
Estimated depth of the sources calculated as 4.5 times the mean
distance between first neighboring sources. This attribute is set to
None if ``points`` is passed.
window_size_ : float or None
Size of the overlapping windows used in gradient-boosting equivalent
point sources. It will be set to None if ``window_size = "default"``
Expand All @@ -105,7 +114,7 @@ def __init__(
self,
damping=None,
points=None,
depth=500,
depth: float | str = "default",
block_size=None,
window_size="default",
parallel=True,
Expand Down Expand Up @@ -221,6 +230,7 @@ def fit(self, coordinates, data, weights=None):
p.astype(self.dtype) for p in self._build_points(coordinates)
)
else:
self.depth_ = None # set depth_ to None so we don't leave it unset
self.points_ = tuple(
p.astype(self.dtype) for p in vdb.n_1d_arrays(self.points, 3)
)
Expand Down
28 changes: 26 additions & 2 deletions harmonica/tests/test_gradient_boosted_eqs.py
Original file line number Diff line number Diff line change
Expand Up @@ -115,7 +115,7 @@ def test_custom_points(region, coordinates_small, data_small):
i.ravel()
for i in vd.grid_coordinates(region=region, shape=(3, 3), extra_coords=-550)
)
eqs = EquivalentSourcesGB(points=points_custom, window_size=500)
eqs = EquivalentSourcesGB(points=points_custom, window_size=500, depth=500)
eqs.fit(coordinates_small, data_small)
# Check that the proper source locations were set
npt.assert_allclose(points_custom, eqs.points_, rtol=1e-5)
Expand Down Expand Up @@ -171,7 +171,7 @@ def test_gradient_boosted_eqs_single_window(region, points, masses, coordinates,
"""
Test GB eq-sources with a single window that covers the whole region
"""
eqs = EquivalentSourcesGB(window_size=region[1] - region[0])
eqs = EquivalentSourcesGB(depth=500, window_size=region[1] - region[0])
eqs.fit(coordinates, data)
npt.assert_allclose(data, eqs.predict(coordinates), rtol=1e-5)
# Gridding onto a denser grid should be reasonably accurate when compared
Expand Down Expand Up @@ -410,3 +410,27 @@ def test_window_size():
def test_invalid_window_size():
with pytest.raises(ValueError, match="Found invalid 'window_size' value equal to"):
EquivalentSourcesGB(window_size="Chuckie took my soul!")


def test_default_depth(coordinates, data):
"""
Test if the depth of sources is correctly set by the default strategy
"""
# Get distance to first neighbour in the grid
easting, northing = coordinates[:2]
d_easting = easting[1, 1] - easting[0, 0]
d_northing = northing[1, 1] - northing[0, 0]
first_neighbour_distance = min(d_easting, d_northing)
# Fit the equivalent sources with default `depth`
eqs = EquivalentSourcesGB().fit(coordinates, data)
npt.assert_allclose(eqs.depth_, first_neighbour_distance * 4.5)


def test_invalid_depth():
"""
Test if error is raised after passing invalid value for depth.
"""
invalid_depth = "this is not a valid one"
msg = f"Found invalid 'depth' value equal to '{invalid_depth}'"
with pytest.raises(ValueError, match=msg):
EquivalentSourcesGB(depth=invalid_depth)

0 comments on commit 5091fcb

Please sign in to comment.