Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
24 changes: 20 additions & 4 deletions simpeg_drivers/options.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@

import numpy as np
from geoapps_utils.base import Options
from geoapps_utils.utils.importing import GeoAppsError
from geoh5py.data import (
BooleanData,
DataAssociationEnum,
Expand Down Expand Up @@ -458,7 +459,7 @@ def property_group_data(self, property_group: PropertyGroup):
k for k in self.data_object.property_groups if k.uid == property_group.uid
)
data = {
freq: self.geoh5.get_entity(p)[0].values
freq: self.data_object.get_entity(p)[0].values
for freq, p in zip(frequencies, group.properties, strict=False)
}
return data
Expand Down Expand Up @@ -615,20 +616,35 @@ def data(self) -> dict[str, dict[float, np.ndarray | None]]:

@property
def uncertainties(self) -> dict[str, dict[float, np.ndarray | None]]:
"""Return dictionary of unceratinty components and associated values."""
"""Return dictionary of uncertainty components and associated values."""
out = {}
flags = []
for k in self.active_components:
out[k] = self.component_uncertainty(k)

for data in out[k].values():
if np.any(np.isnan(data)) or np.any(data < 0):
flags.append(f"{k} component")
break

if flags:
summary = (
"Issues encountered with uncertainties having NDV or negative values:\n\n - "
+ "\n - ".join(flags)
)
Comment thread
domfournier marked this conversation as resolved.
summary += "\n\nPlease review the input values."
raise GeoAppsError(summary)

return out

def component_data(self, component: str) -> np.ndarray | None:
def component_data(self, component: str) -> dict:
"""Return data values associated with the component."""
data = getattr(self, "_".join([component, "channel"]), None)
if isinstance(data, NumericData):
data = data.values
return {None: data}

def component_uncertainty(self, component: str) -> np.ndarray | None:
def component_uncertainty(self, component: str) -> dict:
"""
Return uncertainty values associated with the component.

Expand Down
26 changes: 26 additions & 0 deletions tests/run_tests/driver_grav_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -206,6 +206,32 @@ def test_array_too_large_run(
driver.run()


def test_bad_uncertainties(
tmp_path: Path,
):
workpath = tmp_path / f"{__name__}.ui.geoh5"

with Workspace(workpath) as geoh5:
components = SyntheticsComponents(geoh5)
gz = components.survey.add_data(
{"data": {"values": np.random.randn(components.survey.n_vertices)}}
)

# Run the inverse
params = GravityInversionOptions.build(
geoh5=geoh5,
mesh=components.mesh,
topography_object=components.topography,
data_object=gz.parent,
gz_channel=gz,
gz_uncertainty=-1e-4,
starting_model=1e-4,
)

with raises(GeoAppsError, match="Issues encountered with uncertainties"):
_ = params.uncertainties


if __name__ == "__main__":
# Full run
test_gravity_fwr_run(
Expand Down
45 changes: 44 additions & 1 deletion tests/run_tests/driver_mt_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,9 @@
from pathlib import Path

import numpy as np
from geoapps_utils.utils.importing import GeoAppsError
from geoh5py.workspace import Workspace
from pytest import raises

from simpeg_drivers.natural_sources.magnetotellurics import (
MTForwardDriver,
Expand Down Expand Up @@ -73,7 +75,7 @@ def setup_data(workspace, survey):
}
}
)
uncertainties[f"{cname} uncertainties"].append(uncert.copy(parent=survey))
uncertainties[f"{cname} uncertainties"].append(uncert)

data_groups = survey.add_components_data(data)
uncert_groups = survey.add_components_data(uncertainties)
Expand Down Expand Up @@ -141,6 +143,47 @@ def test_magnetotellurics_fwr_run(
assert geoh5.get_entity("Iteration_0_zyy_real_[0]")[0] is not None


def test_bad_uncertainties(
tmp_path: Path,
):
workpath = (
tmp_path.parent / "test_magnetotellurics_fwr_run0" / "inversion_test.ui.geoh5"
)

Comment thread
MatthieuCMira marked this conversation as resolved.
with Workspace(workpath) as geoh5:
with Workspace.create(tmp_path / f"{__name__}.geoh5") as ws:
components = SyntheticsComponents(geoh5)
survey = components.survey.copy(
parent=ws, copy_children=False, name="bad_uncertainties"
)
mesh = components.mesh
topography = components.topography
data_kwargs = setup_data(geoh5, survey)

for elem in ["uncertainty_zxx_imag_[0]", "uncertainty_zyx_real_[0]"]:
data = survey.get_entity(elem)[0]
vals = data.values
vals[0] = np.nan
data.values = vals

# Run the inverse
params = MTInversionOptions.build(
geoh5=geoh5,
mesh=mesh,
topography_object=topography,
data_object=survey,
starting_model=100.0,
background_conductivity=100.0,
**data_kwargs,
)

with raises(GeoAppsError) as error:
_ = params.uncertainties

assert "zxx_imag" in str(error.value)
assert "zyx_real" in str(error.value)


def test_magnetotellurics_run(tmp_path: Path, max_iterations=1, pytest=True):
# pass
workpath = tmp_path / "inversion_test.ui.geoh5"
Expand Down
Loading