Glassvision/IfcOpenShell
1
0
Fork 0
mirror of https://github.com/IfcOpenShell/IfcOpenShell.git synced 2026-08-10 09:48:32 +00:00
Code Issues Packages Projects Releases Wiki Activity

Typing and update geolocation conversion utils to follow IFC4X3 rules about unit scaling

Browse Source
This commit is contained in:
Dion Moult
2024-06-22 15:52:00 +10:00
parent d0373da0ad
commit a6c3ae8446
2 changed files with 215 additions and 137 deletions
Download Patch File Download Diff File Expand all files Collapse all files
src/ifcopenshell-python/ifcopenshell/util/geolocation.py
+123 -123
View File
@@ -18,6 +18,7 @@
import math
import numpy as np
import numpy.typing as npt
import ifcopenshell
import ifcopenshell.util.unit
import ifcopenshell.util.element
@@ -94,7 +95,7 @@ def xyz2enh(
you are applying your own temporary false origin (such as when federating
models for digital twins of large cities).
For most scenarios you should use ``auto_xyz2enh`` instead.
For most scenarios you should use :func:`auto_xyz2enh` instead.
:param x: The X local engineering coordinate.
:param y: The Y local engineering coordinate.
@@ -131,7 +132,7 @@ def auto_xyz2enh(
The necessary georeferencing map conversion is automatically detected from
the IFC map conversion parameters present in the IFC model. If no map
conversion is present, then the Z coordinate is returned unchanged.
conversion is present, then the coordinates are returned unchanged.
For IFC2X3, the map conversion is detected from the IfcProject's
ePSet_MapConversion. See the "User Guide for Geo-referencing in IFC":
@@ -144,7 +145,6 @@ def auto_xyz2enh(
:param should_return_in_map_units: If true, the result is given in map units.
If false, the result will be converted back into project units.
:return: The global map coordinate eastings, northings, and height.
:rtype: tuple[float]
"""
parameters = get_helmert_transformation_parameters(ifc_file)
if not parameters:
@@ -170,6 +170,7 @@ def auto_enh2xyz(ifc_file, easting, northing, height, is_specified_in_map_units:
:param easting: The global easting map coordinate provided in map units.
:param northing: The global northing map coordinate provided in map units.
:param height: The global height map coordinate provided in map units.
:param is_specified_in_map_units: True if the input eastings, northing, and height are in map units.
:return: The local engineering XYZ coordinates in project length units.
"""
parameters = get_helmert_transformation_parameters(ifc_file)
@@ -270,7 +271,7 @@ def z2e(z: float, orthogonal_height: float = 0.0, scale: float = 1.0, factor_z:
This function is for advanced users as it allows you to specify your own
orthogonal height offset and transformation parameters.
For most scenarios you should use ``auto_z2e`` instead.
For most scenarios you should use :func:`auto_z2e` instead.
:param z: The Z local engineering coordinate provided in project length units.
:param orthogonal_height: The orthogonal height offset to apply.
@@ -310,7 +311,7 @@ def enh2xyz(
you are applying your own temporary false origin (such as when federating
models for digital twins of large cities).
For most scenarios you should use ``auto_enh2xyz`` instead.
For most scenarios you should use :func:`auto_enh2xyz` instead.
:param e: The global easting map coordinate.
:param n: The global northing map coordinate.
@@ -342,7 +343,18 @@ def enh2xyz(
return (x, y, z)
def local2global(matrix, eastings, northings, orthogonal_height, x_axis_abscissa, x_axis_ordinate, scale=None):
def local2global(
matrix: npt.NDArray[np.float64],
eastings: float = 0.0,
northings: float = 0.0,
orthogonal_height: float = 0.0,
x_axis_abscissa: float = 1.0,
x_axis_ordinate: float = 0.0,
scale: float = 1.0,
factor_x: float = 1.0,
factor_y: float = 1.0,
factor_z: float = 1.0,
) -> npt.NDArray[np.float64]:
"""Manually convert a 4x4 matrix from local to global coordinates
This function is for advanced users as it allows you to specify your own
@@ -352,63 +364,20 @@ def local2global(matrix, eastings, northings, orthogonal_height, x_axis_abscissa
you are applying your own temporary false origin (such as when federating
models for digital twins of large cities).
No unit conversion is performed.
For most scenarios you should use :func:`auto_local2global` instead.
:param matrix: A 4x4 numpy matrix representing local coordinates.
:type matrix: np.array
:param eastings: The eastings offset to apply.
:type eastings: float
:param northings: The northings offset to apply.
:type northings: float
:param orthogonal_height: The orthogonal height offset to apply.
:type orthogonal_height: float
:param x_axis_abscissa: The X axis abscissa (i.e. first coordinate) of the
2D vector that points to the local X axis when in map coordinates.
:type x_axis_abscissa: float
:param x_axis_ordinate: The X axis ordinate (i.e. second coordinate) of the
2D vector that points to the local X axis when in map coordinates.
:type x_axis_ordinate: float
:param scale: The combined scale factor to convert from local coordinates
to map coordinates.
:type scale: float
:return: A numpy 4x4 array matrix representing global coordinates.
:rtype: np.array
"""
if scale is None:
scale = 1.0
x = np.array([x_axis_abscissa, x_axis_ordinate, 0])
x /= np.linalg.norm(x)
y = np.cross(np.array([0, 0, 1]), x)
intermediate = (
np.array(
[
[x[0], y[0], 0, 0],
[x[1], y[1], 0, 0],
[x[2], y[2], 1, 0],
[0, 0, 0, 1],
]
)
@ matrix
)
intermediate[0, 3] = (intermediate[0, 3] * scale) + eastings
intermediate[1, 3] = (intermediate[1, 3] * scale) + northings
intermediate[2, 3] = (intermediate[2, 3] * scale) + orthogonal_height
return intermediate
def local2global_ifc4x3(
matrix,
eastings,
northings,
orthogonal_height,
x_axis_abscissa,
x_axis_ordinate,
scale=1.0,
factor_x=1.0,
factor_y=1.0,
factor_z=1.0,
):
# Matrix is a 4x4 matrix typically describing the object placement of an element.
theta = math.atan2(x_axis_ordinate, x_axis_abscissa)
scale_and_factor_matrix = np.array(
[
@@ -436,7 +405,45 @@ def local2global_ifc4x3(
return result
def global2local(matrix, eastings, northings, orthogonal_height, x_axis_abscissa, x_axis_ordinate, scale=None):
def auto_local2global(
ifc_file: ifcopenshell.file, matrix: npt.NDArray[np.float64], should_return_in_map_units: bool = True
) -> npt.NDArray[np.float64]:
"""Convert a local matrix to a global map matrix
The necessary georeferencing map conversion is automatically detected from
the IFC map conversion parameters present in the IFC model. If no map
conversion is present, then the matrix is returned unchanged.
:param ifc_file: The IFC file
:param matrix: A 4x4 numpy matrix representing local coordinates.
:param should_return_in_map_units: If true, the result is given in map units.
If false, the result will be converted back into project units.
:return: A numpy 4x4 array matrix representing global coordinates.
"""
parameters = get_helmert_transformation_parameters(ifc_file)
if not parameters:
return matrix.copy()
result = local2global(matrix, *parameters)
if should_return_in_map_units:
return result
result[0][3] /= parameters.scale
result[1][3] /= parameters.scale
result[2][3] /= parameters.scale
return result
def global2local(
matrix: npt.NDArray[np.float64],
eastings: float = 0.0,
northings: float = 0.0,
orthogonal_height: float = 0.0,
x_axis_abscissa: float = 1.0,
x_axis_ordinate: float = 0.0,
scale: float = 1.0,
factor_x: float = 1.0,
factor_y: float = 1.0,
factor_z: float = 1.0,
) -> npt.NDArray[np.float64]:
"""Manually convert a 4x4 matrix from global to local coordinates
This function is for advanced users as it allows you to specify your own
@@ -446,76 +453,91 @@ def global2local(matrix, eastings, northings, orthogonal_height, x_axis_abscissa
you are applying your own temporary false origin (such as when federating
models for digital twins of large cities).
No unit conversion is performed.
:param matrix: A 4x4 numpy matrix representing global coordinates.
:type matrix: np.array
:param eastings: The eastings offset to apply.
:type eastings: float
:param northings: The northings offset to apply.
:type northings: float
:param orthogonal_height: The orthogonal height offset to apply.
:type orthogonal_height: float
:param x_axis_abscissa: The X axis abscissa (i.e. first coordinate) of the
2D vector that points to the local X axis when in map coordinates.
:type x_axis_abscissa: float
:param x_axis_ordinate: The X axis ordinate (i.e. second coordinate) of the
2D vector that points to the local X axis when in map coordinates.
:type x_axis_ordinate: float
:param scale: The combined scale factor to convert from local coordinates
to map coordinates.
:type scale: float
:return: A numpy 4x4 array matrix representing local coordinates.
:rtype: np.array
"""
if scale is None:
scale = 1.0
x = np.array([x_axis_abscissa, x_axis_ordinate, 0])
x /= np.linalg.norm(x)
y = np.cross(np.array([0, 0, 1]), x)
result = matrix.copy()
result[0, 3] = (result[0, 3] - eastings) / scale
result[1, 3] = (result[1, 3] - northings) / scale
result[2, 3] = (result[2, 3] - orthogonal_height) / scale
return (
np.linalg.inv(
np.array(
[
[x[0], y[0], 0, 0],
[x[1], y[1], 0, 0],
[x[2], y[2], 1, 0],
[0, 0, 0, 1],
]
)
)
@ result
theta = math.atan2(x_axis_ordinate, x_axis_abscissa)
scale_and_factor_matrix = np.array(
[
[scale * factor_x, 0, 0, 0],
[0, scale * factor_y, 0, 0],
[0, 0, scale * factor_z, 0],
[0, 0, 0, 1],
]
)
rotation_matrix = np.array(
[
[math.cos(theta), -math.sin(theta), 0, 0],
[math.sin(theta), math.cos(theta), 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1],
]
)
result = matrix.copy()
result[0][3] -= eastings
result[1][3] -= northings
result[2][3] -= orthogonal_height
result = np.linalg.inv(scale_and_factor_matrix) @ np.linalg.inv(rotation_matrix) @ result
result[:, 0][0:3] /= np.linalg.norm(result[:, 0][0:3])
result[:, 1][0:3] /= np.linalg.norm(result[:, 1][0:3])
result[:, 2][0:3] /= np.linalg.norm(result[:, 2][0:3])
return result
def xaxis2angle(x, y):
def auto_global2local(
ifc_file: ifcopenshell.file, matrix: npt.NDArray[np.float64], is_specified_in_map_units: bool = True
) -> npt.NDArray[np.float64]:
"""Convert a global map matrix to a local matrix
The necessary georeferencing map conversion is automatically detected from
the IFC map conversion parameters present in the IFC model. If no map
conversion is present, then the matrix is returned unchanged.
:param ifc_file: The IFC file
:param matrix: A 4x4 numpy matrix representing local coordinates.
:param should_return_in_map_units: If true, the result is given in map units.
If false, the result will be converted back into project units.
:param is_specified_in_map_units: True if the input matrix is in map units.
:return: A numpy 4x4 array matrix representing global coordinates.
"""
parameters = get_helmert_transformation_parameters(ifc_file)
if not parameters:
return matrix.copy()
if not is_specified_in_map_units:
matrix = matrix.copy()
matrix[0][3] *= parameters.scale
matrix[1][3] *= parameters.scale
matrix[2][3] *= parameters.scale
return global2local(matrix, *parameters)
def xaxis2angle(x: float, y: float) -> float:
"""Converts X axis abscissa and ordinates to an angle in decimal degrees
:param x: The X axis abscissa
:type x: float
:param y: The X axis ordinate
:type y: float
:return: The equivalent angle in decimal degrees from the X axis
:rtype: float
"""
return math.degrees(math.atan2(y, x)) * -1
def yaxis2angle(x, y):
def yaxis2angle(x: float, y: float) -> float:
"""Converts Y axis abscissa and ordinates to an angle in decimal degrees
The Y axis abscissa and ordinate is how IFC stores true north.
:param x: The Y axis abscissa
:type x: float
:param y: The Y axis ordinate
:type y: float
:return: The equivalent angle in decimal degrees from the Y axis
:rtype: float
"""
angle = math.degrees(math.atan2(y, x)) - 90
if angle < -180:
@@ -525,7 +547,7 @@ def yaxis2angle(x, y):
return angle
def get_grid_north(ifc_file):
def get_grid_north(ifc_file: ifcopenshell.file) -> float:
"""Get an angle pointing to map grid north
Anticlockwise is positive.
@@ -539,31 +561,15 @@ def get_grid_north(ifc_file):
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
:param ifc_file: The IFC file
:type ifc_file: ifcopenshell.file
:return: An angle to grid north in decimal degrees
:rtype: float
"""
conversion = None
try:
conversion = ifc_file.by_type("IfcMapConversion")[0]
except:
pass
if conversion:
if not conversion.XAxisAbscissa or not conversion.XAxisOrdinate:
return 0
xaa = conversion.XAxisAbscissa
xao = conversion.XAxisOrdinate
else:
project = ifc_file.by_type("IfcProject")[0]
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion")
if not conversion:
return 0
xaa = conversion.get("XAxisAbscissa", None) or 0
xao = conversion.get("XAxisOrdinate", None) or 0
return xaxis2angle(xaa, xao)
parameters = get_helmert_transformation_parameters(ifc_file)
if not parameters:
return 0
return xaxis2angle(parameters.xaa, parameters.xao)
def get_true_north(ifc_file):
def get_true_north(ifc_file: ifcopenshell.file) -> float:
"""Get an angle pointing to global true north
Anticlockwise is positive.
@@ -571,13 +577,11 @@ def get_true_north(ifc_file):
Always remember that true north is not a constant! (Unless you are working
in polar coordinates) This true north is only a reference value useful for
things like solar analysis on small sites (<1km). If you're after the north
that your surveyor is using, you're probably after ``get_grid_north``
that your surveyor is using, you're probably after :func:`get_grid_north`
instead.
:param ifc_file: The IFC file
:type ifc_file: ifcopenshell.file
:return: An angle to true north in decimal degrees
:rtype: float
"""
try:
for context in ifc_file.by_type("IfcGeometricRepresentationContext", include_subtypes=False):
@@ -588,15 +592,13 @@ def get_true_north(ifc_file):
return 0
def angle2xaxis(angle):
def angle2xaxis(angle: float) -> tuple[float, float]:
"""Converts an angle into an X axis abscissa and ordinate
The inverse of ``xaxis2angle``.
The inverse of :func:`xaxis2angle`.
:param angle: The angle in decimal degrees where anticlockwise is positive.
:type angle: float
:return: A tuple of X axis abscissa and ordinate
:rtype: tuple[float]
"""
angle_rad = math.radians(angle)
x = math.cos(angle_rad)
@@ -604,15 +606,13 @@ def angle2xaxis(angle):
return x, y
def angle2yaxis(angle):
def angle2yaxis(angle: float) -> tuple[float, float]:
"""Converts an angle into an Y axis abscissa and ordinate
The inverse of ``yaxis2angle``.
The inverse of :func:`yaxis2angle`.
:param angle: The angle in decimal degrees where anticlockwise is positive.
:type angle: float
:return: A tuple of Y axis abscissa and ordinate
:rtype: tuple[float]
"""
angle_rad = math.radians(angle)
x = -math.sin(angle_rad)