Skip to content
Merged
Show file tree
Hide file tree
Changes from 1 commit
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
52 changes: 52 additions & 0 deletions geoapps_utils/modelling/plates.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,9 +8,12 @@
# '
# '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''

import warnings

import numpy as np
from geoh5py import Workspace
from geoh5py.objects import Octree, Surface
from geoh5py.objects.maxwell_plate import MaxwellPlate, PlateGeometry, PlatePosition
from geoh5py.shared.utils import fetch_active_workspace
from pydantic import BaseModel, ConfigDict, Field

Expand Down Expand Up @@ -47,6 +50,40 @@ class PlateModel(BaseModel):
direction: float = Field(default=0.0, alias="dip_direction")
dip: float = 0.0

@classmethod
def from_maxwell_plate_geometry(cls, geometry: PlateGeometry):

if geometry.rotation != 0.0:
warnings.warn(
"Plunging plate models are not yet implemented. "
"Ignoring the maxwell plate geometry rotation."
Comment thread
benk-mira marked this conversation as resolved.
Outdated
)

return PlateModel(
strike_length=geometry.width,
dip_length=geometry.length,
width=geometry.thickness,
easting=geometry.position.x,
northing=geometry.position.y,
elevation=geometry.position.z,
direction=geometry.dip_direction,
dip=geometry.dip,
)
Comment thread
benk-mira marked this conversation as resolved.

def to_maxwell_plate_geometry(self) -> PlateGeometry:
return PlateGeometry(
position=PlatePosition(
x=self.easting,
y=self.northing,
z=self.elevation,
),
dip=self.dip,
dip_direction=self.direction,
length=self.dip_length,
width=self.strike_length,
thickness=self.width,
)

@property
def origin(self) -> tuple[float, float, float]:
return (self.easting, self.northing, self.elevation)
Expand All @@ -62,6 +99,21 @@ class Plate:
def __init__(self, params: PlateModel):
self.params = params

@classmethod
def from_maxwell_plate(cls, plate: MaxwellPlate):
if plate.geometry is None:
raise ValueError("Maxwell plate must have its geometry set.")
return Plate(PlateModel.from_maxwell_plate_geometry(plate.geometry))
Comment thread
benk-mira marked this conversation as resolved.
Outdated

def to_maxwell_plate(
self, workspace: Workspace, name: str | None = None
) -> MaxwellPlate:
with fetch_active_workspace(workspace) as ws:
plate = MaxwellPlate.create(
ws, name=name, geometry=self.params.to_maxwell_plate_geometry()
)
return plate

def mask(self, mesh: Octree) -> np.ndarray:
rotations = [
z_rotation_matrix(np.deg2rad(self.params.direction)),
Expand Down
35 changes: 34 additions & 1 deletion tests/plates_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,10 +11,11 @@
from __future__ import annotations

import numpy as np
import pytest
from geoh5py import Workspace
from geoh5py.objects import BlockModel

from geoapps_utils.modelling.plates import PlateModel, inside_plate, make_plate
from geoapps_utils.modelling.plates import Plate, PlateModel, inside_plate, make_plate


def test_inside_plate(tmp_path):
Expand Down Expand Up @@ -200,3 +201,35 @@ def test_plate_alias():
)
assert plate.direction == 90
assert "dip_direction" in plate.model_dump(by_alias=True)


def test_maxwell_plate_integration(tmp_path):

plate = Plate(
PlateModel(
strike_length=100,
dip_length=300,
width=20,
easting=100.0,
northing=0.0,
elevation=0.0,
dip_direction=90,
dip=45,
)
)
with Workspace(tmp_path / "test.geoh5") as workspace:
maxwell_plate = plate.to_maxwell_plate(workspace)
assert maxwell_plate.geometry is not None
maxwell_plate.geometry.rotation = 10.0

with pytest.warns(UserWarning, match="Plunging plate"):
plate = Plate.from_maxwell_plate(maxwell_plate)

assert plate.params.strike_length == 100
assert plate.params.dip_length == 300
assert plate.params.width == 20
assert plate.params.easting == 100
assert plate.params.northing == 0
assert plate.params.elevation == 0
assert plate.params.direction == 90
assert plate.params.dip == 45
Comment thread
benk-mira marked this conversation as resolved.
Outdated
Loading