Skip to content

Update maths.py #229

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 4 commits into from
May 14, 2025
Merged
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
1 change: 1 addition & 0 deletions LoopStructural/utils/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
get_strike_vector,
get_vectors,
strikedip2vector,
plungeazimuth2vector,
azimuthplunge2vector,
normal_vector_to_strike_and_dip,
rotate,
Expand Down
25 changes: 1 addition & 24 deletions LoopStructural/utils/helper.py
Original file line number Diff line number Diff line change
Expand Up @@ -107,30 +107,7 @@ def region(xyz):
return bb, region


# def azimuthplunge2vector(
# plunge: Union[np.ndarray, list], plunge_dir: Union[np.ndarray, list]
# ) -> np.ndarray:
# """Convert plunge and plunge direction to a vector

# Parameters
# ----------
# plunge : Union[np.ndarray, list]
# array or array like of plunge values
# plunge_dir : Union[np.ndarray, list]
# array or array like of plunge direction values

# Returns
# -------
# np.array
# nx3 vector
# """
# plunge = np.deg2rad(plunge)
# plunge_dir = np.deg2rad(plunge_dir)
# vec = np.zeros(3)
# vec[0] = np.sin(plunge_dir) * np.cos(plunge)
# vec[1] = np.cos(plunge_dir) * np.cos(plunge)
# vec[2] = -np.sin(plunge)
# return vec



def create_surface(bounding_box, nstep):
Expand Down
91 changes: 74 additions & 17 deletions LoopStructural/utils/maths.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,20 +28,28 @@ def strikedip2vector(strike: NumericInput, dip: NumericInput) -> np.ndarray:
vec /= np.linalg.norm(vec, axis=1)[:, None]
return vec


def azimuthplunge2vector(
plunge: NumericInput,
azimuth: NumericInput,
degrees: bool = True,
) -> np.ndarray:
raise DeprecationWarning(
"azimuthplunge2vector is deprecated, use plungeazimuth2vector instead"
)

def plungeazimuth2vector(
plunge: NumericInput,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Did you mean to change the order of the variables?

Copy link
Contributor Author

@AxMeNi AxMeNi Apr 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I meant to make them coherent with the name of the function. I have no idea if that would change anything with the actual code though. Do you know if there is a way to access every single call to that function in the whole project?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@AxMeNi, I think I would rather keep it as plunge and plunge_dir because that is the order we refer to them as geologists. But maybe we should change the function name, but we need to add a deprecation to the old naming

plunge_dir: NumericInput,
azimuth: NumericInput,
degrees: bool = True,
) -> np.ndarray:
"""Convert plunge and plunge direction to a vector

Parameters
----------
azimuth : Union[np.ndarray, list]
array or array like of plunge direction values
plunge : Union[np.ndarray, list]
array or array like of plunge values
plunge_dir : Union[np.ndarray, list]
array or array like of plunge direction values

Returns
-------
Expand All @@ -52,16 +60,16 @@ def azimuthplunge2vector(
plunge = np.array([plunge], dtype=float)
else:
plunge = np.array(plunge, dtype=float)
if isinstance(plunge_dir, numbers.Number):
plunge_dir = np.array([plunge_dir], dtype=float)
if isinstance(azimuth, numbers.Number):
azimuth = np.array([azimuth], dtype=float)
else:
plunge_dir = np.array(plunge_dir, dtype=float)
azimuth = np.array(azimuth, dtype=float)
if degrees:
plunge = np.deg2rad(plunge)
plunge_dir = np.deg2rad(plunge_dir)
azimuth = np.deg2rad(azimuth)
vec = np.zeros((len(plunge), 3))
vec[:, 0] = np.sin(plunge_dir) * np.cos(plunge)
vec[:, 1] = np.cos(plunge_dir) * np.cos(plunge)
vec[:, 0] = np.sin(azimuth) * np.cos(plunge)
vec[:, 1] = np.cos(azimuth) * np.cos(plunge)
vec[:, 2] = -np.sin(plunge)
return vec

Expand Down Expand Up @@ -204,19 +212,21 @@ def get_vectors(normal: NumericInput) -> Tuple[np.ndarray, np.ndarray]:


def get_strike_vector(strike: NumericInput, degrees: bool = True) -> np.ndarray:
"""Return the vector aligned with the strike direction
"""Return strike direction vector(s) from strike angle(s).

Parameters
----------
strike : np.ndarray
strike direction
strike : NumericInput
Single strike angle or array-like of strike angles, measured clockwise from North.
degrees : bool, optional
whether to return in degrees or radians, by default True
Whether the input angles are in degrees. If False, angles are assumed to be in radians.
Default is True.

Returns
-------
np.ndarray
vector aligned with strike direction
Array of shape (3, n) where each column is a 3D unit vector (x, y, z) representing
the horizontal strike direction. The z-component is always 0.

"""
if isinstance(strike, numbers.Number):
Expand All @@ -236,6 +246,21 @@ def get_strike_vector(strike: NumericInput, degrees: bool = True) -> np.ndarray:


