more work on transitions for MEP

- support transitions from and to circle profiles
- reworked transition length algorithm now it should be more accurate
- added support for creating transitions between profiles that are parallel but not collinear
This commit is contained in:
Andrej730
2023-08-21 17:22:20 +05:00
parent 941ede117c
commit 2d3295ec54
4 changed files with 405 additions and 80 deletions
@@ -29,6 +29,7 @@ import ifcopenshell.util.system
import ifcopenshell.util.element
import ifcopenshell.util.representation
import mathutils.geometry
import numpy as np
import blenderbim.bim.handler
import blenderbim.core.type
import blenderbim.core.root
@@ -347,6 +348,8 @@ class MEPGenerator:
There lies the problem that it won't be
able to identify the fittings that were not yet connected to any segments yet.
"""
# TODO: check angle, start, end and offset for transitions
if not isinstance(segment_or_segments, collections.abc.Iterable):
segments = [segment_or_segments]
ports = [port_or_ports]
@@ -436,7 +439,7 @@ class MEPGenerator:
if element_type is None:
skipped_the_occurrence = True
break
fitting_data.append((element_type, port.PredefinedType, port.SystemType))
# if we skipped the occurrence we still can other occurrences
@@ -624,10 +627,13 @@ class MEPAddTransition(bpy.types.Operator, tool.Ifc.Operator):
start_axis = tool.Model.get_flow_segment_axis(start_object)
end_axis = tool.Model.get_flow_segment_axis(end_object)
start_object_rotation = start_object.matrix_world.to_quaternion()
start_object_z_basis = start_object_rotation.to_matrix().col[2] # z basis vector
keep_only_z_axis = lambda p_ws: p_ws.dot(start_object_z_basis) * start_object_z_basis
# TODO: support cases when segments are partially or completely overlapping each other
if not tool.Cad.are_edges_collinear(start_axis, end_axis):
self.report({"ERROR"}, f"Failed to add transition - non collinear segments are not yet supported.")
if not tool.Cad.are_edges_parallel(start_axis, end_axis):
self.report({"ERROR"}, f"Failed to add transition - segments are not parallel.")
return {"CANCELLED"}
start_segment_data = MEPGenerator().get_segment_data(start_element)
@@ -647,42 +653,78 @@ class MEPAddTransition(bpy.types.Operator, tool.Ifc.Operator):
(start_segment_data["start_point"], start_segment_data["end_point"]),
(end_segment_data["start_point"], end_segment_data["end_point"]),
)
# figure profile offset
base_transition_dir = keep_only_z_axis(end_point - start_point).normalized()
flip_profile_offset = base_transition_dir.dot(start_object_z_basis) < 0
if tool.Cad.are_edges_collinear(start_axis, end_axis):
profile_offset = None
else:
to_start_object_space = start_object_rotation.inverted()
profile_offset = (
(to_start_object_space @ end_object.location) - (to_start_object_space @ start_object.location)
).xy
if tool.Cad.is_x(profile_offset.length_squared, 0):
profile_offset = None
else:
profile_offset = profile_offset / si_conversion
if flip_profile_offset:
profile_offset *= V(1, -1)
# world space profile offset
profile_offset_ws = (
start_object_rotation @ (profile_offset * si_conversion).to_3d() if profile_offset else V(0, 0, 0)
)
# will need entire_length to check that transition length fill fit
first_segment_start, second_segment_end = [
p for p in (
start_segment_data["start_point"],
start_segment_data["end_point"],
end_segment_data["start_point"],
end_segment_data["end_point"])
p
for p in (
start_segment_data["start_point"],
start_segment_data["end_point"],
end_segment_data["start_point"],
end_segment_data["end_point"],
)
if p not in (start_point, end_point)
]
entire_length = (first_segment_start - second_segment_end).length
transition_dir = (end_point - start_point).normalized()
# can't rely on (end_point-start_point) here because
# transition might change the segments length and therefore direction will be changed
segments_dir = (start_point - first_segment_start).normalized()
start_port = points_ports_map[start_point]
end_port = points_ports_map[end_point]
# add transition representation
builder = ShapeBuilder(ifc_file)
rep, transition_data = builder.mep_transition_shape(
start_element, end_element, self.start_length / si_conversion, self.end_length / si_conversion
start_element,
end_element,
self.start_length / si_conversion,
self.end_length / si_conversion,
profile_offset=profile_offset,
)
if not rep:
self.report({"ERROR"}, f"Failed to add transition - this kind of profiles is not yet supported.")
return {"CANCELLED"}
# TODO: test it
full_transition_length = transition_data["full_transition_length"] * si_conversion
if full_transition_length >= entire_length:
self.report({"ERROR"}, f"Failed to add transition - transition length is larger the segments and the distance between them.")
# TODO: handle the case without creating representation in the first place?
self.report(
{"ERROR"},
f"Failed to add transition - transition length is larger the segments and the distance between them.\n"
+ f"Transition length: {full_transition_length:.2f}m, segments length: {entire_length:.2f}m",
)
# TODO: handle the case without creating a representation in the first place?
ifcopenshell.api.run("geometry.remove_representation", ifc_file, representation=rep)
return {"CANCELLED"}
middle_point = (start_point + end_point) / 2
start_segment_extend_point = middle_point - transition_dir * full_transition_length / 2
end_segment_extend_point = middle_point + transition_dir * full_transition_length / 2
middle_point = keep_only_z_axis((start_point + end_point) / 2 - start_point) + start_point
start_segment_extend_point = middle_point - segments_dir * full_transition_length / 2
end_segment_extend_point = middle_point + segments_dir * full_transition_length / 2 + profile_offset_ws
transition_dir = keep_only_z_axis(end_segment_extend_point - start_segment_extend_point).normalized()
DumbProfileJoiner().join_E(start_object, start_segment_extend_point)
DumbProfileJoiner().join_E(end_object, end_segment_extend_point)
@@ -691,6 +733,10 @@ class MEPAddTransition(bpy.types.Operator, tool.Ifc.Operator):
)
transition_type = fitting_data["fitting_type"] if fitting_data else None
if transition_type:
# TODO: handle the case without creating a representation in the first place?
ifcopenshell.api.run("geometry.remove_representation", ifc_file, representation=rep)
start_port_match = fitting_data["start_port_match"] if fitting_data else True
if not transition_type:
@@ -722,8 +768,9 @@ class MEPAddTransition(bpy.types.Operator, tool.Ifc.Operator):
# adjust transition segment rotation and location
transition_obj.matrix_world = start_object.matrix_world
context.view_layer.update()
transition_obj_dir = tool.Cad.get_edge_direction(tool.Model.get_flow_segment_axis(transition_obj))
direction_match = tool.Cad.are_vectors_equal(transition_obj_dir, transition_dir)
direction_match = tool.Cad.are_vectors_equal(transition_dir, transition_obj_dir)
# if there are no mismatches or everything matches up we don't need to flip the transition
if start_port_match != direction_match:
+7 -11
View File
@@ -236,20 +236,16 @@ class Cad:
return (edge[1] - edge[0]).normalized()
@classmethod
def are_edges_collinear(cls, edge1, edge2):
def is_point_on_line(p, edge):
a1, a2 = edge
# comparing slopes between PA1 and A2A1
# using cross multiplication to avoid division by zero
return cls.is_x((p.y - a1.y) * (a2.x - a1.x), (a2.y - a1.y) * (p.x - a1.x))
def are_edges_parallel(cls, edge1, edge2):
edge1_dir = edge1[1] - edge1[0]
edge2_dir = edge2[1] - edge2[0]
return cls.is_x(edge1_dir.cross(edge2_dir).length_squared, 0)
if cls.is_x(edge1_dir.cross(edge2_dir).length_squared, 0): # check they are parallel
if is_point_on_line(edge1[0], edge2) or is_point_on_line(edge1[1], edge2):
return True
return False
@classmethod
def are_edges_collinear(cls, edge1, edge2):
if not cls.are_edges_parallel(edge1, edge2):
return False
return cls.are_edges_parallel((edge2[0], edge1[0]), edge2)
@classmethod
def closest_points(cls, edge1, edge2) -> bool:
+58
View File
@@ -0,0 +1,58 @@
# BlenderBIM Add-on - OpenBIM Blender Add-on
# Copyright (C) 2023 Dion Moult <dion@thinkmoult.com>, @Andrej730
#
# This file is part of BlenderBIM Add-on.
#
# BlenderBIM Add-on is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# BlenderBIM Add-on is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with BlenderBIM Add-on. If not, see <http://www.gnu.org/licenses/>.
from test.bim.bootstrap import NewFile
from blenderbim.tool.cad import Cad as subject
from mathutils import Vector
V = lambda *x: Vector([float(i) for i in x])
class TestAreEdgesCollinear(NewFile):
def test_run(self):
# fmt: off
# Parallel edges but not collinear (different z-coordinates)
assert not subject.are_edges_collinear(
(V(-1,0,-1), V(1,0,-1)),
(V(-1,0,1), V(1,0,1))
)
# One edge is just a point and the other is a line segment.
assert not subject.are_edges_collinear(
(V(1,-1,0), V(1,-1,0)),
(V(-1,1,0), V(1,1,0))
)
# Both edges are collinear and overlap.
assert subject.are_edges_collinear(
(V(0,0,0), V(2,2,2)),
(V(1,1,1), V(3,3,3))
)
# Both edges are collinear but don't overlap.
assert subject.are_edges_collinear(
(V(0,0,0), V(1,1,1)),
(V(2,2,2), V(3,3,3))
)
# Edges are not parallel and not collinear.
assert not subject.are_edges_collinear(
(V(0,0,0), V(1,1,1)),
(V(0,1,0), V(1,0,1))
)
# fmt: on
@@ -19,12 +19,17 @@
import collections
import ifcopenshell
import ifcopenshell.api
from math import cos, sin, pi, tan, radians
from math import cos, sin, pi, tan, radians, degrees, atan, sqrt
from mathutils import Vector, Matrix
from itertools import chain
V = lambda *x: Vector([float(i) for i in x])
sign = lambda x: x and (1, -1)[x < 0]
PRECISION = 1.0e-5
is_x = lambda value, x: (x + PRECISION) > value > (x - PRECISION)
round_to_precision = lambda x, si_conversion: round(x * si_conversion, 5) / si_conversion
round_vector_to_precision = lambda v, si_conversion: Vector([round_to_precision(i, si_conversion) for i in v])
# Note: using ShapeBuilder try not to reuse IFC elements in the process
# otherwise you might run into situation where builder.mirror or other operation
@@ -840,7 +845,9 @@ class ShapeBuilder:
return face_set
def mep_transition_shape(self, start_segment, end_segment, start_length, end_length, angle=30.0):
def mep_transition_shape(
self, start_segment, end_segment, start_length, end_length, angle=30.0, profile_offset=None
):
"""
returns tuple of Model/Body/MODEL_VIEW IfcRepresentation and transition shape data
"""
@@ -853,68 +860,285 @@ class ShapeBuilder:
if material and material.is_a("IfcMaterialProfileSet") and len(material.MaterialProfiles) == 1:
return material.MaterialProfiles[0].Profile
def get_circle_points(radius, segments=16):
"""starting from (R,0), going counter-clockwise"""
angle_d = 2 * pi / segments
verts = []
for i in range(segments):
angle = angle_d * i
verts.append(V(cos(angle), sin(angle), 0) * radius)
return verts
def get_rectangle_points(dim):
"""Starting from (+X/2, +Y/2) going counter-clockwise"""
dim = dim / 2
points = [
dim * V(1, 1, 0),
dim * V(-1, 1, 0),
dim * V(-1, -1, 0),
dim * V(1, -1, 0),
]
return points
# TODO: support more profiles
def get_dim(profile, depth):
if profile.is_a("IfcRectangleProfileDef"):
return V(profile.XDim / 2, profile.YDim / 2, depth)
elif profile.is_a("IfcCircleProfileDef"):
return V(profile.Radius, profile.Radius, depth)
return None
def get_profile_faceset(points, length, offset=None):
# prevent mutating arguments, deepcopy doesn't work
start_points = [p.copy() if not offset else (p + offset) for p in points]
end_points = [p.copy() for p in start_points]
for p in end_points:
p.z += length
points = start_points + end_points
faces = []
n_verts = len(start_points)
last_vert_i = n_verts - 1
for i in range(last_vert_i):
face = (i, i + 1, n_verts + i + 1, n_verts + i)
faces.append(face)
faces.append((last_vert_i, 0, n_verts + 0, n_verts + last_vert_i)) # close the loop
# if there is offset we put a cap at the end
# otherwise at the start
if offset:
faces.append(tuple(range(n_verts, n_verts * 2)))
else:
faces.append(tuple(reversed(range(n_verts))))
face_set = self.polygonal_face_set(points, faces)
return face_set
start_profile = get_profile(start_segment)
end_profile = get_profile(end_segment)
# TODO: support more profiles
if not start_profile.is_a("IfcRectangleProfileDef") or not end_profile.is_a("IfcRectangleProfileDef"):
# Non rectangular profiles are not yet supported
start_half_dim = get_dim(start_profile, start_length)
end_half_dim = get_dim(end_profile, end_length)
# if profile types are not supported
if not start_half_dim or not end_half_dim:
return None, None
start_half_dim = V(start_profile.XDim / 2, start_profile.YDim / 2, start_length)
end_half_dim = V(end_profile.XDim / 2, end_profile.YDim / 2, end_length)
transition_items = []
end_extrusion_offset = V(0, 0, start_length)
start_offset = V(0, 0, start_length)
end_extrusion_offset = start_offset.copy()
si_conversion = ifcopenshell.util.unit.calculate_unit_scale(self.file)
# TODO: support offseted profiles
def get_transition_length(start_half_dim, end_half_dim, angle, profile_offset=None):
# NOTE: transition_length == 0 when profiles have the same dimensions
# holy grail of the transition length:
def get_transition_legth(start_half_dim, end_half_dim, angle):
diff = start_half_dim.xy - end_half_dim.xy
diff = Vector([abs(i) for i in diff])
c = diff.x * tan(radians(90 - angle / 2))
a = diff.y
b = (c**2 - a**2) ** 0.5
return b
transition_length = get_transition_legth(start_half_dim, end_half_dim, angle)
def calculate_transition(diff, profile_offset, end_profile=False, angle=None, length=None):
"""will return transition length based on the profile dimension differences and offset.
If `length` is provided will return transition angle"""
# offsets tend to have bunch of float point garbage
# that can result in errors when we're calculating value for square root below
offset = V(0, 0) if profile_offset is None else round_vector_to_precision(profile_offset, si_conversion)
if end_profile:
diff, offset = diff.yx, offset.yx
a = diff.x + offset.x
b = diff.x - offset.x
if length is None:
if diff.x == 0:
return 0
t = tan(radians(angle))
l1 = (a + b + sqrt(a**2 + 4 * a * b * t**2 + 2 * a * b + b**2)) / (2 * t)
length = sqrt(l1**2 - offset.y**2)
# TODO: remove after debug, move somewhere to tests?
if True:
A = (end_profile if end_profile else start_half_dim) * V(1, 0, 0)
end_profile_offset = offset.to_3d() + V(0, 0, length)
D = (start_half_dim if end_profile else end_half_dim) * V(1, 0, 0)
B, C = -A, -D
C += end_profile_offset
D += end_profile_offset
tested_angle = degrees((A - D).angle(B - C))
print(f"length = {length}, requested angle = {angle}, tested angle = {tested_angle}")
return length
elif angle is None:
# TODO: need to handle angle differently for that case
# it occurs when diff == 0
if length == 0:
return 0
l1 = sqrt(length**2 + offset.y**2)
t = -l1 * (a + b) / (a * b - l1**2)
angle = atan(t)
return angle
transition_lengths = [
calculate_transition(diff, profile_offset, angle=angle),
calculate_transition(diff, profile_offset, angle=angle, end_profile=True),
]
other_side_angles = [
calculate_transition(diff, profile_offset, length=transition_lengths[0]),
calculate_transition(diff, profile_offset, length=transition_lengths[1], end_profile=True),
]
# NOTE: debug values
print(f"offset = {profile_offset}")
print(f"diff = {diff}")
print(f"lengths = {transition_lengths}")
print(f"other angles = {other_side_angles}")
print(f"measurable angles = {[(180 - deg)/2 for deg in other_side_angles]}")
# need to make sure that the worst angle (maximum angle)
# for this transition angle is `angle`
for transition_length, other_side_angle in zip(transition_lengths, other_side_angles):
if other_side_angle < angle or is_x(other_side_angle, angle):
print(f"final length = {transition_length}") # TODO: remove after debug
return transition_length
transition_length = get_transition_length(start_half_dim, end_half_dim, angle, profile_offset)
if transition_length is None:
return None, None
faces = []
if transition_length != 0:
end_extrusion_offset.z += transition_length
end_extrusion_offset.z += transition_length
if profile_offset:
end_extrusion_offset.xy += profile_offset
if start_profile.is_a("IfcRectangleProfileDef") and end_profile.is_a("IfcRectangleProfileDef"):
# no transitions for exactly the same profiles
if transition_length == 0:
return None, None
faces += [(3, 4, 7, 0), (11, 8, 15, 12), (3, 11, 12, 4), (7, 15, 8, 0)]
# NOTE: clockwise order for correct face orientation
faces += [
# start extrusion
(0, 1, 2, 3),
(8, 11, 10, 9),
(0, 8, 9, 1),
(1, 9, 10, 2),
(2, 10, 11, 3),
# end extrusion
(4, 5, 6, 7),
(12, 15, 14, 13),
(4, 12, 13, 5),
(5, 13, 14, 6),
(6, 14, 15, 7),
]
points = [
start_half_dim * V(-1, -1, 1),
start_half_dim * V(-1, -1, 0),
start_half_dim * V(1, -1, 0),
start_half_dim * V(1, -1, 1),
end_half_dim * V(1, -1, 0) + end_extrusion_offset,
end_half_dim * V(1, -1, 1) + end_extrusion_offset,
end_half_dim * V(-1, -1, 1) + end_extrusion_offset,
end_half_dim * V(-1, -1, 0) + end_extrusion_offset,
start_half_dim * V(-1, 1, 1),
start_half_dim * V(-1, 1, 0),
start_half_dim * V(1, 1, 0),
start_half_dim * V(1, 1, 1),
end_half_dim * V(1, 1, 0) + end_extrusion_offset,
end_half_dim * V(1, 1, 1) + end_extrusion_offset,
end_half_dim * V(-1, 1, 1) + end_extrusion_offset,
end_half_dim * V(-1, 1, 0) + end_extrusion_offset,
]
# NOTE: clockwise order for correct face orientation
faces += [
# start extrusion
(0, 1, 2, 3),
(8, 11, 10, 9),
(0, 8, 9, 1),
(1, 9, 10, 2),
(2, 10, 11, 3),
# end extrusion
(4, 5, 6, 7),
(12, 15, 14, 13),
(4, 12, 13, 5),
(5, 13, 14, 6),
(6, 14, 15, 7),
]
points = [
start_half_dim * V(-1, -1, 1),
start_half_dim * V(-1, -1, 0),
start_half_dim * V(1, -1, 0),
start_half_dim * V(1, -1, 1),
end_half_dim * V(1, -1, 0) + end_extrusion_offset,
end_half_dim * V(1, -1, 1) + end_extrusion_offset,
end_half_dim * V(-1, -1, 1) + end_extrusion_offset,
end_half_dim * V(-1, -1, 0) + end_extrusion_offset,
start_half_dim * V(-1, 1, 1),
start_half_dim * V(-1, 1, 0),
start_half_dim * V(1, 1, 0),
start_half_dim * V(1, 1, 1),
end_half_dim * V(1, 1, 0) + end_extrusion_offset,
end_half_dim * V(1, 1, 1) + end_extrusion_offset,
end_half_dim * V(-1, 1, 1) + end_extrusion_offset,
end_half_dim * V(-1, 1, 0) + end_extrusion_offset,
]
elif start_profile.is_a("IfcCircleProfileDef") and end_profile.is_a("IfcCircleProfileDef"):
# no transitions for exactly the same profiles
if transition_length == 0:
return None, None
n_segments = 16
first_profile_points = get_circle_points(start_profile.Radius, n_segments)
second_profile_points = get_circle_points(end_profile.Radius, n_segments)
faces = []
for i in range(n_segments):
# For wrapping around the circle
next_i = (i + 1) % n_segments
face = [i, next_i, next_i + n_segments, i + n_segments]
faces.append(face)
transition_items.append(get_profile_faceset(first_profile_points, start_length))
transition_items.append(get_profile_faceset(second_profile_points, end_length, end_extrusion_offset))
first_profile_points = [p + start_offset for p in first_profile_points]
second_profile_points = [p + end_extrusion_offset for p in second_profile_points]
points = first_profile_points + second_profile_points
else: # one is circular, another one is rectangular
# support transition from rectangle to circle of the same dimensions
if transition_length == 0:
transition_length = (start_length + end_length) / 2
end_extrusion_offset.z += transition_length
starting_with_circle = start_profile.is_a("IfcCircleProfileDef")
if starting_with_circle:
circle_profile, rect_profile = start_profile, end_profile
else:
circle_profile, rect_profile = end_profile, start_profile
circle_points = get_circle_points(circle_profile.Radius)
rect_points = get_rectangle_points(V(rect_profile.XDim, rect_profile.YDim, 0))
if starting_with_circle:
start_points, end_points = circle_points, rect_points
else:
start_points, end_points = rect_points, circle_points
transition_items.append(get_profile_faceset(start_points, start_length))
transition_items.append(get_profile_faceset(end_points, end_length, end_extrusion_offset))
# offset verts
if starting_with_circle:
circle_points = [p + start_offset for p in circle_points]
rect_points = [p + end_extrusion_offset for p in rect_points]
else:
rect_points = [p + start_offset for p in rect_points]
circle_points = [p + end_extrusion_offset for p in circle_points]
# circle verts are 0-15, rect verts are 16-19
points = circle_points + rect_points
transition_faces = [
(0, 19, 16), # base
(0, 16, 1),
(1, 16, 2),
(2, 16, 3),
(3, 16, 4),
(4, 16, 17), # base
(4, 17, 5),
(5, 17, 6),
(6, 17, 7),
(7, 17, 8),
(8, 17, 18), # base
(8, 18, 9),
(9, 18, 10),
(10, 18, 11),
(11, 18, 12),
(12, 18, 19), # base
(12, 19, 13),
(13, 19, 14),
(14, 19, 15),
(15, 19, 0),
]
# revert them in case it's starting with circle profile to keep the face orientation
if starting_with_circle:
transition_faces = [f[::-1] for f in transition_faces]
faces += transition_faces
face_set = self.polygonal_face_set(points, faces)
transition_items.append(face_set)