mirror of
https://github.com/IfcOpenShell/IfcOpenShell.git
synced 2026-09-21 14:13:41 +00:00
Update XYZ <-> ENH geolocation conversion after the clarifications made on how to handle scale in IFC4X3
This commit is contained in:
@@ -20,36 +20,42 @@ import math
|
|||||||
import numpy as np
|
import numpy as np
|
||||||
import ifcopenshell
|
import ifcopenshell
|
||||||
import ifcopenshell.util.unit
|
import ifcopenshell.util.unit
|
||||||
|
import ifcopenshell.util.element
|
||||||
|
from typing import NamedTuple, Optional, Union
|
||||||
|
|
||||||
|
|
||||||
def dms2dd(degrees, minutes, seconds, ms=0):
|
class HelmertTransformation(NamedTuple):
|
||||||
|
e: float
|
||||||
|
n: float
|
||||||
|
h: float
|
||||||
|
xaa: float
|
||||||
|
xao: float
|
||||||
|
scale: float
|
||||||
|
factor_x: float
|
||||||
|
factor_y: float
|
||||||
|
factor_z: float
|
||||||
|
|
||||||
|
|
||||||
|
def dms2dd(degrees: int, minutes: int, seconds: int, ms: int = 0) -> float:
|
||||||
"""Convert degrees, minutes, and (milli)seconds to decimal degrees
|
"""Convert degrees, minutes, and (milli)seconds to decimal degrees
|
||||||
|
|
||||||
:param degrees: The degrees component
|
:param degrees: The degrees component
|
||||||
:type degrees: int
|
|
||||||
:param minutes: The minutes component
|
:param minutes: The minutes component
|
||||||
:type minutes: int
|
|
||||||
:param seconds: The seconds component
|
:param seconds: The seconds component
|
||||||
:type seconds: int
|
|
||||||
:param ms: The milliseconds component
|
:param ms: The milliseconds component
|
||||||
:type ms: int
|
|
||||||
:return: The angle in decimal degrees.
|
:return: The angle in decimal degrees.
|
||||||
:rtype: float
|
|
||||||
"""
|
"""
|
||||||
dd = float(degrees) + float(minutes) / 60.0 + float(seconds) / (3600.0) + float(ms / 3600000000.0)
|
dd = float(degrees) + float(minutes) / 60.0 + float(seconds) / (3600.0) + float(ms / 3600000000.0)
|
||||||
return dd
|
return dd
|
||||||
|
|
||||||
|
|
||||||
def dd2dms(dd, use_ms=False):
|
def dd2dms(dd: float, use_ms: bool = False) -> Union[tuple[float, float, float, float], tuple[float, float, float]]:
|
||||||
"""Convert decimal degrees to degrees, minutes, and (milli)seconds format
|
"""Convert decimal degrees to degrees, minutes, and (milli)seconds format
|
||||||
|
|
||||||
:param dd: The decimal degrees
|
:param dd: The decimal degrees
|
||||||
:type dd: float
|
|
||||||
:param use_ms: True if to include milliseconds and false otherwise. Defaults to false.
|
:param use_ms: True if to include milliseconds and false otherwise. Defaults to false.
|
||||||
:type use_ms: bool
|
|
||||||
:return: The angle in a tuple of either 3 or 4 values, being degrees,
|
:return: The angle in a tuple of either 3 or 4 values, being degrees,
|
||||||
minutes, seconds, and optionally milliseconds.
|
minutes, seconds, and optionally milliseconds.
|
||||||
:rtype: tuple[float]
|
|
||||||
"""
|
"""
|
||||||
dd = float(dd)
|
dd = float(dd)
|
||||||
sign = 1 if dd >= 0 else -1
|
sign = 1 if dd >= 0 else -1
|
||||||
@@ -65,7 +71,20 @@ def dd2dms(dd, use_ms=False):
|
|||||||
return (int(degrees) * sign, int(minutes) * sign, int(seconds) * sign)
|
return (int(degrees) * sign, int(minutes) * sign, int(seconds) * sign)
|
||||||
|
|
||||||
|
|
||||||
def xyz2enh(x, y, z, eastings, northings, orthogonal_height, x_axis_abscissa, x_axis_ordinate, scale=None):
|
def xyz2enh(
|
||||||
|
x: float,
|
||||||
|
y: float,
|
||||||
|
z: float,
|
||||||
|
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,
|
||||||
|
) -> tuple[float, float, float]:
|
||||||
"""Manually convert local XYZ coordinates to map eastings, northings, and height
|
"""Manually convert local XYZ coordinates to map eastings, northings, and height
|
||||||
|
|
||||||
This function is for advanced users as it allows you to specify your own
|
This function is for advanced users as it allows you to specify your own
|
||||||
@@ -75,59 +94,29 @@ def xyz2enh(x, y, z, eastings, northings, orthogonal_height, x_axis_abscissa, x_
|
|||||||
you are applying your own temporary false origin (such as when federating
|
you are applying your own temporary false origin (such as when federating
|
||||||
models for digital twins of large cities).
|
models for digital twins of large cities).
|
||||||
|
|
||||||
No unit conversion is performed.
|
|
||||||
|
|
||||||
For most scenarios you should use ``auto_xyz2enh`` instead.
|
For most scenarios you should use ``auto_xyz2enh`` instead.
|
||||||
|
|
||||||
:param x: The X local engineering coordinate.
|
:param x: The X local engineering coordinate.
|
||||||
:type x: float
|
|
||||||
:param y: The Y local engineering coordinate.
|
:param y: The Y local engineering coordinate.
|
||||||
:type y: float
|
|
||||||
:param z: The Z local engineering coordinate.
|
:param z: The Z local engineering coordinate.
|
||||||
:type z: float
|
|
||||||
:param eastings: The eastings offset to apply.
|
:param eastings: The eastings offset to apply.
|
||||||
:type eastings: float
|
|
||||||
:param northings: The northings offset to apply.
|
:param northings: The northings offset to apply.
|
||||||
:type northings: float
|
|
||||||
:param orthogonal_height: The orthogonal height offset to apply.
|
: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
|
: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.
|
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
|
: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.
|
2D vector that points to the local X axis when in map coordinates.
|
||||||
:type x_axis_ordinate: float
|
:param scale: The unit scale such that local ordinate * scale = map
|
||||||
:param scale: The combined scale factor to convert from local coordinates
|
ordinate. E.g. if your project is in millimeters but your CRS is in
|
||||||
to map coordinates.
|
meters, your scale should be 0.001.
|
||||||
:type scale: float
|
:param factor_x: The combined scale factor for the X value to convert from
|
||||||
|
local coordinates to map coordinates. Your surveyor will typically know
|
||||||
|
this number and approximate it as a constant on a small site. Typically
|
||||||
|
factor_x and factor_y will be identical, and factor_z will be 1.
|
||||||
|
:param factor_y: Same but for the Y value.
|
||||||
|
:param factor_z: Same but for the Z value.
|
||||||
:return: A tuple of three ordinates representing the easting, northing and height.
|
:return: A tuple of three ordinates representing the easting, northing and height.
|
||||||
:rtype: tuple[float]
|
|
||||||
"""
|
"""
|
||||||
if scale is None:
|
|
||||||
scale = 1.0
|
|
||||||
rotation = math.atan2(x_axis_ordinate, x_axis_abscissa)
|
|
||||||
a = scale * math.cos(rotation)
|
|
||||||
b = scale * math.sin(rotation)
|
|
||||||
eastings = (a * x) - (b * y) + eastings
|
|
||||||
northings = (b * x) + (a * y) + northings
|
|
||||||
height = z + orthogonal_height
|
|
||||||
return (eastings, northings, height)
|
|
||||||
|
|
||||||
|
|
||||||
def xyz2enh_ifc4x3(
|
|
||||||
x,
|
|
||||||
y,
|
|
||||||
z,
|
|
||||||
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,
|
|
||||||
):
|
|
||||||
theta = math.atan2(x_axis_ordinate, x_axis_abscissa)
|
theta = math.atan2(x_axis_ordinate, x_axis_abscissa)
|
||||||
eastings = (scale * factor_x * math.cos(theta) * x) - (scale * factor_y * math.sin(theta) * y) + eastings
|
eastings = (scale * factor_x * math.cos(theta) * x) - (scale * factor_y * math.sin(theta) * y) + eastings
|
||||||
northings = (scale * factor_x * math.sin(theta) * x) + (scale * factor_y * math.cos(theta) * y) + northings
|
northings = (scale * factor_x * math.sin(theta) * x) + (scale * factor_y * math.cos(theta) * y) + northings
|
||||||
@@ -135,7 +124,9 @@ def xyz2enh_ifc4x3(
|
|||||||
return (eastings, northings, height)
|
return (eastings, northings, height)
|
||||||
|
|
||||||
|
|
||||||
def auto_xyz2enh(ifc_file, x, y, z):
|
def auto_xyz2enh(
|
||||||
|
ifc_file: ifcopenshell.file, x: float, y: float, z: float, should_return_in_map_units: bool = True
|
||||||
|
) -> tuple[float, float, float]:
|
||||||
"""Convert from local XYZ coordinates to global map coordinate eastings, northings, and heights
|
"""Convert from local XYZ coordinates to global map coordinate eastings, northings, and heights
|
||||||
|
|
||||||
The necessary georeferencing map conversion is automatically detected from
|
The necessary georeferencing map conversion is automatically detected from
|
||||||
@@ -147,63 +138,24 @@ def auto_xyz2enh(ifc_file, x, y, z):
|
|||||||
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
||||||
|
|
||||||
:param ifc_file: The IFC file
|
:param ifc_file: The IFC file
|
||||||
:type ifc_file: ifcopenshell.file
|
|
||||||
:param x: The X local engineering coordinate provided in project length units.
|
:param x: The X local engineering coordinate provided in project length units.
|
||||||
:type x: float
|
|
||||||
:param y: The Y local engineering coordinate provided in project length units.
|
:param y: The Y local engineering coordinate provided in project length units.
|
||||||
:type y: float
|
|
||||||
:param z: The Z local engineering coordinate provided in project length units.
|
:param z: The Z local engineering coordinate provided in project length units.
|
||||||
:type z: float
|
:param should_return_in_map_units: If true, the result is given in map units.
|
||||||
:return: The global map coordinate eastings, northings, and height 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]
|
:rtype: tuple[float]
|
||||||
"""
|
"""
|
||||||
conversion = None
|
parameters = get_helmert_transformation_parameters(ifc_file)
|
||||||
try:
|
if not parameters:
|
||||||
conversion = ifc_file.by_type("IfcMapConversion")
|
return x, y, z
|
||||||
except:
|
enh = xyz2enh(x, y, z, *parameters)
|
||||||
pass
|
if should_return_in_map_units:
|
||||||
|
return enh
|
||||||
if conversion:
|
return enh[0] / parameters.scale, enh[1] / parameters.scale, enh[2] / parameters.scale
|
||||||
conversion = conversion[0]
|
|
||||||
e = conversion.Eastings or 0
|
|
||||||
n = conversion.Northings or 0
|
|
||||||
h = conversion.OrthogonalHeight or 0
|
|
||||||
xaa = conversion.XAxisAbscissa or 0
|
|
||||||
xao = conversion.XAxisOrdinate or 0
|
|
||||||
scale = conversion.Scale or 1
|
|
||||||
map_unit = conversion.TargetCRS.MapUnit
|
|
||||||
else:
|
|
||||||
project = ifc_file.by_type("IfcProject")[0]
|
|
||||||
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion")
|
|
||||||
if not conversion:
|
|
||||||
return (x, y, z)
|
|
||||||
|
|
||||||
e = conversion.get("Eastings", None) or 0
|
|
||||||
n = conversion.get("Northings", None) or 0
|
|
||||||
h = conversion.get("OrthogonalHeight", None) or 0
|
|
||||||
xaa = conversion.get("XAxisAbscissa", None) or 0
|
|
||||||
xao = conversion.get("XAxisOrdinate", None) or 0
|
|
||||||
scale = conversion.get("Scale", None) or 1
|
|
||||||
map_unit = None
|
|
||||||
|
|
||||||
if not xaa and not xao:
|
|
||||||
xaa = 1.0
|
|
||||||
xao = 0.0
|
|
||||||
|
|
||||||
if map_unit:
|
|
||||||
# Warning! This definition has changed in IFC4X3 such that map_unit no
|
|
||||||
# longer affects unit conversion, only the Scale attribute affects unit
|
|
||||||
# conversion. TODO: consolidate once IFC4X3 confirmed.
|
|
||||||
project_unit = ifcopenshell.util.unit.get_project_unit(ifc_file, "LENGTHUNIT")
|
|
||||||
map_prefix = getattr(map_unit, "Prefix", None)
|
|
||||||
project_prefix = getattr(project_unit, "Prefix", None)
|
|
||||||
e = ifcopenshell.util.unit.convert(e, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
n = ifcopenshell.util.unit.convert(n, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
h = ifcopenshell.util.unit.convert(h, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
return xyz2enh(x, y, z, e, n, h, xaa, xao, scale)
|
|
||||||
|
|
||||||
|
|
||||||
def auto_enh2xyz(ifc_file, easting, northing, height):
|
def auto_enh2xyz(ifc_file, easting, northing, height, is_specified_in_map_units: bool = True):
|
||||||
"""Convert from global map coordinate eastings, northings, and heights to local XYZ coordinates
|
"""Convert from global map coordinate eastings, northings, and heights to local XYZ coordinates
|
||||||
|
|
||||||
The necessary georeferencing map conversion is automatically detected from
|
The necessary georeferencing map conversion is automatically detected from
|
||||||
@@ -215,63 +167,80 @@ def auto_enh2xyz(ifc_file, easting, northing, height):
|
|||||||
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
||||||
|
|
||||||
:param ifc_file: The IFC file
|
:param ifc_file: The IFC file
|
||||||
:type ifc_file: ifcopenshell.file
|
|
||||||
:param easting: The global easting map coordinate provided in map units.
|
:param easting: The global easting map coordinate provided in map units.
|
||||||
:type easting: float
|
|
||||||
:param northing: The global northing map coordinate provided in map units.
|
:param northing: The global northing map coordinate provided in map units.
|
||||||
:type northing: float
|
|
||||||
:param height: The global height map coordinate provided in map units.
|
:param height: The global height map coordinate provided in map units.
|
||||||
:type height: float
|
|
||||||
:return: The local engineering XYZ coordinates in project length units.
|
:return: The local engineering XYZ coordinates in project length units.
|
||||||
:rtype: tuple[float]
|
|
||||||
"""
|
"""
|
||||||
conversion = None
|
parameters = get_helmert_transformation_parameters(ifc_file)
|
||||||
try:
|
if not parameters:
|
||||||
conversion = ifc_file.by_type("IfcMapConversion")
|
return easting, northing, height
|
||||||
except:
|
if not is_specified_in_map_units:
|
||||||
pass
|
easting *= parameters.scale
|
||||||
|
northing *= parameters.scale
|
||||||
|
height *= parameters.scale
|
||||||
|
return enh2xyz(easting, northing, height, *parameters)
|
||||||
|
|
||||||
if conversion:
|
|
||||||
conversion = conversion[0]
|
def get_helmert_transformation_parameters(ifc_file: ifcopenshell.file) -> Optional[HelmertTransformation]:
|
||||||
e = conversion.Eastings or 0
|
"""Retrieves the parameters of a helmert transformation that represents a
|
||||||
n = conversion.Northings or 0
|
coordinate operation
|
||||||
h = conversion.OrthogonalHeight or 0
|
|
||||||
xaa = conversion.XAxisAbscissa or 0
|
This coordinate operation is typically what is used to convert between
|
||||||
xao = conversion.XAxisOrdinate or 0
|
local engineering coordinates and map coordinates.
|
||||||
scale = conversion.Scale or 1
|
|
||||||
map_unit = conversion.TargetCRS.MapUnit
|
:param ifc_file: The IFC model, typically containing an
|
||||||
else:
|
IfcCoordinateOperation such as an IfcMapConversion.
|
||||||
|
:return: The parameters of the transformation.
|
||||||
|
"""
|
||||||
|
if ifc_file.schema == "IFC2X3":
|
||||||
project = ifc_file.by_type("IfcProject")[0]
|
project = ifc_file.by_type("IfcProject")[0]
|
||||||
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion")
|
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion")
|
||||||
if not conversion:
|
if not conversion:
|
||||||
return (easting, northing, height)
|
return
|
||||||
|
|
||||||
e = conversion.get("Eastings", None) or 0
|
e = conversion.get("Eastings", None) or 0
|
||||||
n = conversion.get("Northings", None) or 0
|
n = conversion.get("Northings", None) or 0
|
||||||
h = conversion.get("OrthogonalHeight", None) or 0
|
h = conversion.get("OrthogonalHeight", None) or 0
|
||||||
xaa = conversion.get("XAxisAbscissa", None) or 0
|
xaa = conversion.get("XAxisAbscissa", None) or 0
|
||||||
xao = conversion.get("XAxisOrdinate", None) or 0
|
xao = conversion.get("XAxisOrdinate", None) or 0
|
||||||
scale = conversion.get("Scale", None) or 1
|
scale = conversion.get("Scale", None) or 1
|
||||||
map_unit = None
|
factor_x = factor_y = factor_z = 1
|
||||||
|
else:
|
||||||
|
conversion = ifc_file.by_type("IfcCoordinateOperation")
|
||||||
|
if not conversion:
|
||||||
|
return
|
||||||
|
conversion = conversion[0]
|
||||||
|
|
||||||
|
if conversion.is_a("IfcMapConversion"):
|
||||||
|
e = conversion.Eastings or 0
|
||||||
|
n = conversion.Northings or 0
|
||||||
|
h = conversion.OrthogonalHeight or 0
|
||||||
|
xaa = conversion.XAxisAbscissa or 0
|
||||||
|
xao = conversion.XAxisOrdinate or 0
|
||||||
|
scale = conversion.Scale or 1
|
||||||
|
if conversion.is_a() == "IfcMapConversionScaled":
|
||||||
|
factor_x = conversion.FactorX
|
||||||
|
factor_y = conversion.FactorY
|
||||||
|
factor_z = conversion.FactorZ
|
||||||
|
else:
|
||||||
|
factor_x = factor_y = factor_z = 1
|
||||||
|
elif conversion.is_a() == "IfcRigidOperation":
|
||||||
|
# TODO
|
||||||
|
e = conversion.FirstCoordinate
|
||||||
|
n = conversion.SecondCoordinate
|
||||||
|
h = conversion.Height or 0
|
||||||
|
xaa = 1.0
|
||||||
|
xao = 0.0
|
||||||
|
factor_x = factor_y = factor_z = 1
|
||||||
|
|
||||||
if not xaa and not xao:
|
if not xaa and not xao:
|
||||||
xaa = 1.0
|
xaa = 1.0
|
||||||
xao = 0.0
|
xao = 0.0
|
||||||
|
|
||||||
if map_unit:
|
return HelmertTransformation(e, n, h, xaa, xao, scale, factor_x, factor_y, factor_z)
|
||||||
# Warning! This definition has changed in IFC4X3 such that map_unit no
|
|
||||||
# longer affects unit conversion, only the Scale attribute affects unit
|
|
||||||
# conversion. TODO: consolidate once IFC4X3 confirmed.
|
|
||||||
project_unit = ifcopenshell.util.unit.get_project_unit(ifc_file, "LENGTHUNIT")
|
|
||||||
map_prefix = getattr(map_unit, "Prefix", None)
|
|
||||||
project_prefix = getattr(project_unit, "Prefix", None)
|
|
||||||
e = ifcopenshell.util.unit.convert(e, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
n = ifcopenshell.util.unit.convert(n, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
h = ifcopenshell.util.unit.convert(h, map_prefix, map_unit.Name, project_prefix, project_unit.Name)
|
|
||||||
return enh2xyz(easting, northing, height, e, n, h, xaa, xao, scale)
|
|
||||||
|
|
||||||
|
|
||||||
def auto_z2e(ifc_file, z):
|
def auto_z2e(ifc_file: ifcopenshell.file, z: float, should_return_in_map_units: bool = True) -> float:
|
||||||
"""Convert a Z coordinate to an elevation using model georeferencing data
|
"""Convert a Z coordinate to an elevation using model georeferencing data
|
||||||
|
|
||||||
The necessary georeferencing map conversion is automatically detected from
|
The necessary georeferencing map conversion is automatically detected from
|
||||||
@@ -283,66 +252,56 @@ def auto_z2e(ifc_file, z):
|
|||||||
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
https://www.buildingsmart.org/standards/bsi-standards/standards-library/
|
||||||
|
|
||||||
:param ifc_file: The IFC file
|
:param ifc_file: The IFC file
|
||||||
:type ifc_file: ifcopenshell.file
|
|
||||||
:param z: The Z local engineering coordinate provided in project length units.
|
:param z: The Z local engineering coordinate provided in project length units.
|
||||||
:type z: float
|
|
||||||
:return: The elevation in project length units.
|
:return: The elevation in project length units.
|
||||||
:rtype: float
|
|
||||||
"""
|
"""
|
||||||
conversion = None
|
parameters = get_helmert_transformation_parameters(ifc_file)
|
||||||
try:
|
if not parameters:
|
||||||
conversion = ifc_file.by_type("IfcMapConversion")
|
return z
|
||||||
except:
|
e = z2e(z, parameters.h, parameters.scale, parameters.factor_z)
|
||||||
pass
|
if should_return_in_map_units:
|
||||||
|
return e
|
||||||
if conversion and not conversion[0].OrthogonalHeight:
|
return e / parameters.scale
|
||||||
conversion = conversion[0]
|
|
||||||
h = conversion.OrthogonalHeight
|
|
||||||
map_unit = conversion.TargetCRS.MapUnit
|
|
||||||
else:
|
|
||||||
project = ifc_file.by_type("IfcProject")[0]
|
|
||||||
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion")
|
|
||||||
if not conversion:
|
|
||||||
return z
|
|
||||||
|
|
||||||
h = conversion.get("OrthogonalHeight", None) or 0
|
|
||||||
map_unit = None
|
|
||||||
|
|
||||||
if map_unit:
|
|
||||||
# Warning! This definition has changed in IFC4X3 such that map_unit no
|
|
||||||
# longer affects unit conversion, only the Scale attribute affects unit
|
|
||||||
# conversion. TODO: consolidate once IFC4X3 confirmed.
|
|
||||||
project_unit = ifcopenshell.util.unit.get_project_unit(ifc_file, "LENGTHUNIT")
|
|
||||||
h = ifcopenshell.util.unit.convert(
|
|
||||||
h,
|
|
||||||
getattr(map_unit, "Prefix", None),
|
|
||||||
map_unit.Name,
|
|
||||||
getattr(project_unit, "Prefix", None),
|
|
||||||
project_unit.Name,
|
|
||||||
)
|
|
||||||
return z2e(z, h)
|
|
||||||
|
|
||||||
|
|
||||||
def z2e(z, h):
|
def z2e(z: float, orthogonal_height: float = 0.0, scale: float = 1.0, factor_z: float = 1.0) -> float:
|
||||||
"""Manually convert a Z coordinate to an elevation
|
"""Manually convert a Z coordinate to a map elevation
|
||||||
|
|
||||||
This function is for advanced users as it allows you to specify your own
|
This function is for advanced users as it allows you to specify your own
|
||||||
orthogonal height offset.
|
orthogonal height offset and transformation parameters.
|
||||||
|
|
||||||
For most scenarios you should use ``auto_z2e`` instead.
|
For most scenarios you should use ``auto_z2e`` instead.
|
||||||
|
|
||||||
:param z: The Z local engineering coordinate provided in project length units.
|
:param z: The Z local engineering coordinate provided in project length units.
|
||||||
:type z: float
|
:param orthogonal_height: The orthogonal height offset to apply.
|
||||||
:param h: The orthogonal height offset in project length units.
|
:param scale: The unit scale such that local ordinate * scale = map
|
||||||
:type h: float
|
ordinate. E.g. if your project is in millimeters but your CRS is in
|
||||||
:return: The elevation in project length units.
|
meters, your scale should be 0.001.
|
||||||
:rtype: float
|
:param factor_x: The combined scale factor for the Z value to convert from
|
||||||
|
local coordinates to map coordinates. Your surveyor will typically know
|
||||||
|
this number and approximate it as a constant on a small site. This is
|
||||||
|
typically just 1.0, as average combined scale factors usually only
|
||||||
|
affect the XY axes.
|
||||||
|
:return: The elevation in map units.
|
||||||
"""
|
"""
|
||||||
return z + h
|
return (scale * factor_z * z) + orthogonal_height
|
||||||
|
|
||||||
|
|
||||||
def enh2xyz(e, n, h, eastings, northings, orthogonal_height, x_axis_abscissa, x_axis_ordinate, scale=None):
|
def enh2xyz(
|
||||||
"""Manually convert map eastings, northings, and height to local XYZ coordinates
|
e: float,
|
||||||
|
n: float,
|
||||||
|
h: float,
|
||||||
|
eastings: float = 0.0,
|
||||||
|
northings: float = 0.0,
|
||||||
|
orthogonal_height: float = 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,
|
||||||
|
) -> tuple[float, float, float]:
|
||||||
|
"""Manually convert map eastings, northings, and height to local XYZ coordinates
|
||||||
|
|
||||||
This function is for advanced users as it allows you to specify your own
|
This function is for advanced users as it allows you to specify your own
|
||||||
helmert transformation parameters (i.e. those typically stored in
|
helmert transformation parameters (i.e. those typically stored in
|
||||||
@@ -351,42 +310,35 @@ def enh2xyz(e, n, h, eastings, northings, orthogonal_height, x_axis_abscissa, x_
|
|||||||
you are applying your own temporary false origin (such as when federating
|
you are applying your own temporary false origin (such as when federating
|
||||||
models for digital twins of large cities).
|
models for digital twins of large cities).
|
||||||
|
|
||||||
No unit conversion is performed.
|
|
||||||
|
|
||||||
For most scenarios you should use ``auto_enh2xyz`` instead.
|
For most scenarios you should use ``auto_enh2xyz`` instead.
|
||||||
|
|
||||||
:param e: The global easting map coordinate.
|
:param e: The global easting map coordinate.
|
||||||
:type e: float
|
|
||||||
:param n: The global northing map coordinate.
|
:param n: The global northing map coordinate.
|
||||||
:type n: float
|
|
||||||
:param h: The global height map coordinate.
|
:param h: The global height map coordinate.
|
||||||
:type h: float
|
|
||||||
:param eastings: The eastings offset to apply.
|
:param eastings: The eastings offset to apply.
|
||||||
:type eastings: float
|
|
||||||
:param northings: The northings offset to apply.
|
:param northings: The northings offset to apply.
|
||||||
:type northings: float
|
|
||||||
:param orthogonal_height: The orthogonal height offset to apply.
|
: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
|
: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.
|
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
|
: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.
|
2D vector that points to the local X axis when in map coordinates.
|
||||||
:type x_axis_ordinate: float
|
:param scale: The unit scale such that local ordinate * scale = map
|
||||||
:param scale: The combined scale factor to convert from local coordinates
|
ordinate. E.g. if your project is in millimeters but your CRS is in
|
||||||
to map coordinates.
|
meters, your scale should be 0.001.
|
||||||
:type scale: float
|
:param factor_x: The combined scale factor for the X value to convert from
|
||||||
|
local coordinates to map coordinates. Your surveyor will typically know
|
||||||
|
this number and approximate it as a constant on a small site. Typically
|
||||||
|
factor_x and factor_y will be identical, and factor_z will be 1.
|
||||||
|
:param factor_y: Same but for the Y value.
|
||||||
|
:param factor_z: Same but for the Z value.
|
||||||
:return: A tuple of three ordinates representing XYZ.
|
:return: A tuple of three ordinates representing XYZ.
|
||||||
:rtype: tuple[float]
|
|
||||||
"""
|
"""
|
||||||
if scale is None:
|
theta = math.atan2(x_axis_ordinate, x_axis_abscissa)
|
||||||
scale = 1.0
|
sint = math.sin(theta)
|
||||||
rotation = math.atan2(x_axis_ordinate, x_axis_abscissa)
|
cost = math.cos(theta)
|
||||||
a = scale * math.cos(rotation)
|
x = (((e - eastings) * cost) + ((n - northings) * sint)) / (scale * factor_x)
|
||||||
b = scale * math.sin(rotation)
|
y = (((eastings - e) * sint) + ((n - northings) * cost)) / (scale * factor_y)
|
||||||
x = ((b * n) - (b * northings) - (a * eastings) + (a * e)) / ((a * a) + (b * b))
|
z = ((h - orthogonal_height) / scale) / factor_z
|
||||||
y = ((a * n) - (a * northings) + (b * eastings) - (b * e)) / ((a * a) + (b * b))
|
|
||||||
z = h - orthogonal_height
|
|
||||||
return (x, y, z)
|
return (x, y, z)
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
@@ -31,10 +31,7 @@ class TestEditGeoreferencing(test.bootstrap.IFC4):
|
|||||||
ifcopenshell.api.georeference.edit_georeferencing(
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
self.file,
|
self.file,
|
||||||
projected_crs={"Name": "EPSG:7856"},
|
projected_crs={"Name": "EPSG:7856"},
|
||||||
map_conversion={
|
map_conversion={"Eastings": 123.45, "Northings": 234.56},
|
||||||
"Eastings": 123.45,
|
|
||||||
"Northings": 234.56,
|
|
||||||
},
|
|
||||||
)
|
)
|
||||||
crs = self.file.by_type("IfcProjectedCRS")[0]
|
crs = self.file.by_type("IfcProjectedCRS")[0]
|
||||||
assert crs.Name == "EPSG:7856"
|
assert crs.Name == "EPSG:7856"
|
||||||
@@ -47,9 +44,9 @@ class TestEditGeoreferencing(test.bootstrap.IFC4):
|
|||||||
model = ifcopenshell.api.context.add_context(self.file, "Model")
|
model = ifcopenshell.api.context.add_context(self.file, "Model")
|
||||||
plan = ifcopenshell.api.context.add_context(self.file, "Plan")
|
plan = ifcopenshell.api.context.add_context(self.file, "Plan")
|
||||||
ifcopenshell.api.georeference.add_georeferencing(self.file)
|
ifcopenshell.api.georeference.add_georeferencing(self.file)
|
||||||
ifcopenshell.api.georeference.edit_georeferencing(self.file, true_north=[0., 1.])
|
ifcopenshell.api.georeference.edit_georeferencing(self.file, true_north=[0.0, 1.0])
|
||||||
assert model.TrueNorth[0] == (0., 1., 0.)
|
assert model.TrueNorth[0] == (0.0, 1.0, 0.0)
|
||||||
assert plan.TrueNorth[0] == (0., 1.)
|
assert plan.TrueNorth[0] == (0.0, 1.0)
|
||||||
|
|
||||||
|
|
||||||
class TestEditGeoreferencingIFC2X3(test.bootstrap.IFC2X3):
|
class TestEditGeoreferencingIFC2X3(test.bootstrap.IFC2X3):
|
||||||
@@ -59,10 +56,7 @@ class TestEditGeoreferencingIFC2X3(test.bootstrap.IFC2X3):
|
|||||||
ifcopenshell.api.georeference.edit_georeferencing(
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
self.file,
|
self.file,
|
||||||
projected_crs={"Name": "EPSG:7856"},
|
projected_crs={"Name": "EPSG:7856"},
|
||||||
map_conversion={
|
map_conversion={"Eastings": 123.45, "Northings": 234.56},
|
||||||
"Eastings": 123.45,
|
|
||||||
"Northings": 234.56,
|
|
||||||
},
|
|
||||||
)
|
)
|
||||||
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion", verbose=True)
|
conversion = ifcopenshell.util.element.get_pset(project, "ePSet_MapConversion", verbose=True)
|
||||||
crs = ifcopenshell.util.element.get_pset(project, "ePSet_ProjectedCRS", verbose=True)
|
crs = ifcopenshell.util.element.get_pset(project, "ePSet_ProjectedCRS", verbose=True)
|
||||||
|
|||||||
@@ -19,6 +19,9 @@
|
|||||||
import pytest
|
import pytest
|
||||||
import numpy as np
|
import numpy as np
|
||||||
import test.bootstrap
|
import test.bootstrap
|
||||||
|
import ifcopenshell.api.root
|
||||||
|
import ifcopenshell.api.context
|
||||||
|
import ifcopenshell.api.georeference
|
||||||
import ifcopenshell.util.geolocation as subject
|
import ifcopenshell.util.geolocation as subject
|
||||||
|
|
||||||
|
|
||||||
@@ -29,18 +32,109 @@ class TestXYZ2ENH(test.bootstrap.IFC4):
|
|||||||
assert subject.xyz2enh(0, 0, 0, 1, 2, 3, 0, 1) == (1, 2, 3)
|
assert subject.xyz2enh(0, 0, 0, 1, 2, 3, 0, 1) == (1, 2, 3)
|
||||||
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 1, 0), (2, 3, 3))
|
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 1, 0), (2, 3, 3))
|
||||||
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 1, 0, 2), (3, 4, 3))
|
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 1, 0, 2), (3, 4, 3))
|
||||||
|
assert np.allclose(subject.xyz2enh(1, 1, 1, 1, 2, 3, 1, 0, 2, 2, 3, 4), (5, 8, 11))
|
||||||
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 0, 1), (0, 3, 3))
|
assert np.allclose(subject.xyz2enh(1, 1, 0, 1, 2, 3, 0, 1), (0, 3, 3))
|
||||||
|
|
||||||
|
|
||||||
class TestXYZ2ENHIfc4X3(test.bootstrap.IFC4):
|
class TestENH2XYZ(test.bootstrap.IFC4):
|
||||||
def test_converting_from_a_local_xyz_point_to_a_global_easting_northing_height(self):
|
def test_converting_from_a_global_easting_northing_height_to_a_local_xyz_point(self):
|
||||||
assert subject.xyz2enh_ifc4x3(0, 0, 0, 0, 0, 0, 1, 0) == (0, 0, 0)
|
assert subject.enh2xyz(0, 0, 0, 0, 0, 0, 1, 0) == (0, 0, 0)
|
||||||
assert subject.xyz2enh_ifc4x3(0, 0, 0, 1, 2, 3, 1, 0) == (1, 2, 3)
|
assert subject.enh2xyz(1, 2, 3, 1, 2, 3, 1, 0) == (0, 0, 0)
|
||||||
assert subject.xyz2enh_ifc4x3(0, 0, 0, 1, 2, 3, 0, 1) == (1, 2, 3)
|
assert subject.enh2xyz(1, 2, 3, 1, 2, 3, 0, 1) == (0, 0, 0)
|
||||||
assert np.allclose(subject.xyz2enh_ifc4x3(1, 1, 0, 1, 2, 3, 1, 0), (2, 3, 3))
|
assert np.allclose(subject.enh2xyz(2, 3, 3, 1, 2, 3, 1, 0), (1, 1, 0))
|
||||||
assert np.allclose(subject.xyz2enh_ifc4x3(1, 1, 0, 1, 2, 3, 1, 0, 2), (3, 4, 3))
|
assert np.allclose(subject.enh2xyz(3, 4, 3, 1, 2, 3, 1, 0, 2), (1, 1, 0))
|
||||||
assert np.allclose(subject.xyz2enh_ifc4x3(1, 1, 1, 1, 2, 3, 1, 0, 2, 2, 3, 4), (5, 8, 11))
|
assert np.allclose(subject.enh2xyz(5, 8, 11, 1, 2, 3, 1, 0, 2, 2, 3, 4), (1, 1, 1))
|
||||||
assert np.allclose(subject.xyz2enh_ifc4x3(1, 1, 0, 1, 2, 3, 0, 1), (0, 3, 3))
|
assert np.allclose(subject.enh2xyz(0, 3, 3, 1, 2, 3, 0, 1), (1, 1, 0))
|
||||||
|
|
||||||
|
|
||||||
|
class TestZ2E(test.bootstrap.IFC4):
|
||||||
|
def test_converting_from_a_local_z_to_a_global_elevation(self):
|
||||||
|
assert subject.z2e(0) == 0
|
||||||
|
assert subject.z2e(0, 0, 1, 1) == 0
|
||||||
|
assert subject.z2e(0, 2) == 2
|
||||||
|
assert subject.z2e(1, 2) == 3
|
||||||
|
assert subject.z2e(1000, 2, 0.001) == 3
|
||||||
|
assert np.isclose(subject.z2e(1000, 2, 0.001, 0.9), 2.9)
|
||||||
|
|
||||||
|
|
||||||
|
class TestAutoXYZ2ENH(test.bootstrap.IFC4):
|
||||||
|
def test_no_georeferencing(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
assert subject.auto_xyz2enh(self.file, 0, 0, 0) == (0, 0, 0)
|
||||||
|
assert subject.auto_xyz2enh(self.file, 1, 2, 3) == (1, 2, 3)
|
||||||
|
|
||||||
|
def test_map_conversion(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
ifcopenshell.api.context.add_context(self.file, "Model")
|
||||||
|
ifcopenshell.api.georeference.add_georeferencing(self.file)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file,
|
||||||
|
projected_crs={"Name": "EPSG:7856"},
|
||||||
|
map_conversion={"Eastings": 1, "Northings": 2, "OrthogonalHeight": 3},
|
||||||
|
)
|
||||||
|
assert subject.auto_xyz2enh(self.file, 0, 0, 0) == (1, 2, 3)
|
||||||
|
assert subject.auto_xyz2enh(self.file, 1, 3, 5) == (2, 5, 8)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(self.file, map_conversion={"Scale": 0.001})
|
||||||
|
assert subject.auto_xyz2enh(self.file, 1000, 1000, 0) == (2, 3, 3)
|
||||||
|
assert subject.auto_xyz2enh(self.file, 1000, 1000, 0, should_return_in_map_units=False) == (2000, 3000, 3000)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file, map_conversion={"XAxisAbscissa": 0, "XAxisOrdinate": 1}
|
||||||
|
)
|
||||||
|
assert np.allclose(subject.auto_xyz2enh(self.file, 1000, 1000, 0), (0, 3, 3))
|
||||||
|
|
||||||
|
|
||||||
|
class TestAutoENH2XYZ(test.bootstrap.IFC4):
|
||||||
|
def test_no_georeferencing(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
assert subject.auto_enh2xyz(self.file, 0, 0, 0) == (0, 0, 0)
|
||||||
|
assert subject.auto_enh2xyz(self.file, 1, 2, 3) == (1, 2, 3)
|
||||||
|
|
||||||
|
def test_map_conversion(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
ifcopenshell.api.context.add_context(self.file, "Model")
|
||||||
|
ifcopenshell.api.georeference.add_georeferencing(self.file)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file,
|
||||||
|
projected_crs={"Name": "EPSG:7856"},
|
||||||
|
map_conversion={"Eastings": 1, "Northings": 2, "OrthogonalHeight": 3},
|
||||||
|
)
|
||||||
|
assert subject.auto_enh2xyz(self.file, 1, 2, 3) == (0, 0, 0)
|
||||||
|
assert subject.auto_enh2xyz(self.file, 2, 5, 8) == (1, 3, 5)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(self.file, map_conversion={"Scale": 0.001})
|
||||||
|
assert np.allclose(subject.auto_enh2xyz(self.file, 2, 3, 3), (1000, 1000, 0))
|
||||||
|
assert np.allclose(
|
||||||
|
subject.auto_enh2xyz(self.file, 2000, 3000, 3000, is_specified_in_map_units=False), (1000, 1000, 0)
|
||||||
|
)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file, map_conversion={"XAxisAbscissa": 0, "XAxisOrdinate": 1}
|
||||||
|
)
|
||||||
|
assert np.allclose(subject.auto_enh2xyz(self.file, 0, 3, 3), (1000, 1000, 0))
|
||||||
|
|
||||||
|
|
||||||
|
class TestAutoZ2E(test.bootstrap.IFC4):
|
||||||
|
def test_no_georeferencing(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
assert subject.auto_z2e(self.file, 0) == 0
|
||||||
|
assert subject.auto_z2e(self.file, 1) == 1
|
||||||
|
|
||||||
|
def test_map_conversion(self):
|
||||||
|
ifcopenshell.api.root.create_entity(self.file, ifc_class="IfcProject")
|
||||||
|
ifcopenshell.api.context.add_context(self.file, "Model")
|
||||||
|
ifcopenshell.api.georeference.add_georeferencing(self.file)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file,
|
||||||
|
projected_crs={"Name": "EPSG:7856"},
|
||||||
|
map_conversion={"Eastings": 1, "Northings": 2, "OrthogonalHeight": 3},
|
||||||
|
)
|
||||||
|
assert subject.auto_z2e(self.file, 0) == 3
|
||||||
|
assert subject.auto_z2e(self.file, 5) == 8
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(self.file, map_conversion={"Scale": 0.001})
|
||||||
|
assert np.isclose(subject.auto_z2e(self.file, 0), 3)
|
||||||
|
assert np.isclose(subject.auto_z2e(self.file, 0, should_return_in_map_units=False), 3000)
|
||||||
|
ifcopenshell.api.georeference.edit_georeferencing(
|
||||||
|
self.file, map_conversion={"XAxisAbscissa": 0, "XAxisOrdinate": 1}
|
||||||
|
)
|
||||||
|
assert np.isclose(subject.auto_z2e(self.file, 0), 3)
|
||||||
|
|
||||||
|
|
||||||
class TestLocal2Global(test.bootstrap.IFC4):
|
class TestLocal2Global(test.bootstrap.IFC4):
|
||||||
|
|||||||
Reference in New Issue
Block a user