def get_dip_vector(strike, dip):
"""Return the dip vector based on strike and dip angles.

Parameters
----------
strike : float
Strike angle in degrees, measured clockwise from North.
dip : float
Dip angle in degrees, measured from the horizontal plane.

Returns
-------
np.ndarray
Unit vector (length 3) representing the dip direction in 3D space.

"""
v = np.array(
[
-np.cos(np.deg2rad(-strike)) * np.cos(-np.deg2rad(dip)),
Expand All @@ -247,6 +272,23 @@ def get_dip_vector(strike, dip):


def regular_tetraherdron_for_points(xyz, scale_parameter):
"""Generate regular tetrahedrons centered at given 3D points.

Parameters
----------
xyz : np.ndarray
Array of shape (n, 3) representing the coordinates of n points in 3D space,
which will serve as the centers of the generated tetrahedrons.
scale_parameter : float
Scaling factor controlling the size of the regular tetrahedrons.

Returns
-------
np.ndarray
Array of shape (n, 4, 3) representing n regular tetrahedrons, where each
tetrahedron has 4 vertices in 3D space, positioned relative to the corresponding center point.

"""
regular_tetrahedron = np.array(
[
[np.sqrt(8 / 9), 0, -1 / 3],
Expand All @@ -264,8 +306,23 @@ def regular_tetraherdron_for_points(xyz, scale_parameter):


def gradient_from_tetrahedron(tetrahedron, value):
"""
Calculate the gradient from a tetrahedron
"""Compute the gradient of values within tetrahedral elements

Parameters
----------
tetrahedron : np.ndarray
Array of shape (n, 4, 3) representing the coordinates of tetrahedral elements,
where each tetrahedron is defined by 4 vertices in 3D space.
value : np.ndarray
Array of shape (n, 4) representing the scalar values at the 4 vertices
of each tetrahedron.

Returns
-------
np.ndarray
Array of shape (n, 3) representing the gradient vector of the scalar field
inside each tetrahedral element.

"""
tetrahedron = tetrahedron.reshape(-1, 4, 3)
m = np.array(
Expand Down
16 changes: 8 additions & 8 deletions tests/unit/utils/test_conversions.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
from LoopStructural.utils import strikedip2vector, azimuthplunge2vector
from LoopStructural.utils import strikedip2vector, plungeazimuth2vector
import numpy as np


Expand All @@ -19,18 +19,18 @@ def test_strikedip2vector():


# import numpy as np
# from LoopStructural.utils.maths import azimuthplunge2vector
# from LoopStructural.utils.maths import plungeazimuth2vector


def test_azimuthplunge2vector_single_values():
def test_plungeazimuth2vector_single_values():
plunge = 0
plunge_dir = 90
expected_result = np.array([[1, 0, 0]])
result = azimuthplunge2vector(plunge, plunge_dir)
result = plungeazimuth2vector(plunge, plunge_dir)
assert np.allclose(result, expected_result)


def test_azimuthplunge2vector_array_values():
def test_plungeazimuth2vector_array_values():
plunge = [0, 90, 0]
plunge_dir = [90, 90, 0]
expected_result = np.array(
Expand All @@ -40,12 +40,12 @@ def test_azimuthplunge2vector_array_values():
[0, 1, 0],
]
)
result = azimuthplunge2vector(plunge, plunge_dir)
result = plungeazimuth2vector(plunge, plunge_dir)
assert np.allclose(result, expected_result)


def test_azimuthplunge2vector_empty_arrays():
def test_plungeazimuth2vector_empty_arrays():
plunge = []
plunge_dir = []
result = azimuthplunge2vector(plunge, plunge_dir)
result = plungeazimuth2vector(plunge, plunge_dir)
assert result.shape[0] == 0