From a41e7fc29890ca4fd28e7e91cb5864ac55e741ff Mon Sep 17 00:00:00 2001 From: Dion Moult Date: Thu, 26 Sep 2019 07:37:08 +1000 Subject: [PATCH] Add SVG sectioning utility --- src/ifcblenderexport/generate_section.py | 447 +++++++++++++++++++++++ 1 file changed, 447 insertions(+) create mode 100644 src/ifcblenderexport/generate_section.py diff --git a/src/ifcblenderexport/generate_section.py b/src/ifcblenderexport/generate_section.py new file mode 100644 index 0000000000..57785d9873 --- /dev/null +++ b/src/ifcblenderexport/generate_section.py @@ -0,0 +1,447 @@ +import os +import math +import svgwrite +import time + +import OCC.gp +import OCC.Geom +import OCC.Bnd +import OCC.BRepBndLib +import OCC.BRep +import OCC.BRepPrimAPI +import OCC.BRepAlgoAPI +import OCC.BRepBuilderAPI +import OCC.BRepAlgo +import OCC.TopOpeBRepTool +import OCC.TopOpeBRepBuild +import OCC.ShapeExtend +import OCC.GProp +import OCC.BRepGProp +import OCC.GC +import OCC.ShapeAnalysis +import OCC.TopTools +import OCC.TopExp +import OCC.HLRAlgo +import OCC.HLRBRep +import OCC.TopLoc +import OCC.Bnd +import OCC.BRepBndLib +import OCC.BRepTools +import OCC.TopoDS +import OCC.GeomLProp + +from OCC.TopoDS import topods + +import ifcopenshell +import ifcopenshell.geom + +class IfcCutter: + def __init__(self): + self.product_shapes = [] + self.background_elements = [] + self.cut_polygons = [] + self.ifc_file = None + self.unit = None + self.section_box = { + 'top_left_corner': (-3., -2.45, 4.75), + 'projection': (0., 1., 0.), + 'x_axis': (1., 0., 0.), + 'x': 15., + 'y': 7., + 'z': 3.5, + 'shape': None, + 'face': None + } + + def cut(self): + self.ifc_file = ifcopenshell.open('output.ifc') + self.get_units() + self.get_product_shapes() + self.create_section_box() + self.get_cut_polygons() + self.get_background_elements() + self.sort_background_elements() + + def get_units(self): + unit_assignment = self.ifc_file.by_type('IfcUnitAssignment')[0] + for unit in unit_assignment.Units: + if unit.UnitType == 'LENGTHUNIT': + self.unit = unit + break + + def get_product_shapes(self): + settings = ifcopenshell.geom.settings() + settings.set(settings.USE_PYTHON_OPENCASCADE, True) + products = self.ifc_file.by_type('IfcProduct') + for product in products: + if product.is_a('IfcOpeningElement') or product.is_a('IfcSite'): + continue + if product.Representation is not None: + shape = ifcopenshell.geom.create_shape(settings, product).geometry + self.product_shapes.append((product, shape)) + + def sort_background_elements(self): + new_list = sorted(self.background_elements, key=lambda k: k['z']) + self.background_elements = new_list + + def create_section_box(self): + top_left_corner = OCC.gp.gp_Pnt( + self.section_box['top_left_corner'][0], + self.section_box['top_left_corner'][1], + self.section_box['top_left_corner'][2]) + axis = OCC.gp.gp_Ax2( + top_left_corner, + OCC.gp.gp_Dir( + self.section_box['projection'][0], + self.section_box['projection'][1], + self.section_box['projection'][2]), + OCC.gp.gp_Dir( + self.section_box['x_axis'][0], + self.section_box['x_axis'][1], + self.section_box['x_axis'][2]) + ) + section_box = OCC.BRepPrimAPI.BRepPrimAPI_MakeBox( + axis, self.section_box['x'], self.section_box['y'], self.section_box['z'] + ) + self.section_box['shape'] = section_box.Shape() + self.section_box['face'] = section_box.BottomFace() + + source = OCC.gp.gp_Ax3(axis) + destination = OCC.gp.gp_Ax3( + top_left_corner, + OCC.gp.gp_Dir(0, 0, 1), + OCC.gp.gp_Dir(1, 0, 0)) + self.transformation = OCC.gp.gp_Trsf() + self.transformation.SetTransformation(source, destination) + + transformed_face = OCC.BRepBuilderAPI.BRepBuilderAPI_Transform( + self.section_box['face'], self.transformation).Shape() + exp = OCC.TopExp.TopExp_Explorer(transformed_face, OCC.TopAbs.TopAbs_VERTEX) + transformed_vertex = OCC.BRep.BRep_Tool.Pnt(topods.Vertex(exp.Current())) + + self.transformation.SetTranslationPart( + OCC.gp.gp_Vec(-transformed_vertex.X(), -transformed_vertex.Y(), 0)) + + def get_background_elements(self): + for product, shape in self.product_shapes: + print('Processing product {} '.format(product.Name)) + intersection = OCC.BRepAlgoAPI.BRepAlgoAPI_Common(self.section_box['shape'], shape).Shape() + intersection_edges = self.get_booleaned_edges(intersection) + if len(intersection_edges) <= 0: + continue + + transformed_intersection = OCC.BRepBuilderAPI.BRepBuilderAPI_Transform( + intersection, self.transformation) + intersection = transformed_intersection.Shape() + + edge_face_map = OCC.TopTools.TopTools_IndexedDataMapOfShapeListOfShape() + OCC.TopExp.topexp.MapShapesAndAncestors( + intersection, OCC.TopAbs.TopAbs_EDGE, + OCC.TopAbs.TopAbs_FACE, edge_face_map) + + exp = OCC.TopExp.TopExp_Explorer(intersection, OCC.TopAbs.TopAbs_FACE) + while exp.More(): + face = topods.Face(exp.Current()) + zpos, zmax = self.calculate_face_zpos(face) + self.build_new_face(face, zpos, product) + self.get_split_edges(edge_face_map, face, zmax, product) + exp.Next() + + def calculate_face_zpos(self, face): + bbox = OCC.Bnd.Bnd_Box() + OCC.BRepBndLib.brepbndlib_Add(face, bbox) + xmin, ymin, zmin, xmax, ymax, zmax = bbox.Get() + zpos = zmin + ((zmax - zmin)/2) + return zpos, zmax + + def get_split_edges(self, edge_face_map, face, zmax, product): + exp2 = OCC.TopExp.TopExp_Explorer(face, OCC.TopAbs.TopAbs_EDGE) + while exp2.More(): + edge = topods.Edge(exp2.Current()) + adjface = OCC.TopoDS.TopoDS_Face() + getadj = OCC.TopOpeBRepBuild.TopOpeBRepBuild_Tools.GetAdjacentFace(face, edge, edge_face_map, adjface) + if getadj: + edge_angle = math.degrees(self.get_angle_between_faces(face, adjface)) + if edge_angle > 30 and edge_angle < 160: + newedge = self.build_new_edge(edge, zmax+0.01) + if newedge: + self.background_elements.append({ + 'raw': product, + 'geometry': newedge, + 'type': 'line', + 'z': zmax+0.01 + }) + exp2.Next() + + def get_angle_between_faces(self, f1, f2): + return self.convert_dot_product_to_angle( + self.get_dot_product_of_normals( + self.get_normal(f1), self.get_normal(f2))) + + def get_normal(self, face): + surface = OCC.Geom.Handle_Geom_Surface(OCC.BRep.BRep_Tool.Surface(face)) + props = OCC.GeomLProp.GeomLProp_SLProps(surface, 0, 0, 1, .001) + return props.Normal() + + def get_dot_product_of_normals(self, n1, n2): + return n1.X() * n2.X() + n1.Y() * n2.Y() + n1.Z() * n2.Z() + + def convert_dot_product_to_angle(self, dp): + return math.acos(dp) + + def is_same_point(self, p1, p2): + return p1.X() == p2.X() \ + and p1.Y() == p2.Y() \ + and p1.Z() == p2.Z() + + def build_new_edge(self, edge, zpos): + exp = OCC.TopExp.TopExp_Explorer(edge, OCC.TopAbs.TopAbs_VERTEX) + new_vertices = [] + while exp.More(): + current_vertex = topods.Vertex(exp.Current()) + current_point = OCC.BRep.BRep_Tool.Pnt(current_vertex) + current_point.SetZ(zpos) + new_vertices.append(OCC.BRepBuilderAPI.BRepBuilderAPI_MakeVertex(current_point).Vertex()) + exp.Next() + try: + return OCC.BRepBuilderAPI.BRepBuilderAPI_MakeEdge( + new_vertices[0], new_vertices[1] + ).Edge() + except: + return None + + def build_new_face(self, face, zpos, product): + exp = OCC.TopExp.TopExp_Explorer(face, OCC.TopAbs.TopAbs_WIRE) + while exp.More(): + wireexp = OCC.BRepTools.BRepTools_WireExplorer(topods.Wire(exp.Current())) + new_wire_builder = OCC.BRepBuilderAPI.BRepBuilderAPI_MakeWire() + first_vertex = None + previous_vertex = None + while wireexp.More(): + current_vertex = wireexp.CurrentVertex() + current_point = OCC.BRep.BRep_Tool.Pnt(current_vertex) + # Dodgy technique to squash in Z axis + current_point.SetZ(zpos) + current_vertex = OCC.BRepBuilderAPI.BRepBuilderAPI_MakeVertex(current_point).Vertex() + if not first_vertex: + first_vertex = current_vertex + if not previous_vertex: + previous_vertex = current_vertex + else: + try: + new_wire_builder.Add(topods.Edge( + OCC.BRepBuilderAPI.BRepBuilderAPI_MakeEdge( + previous_vertex, current_vertex + ).Edge())) + previous_vertex = current_vertex + except: + pass + wireexp.Next() + + # make last edge + if not wireexp.More(): + try: + new_wire_builder.Add(topods.Edge( + OCC.BRepBuilderAPI.BRepBuilderAPI_MakeEdge( + current_vertex, first_vertex + ).Edge())) + except: + pass + try: + new_wire = new_wire_builder.Wire() + new_face = OCC.BRepBuilderAPI.BRepBuilderAPI_MakeFace(new_wire).Face() + area = self.get_area(new_face) + if area == 0: + self.background_elements.append({ + 'raw': product, + 'geometry': new_wire, + 'geometry_face': new_face, + 'type': 'polyline', + 'z': zpos + }) + else: + self.background_elements.append({ + 'raw': product, + 'geometry': new_wire, + 'geometry_face': new_face, + 'type': 'polygon', + 'z': zpos + }) + except: + print('Could not build face') + exp.Next() + + def get_area(self, shape): + gprops = OCC.GProp.GProp_GProps() + OCC.BRepGProp.brepgprop.SurfaceProperties(shape, gprops) + return gprops.Mass() + + def get_booleaned_edges(self, shape): + edges = [] + exp = OCC.TopExp.TopExp_Explorer(shape, OCC.TopAbs.TopAbs_EDGE) + while exp.More(): + edges.append(topods.Edge(exp.Current())) + exp.Next() + return edges + + def get_cut_polygons(self): + for product, shape in self.product_shapes: + section = OCC.BRepAlgoAPI.BRepAlgoAPI_Section(self.section_box['face'], shape).Shape() + section_edges = self.get_booleaned_edges(section) + if len(section_edges) <= 0: + continue + wires = self.connect_edges_into_wires(section_edges) + for i in range(wires.Length()): + wire_shape = wires.Value(i+1) + + transformed_wire = OCC.BRepBuilderAPI.BRepBuilderAPI_Transform( + wire_shape, self.transformation) + wire_shape = transformed_wire.Shape() + + wire = topods.Wire(wire_shape) + + self.cut_polygons.append({ + 'raw': product, + 'geometry': wire + }) + + def connect_edges_into_wires(self, unconnected_edges): + edges = OCC.TopTools.TopTools_HSequenceOfShape() + edges_handle = OCC.TopTools.Handle_TopTools_HSequenceOfShape(edges) + wires = OCC.TopTools.TopTools_HSequenceOfShape() + wires_handle = OCC.TopTools.Handle_TopTools_HSequenceOfShape(wires) + + for edge in unconnected_edges: + edges.Append(edge) + + OCC.ShapeAnalysis.ShapeAnalysis_FreeBounds.ConnectEdgesToWires(edges_handle, 1e-5, True, wires_handle) + return wires_handle.GetObject() + +class IfcCutterDebug(IfcCutter): + def cut(self): + self.occ_display = ifcopenshell.geom.utils.initialize_display() + super().cut() + + def create_section_box(self): + super().create_section_box() + self.display_everything_with_section_plane() + + def get_cut_polygons(self): + super().get_cut_polygons() + self.display_cut_polygons() + + def get_background_elements(self): + super().get_background_elements() + self.display_background_elements() + + def display_everything_with_section_plane(self): + section_face_display = ifcopenshell.geom.utils.display_shape(self.section_box['face']) + ifcopenshell.geom.utils.set_shape_transparency(section_face_display, 0.5) + section_box_display = ifcopenshell.geom.utils.display_shape(self.section_box['shape']) + ifcopenshell.geom.utils.set_shape_transparency(section_box_display, 0.2) + + transformed_box = OCC.BRepBuilderAPI.BRepBuilderAPI_Transform( + self.section_box['shape'], self.transformation) + box_display = ifcopenshell.geom.utils.display_shape(transformed_box.Shape()) + ifcopenshell.geom.utils.set_shape_transparency(box_display, 0.2) + + for shape in self.product_shapes: + ifcopenshell.geom.utils.display_shape(shape[1]) + input('Debug: showing everything with section plane.') + + def display_cut_polygons(self): + self.occ_display.EraseAll() + for polygon in self.cut_polygons: + ifcopenshell.geom.utils.display_shape(polygon['geometry'], clr='BLACK') + face = OCC.BRepBuilderAPI.BRepBuilderAPI_MakeFace(polygon['geometry']).Face() + face_display = ifcopenshell.geom.utils.display_shape(face) + ifcopenshell.geom.utils.set_shape_transparency(face_display, 0.5) + input('Debug: showing cut polygons.') + + def display_background_elements(self): + self.occ_display.EraseAll() + for element in self.background_elements: + if element['type'] == 'line': + ifcopenshell.geom.utils.display_shape(element['geometry'], clr='PURPLE') + elif element['type'] == 'polyline': + ifcopenshell.geom.utils.display_shape(element['geometry_face'], clr='RED') + elif element['type'] == 'polygon': + ifcopenshell.geom.utils.display_shape(element['geometry_face']) + input('Debug: showing background elements.') + +class SvgWriter(): + def __init__(self, ifc_cutter): + self.ifc_cutter = ifc_cutter + self.scale = 1 / 100 # 1:100 + + def write(self): + self.calculate_scale() + self.svg = svgwrite.Drawing('test.svg', + size=('{}mm'.format(self.ifc_cutter.section_box['x'] * self.scale), + '{}mm'.format(self.ifc_cutter.section_box['y'] * self.scale)), + viewBox=('0 0 {} {}'.format( + self.ifc_cutter.section_box['x'] * self.scale, + self.ifc_cutter.section_box['y'] * self.scale))) + self.draw_background_elements() + self.draw_cut_polygons() + self.svg.save() + + def calculate_scale(self): + # TODO: properly handle units + if self.ifc_cutter.unit.Name == 'METRE': + self.scale *= 1000 + + def draw_background_elements(self): + for element in self.ifc_cutter.background_elements: + if element['type'] == 'polygon': + self.draw_polygon(element, 'white', 0) + elif element['type'] == 'polyline': + self.draw_polyline(element) + elif element['type'] == 'line': + self.draw_line(element) + + def draw_cut_polygons(self): + for polygon in self.ifc_cutter.cut_polygons: + self.draw_polygon(polygon, 'red', 0.5) + + def draw_polyline(self, element): + exp = OCC.BRepTools.BRepTools_WireExplorer(element['geometry']) + points = [] + while exp.More(): + point = OCC.BRep.BRep_Tool.Pnt(exp.CurrentVertex()) + points.append((point.X() * self.scale, -point.Y() * self.scale)) + exp.Next() + self.svg.add(self.svg.polyline(points=points)) + + def draw_line(self, element): + exp = OCC.TopExp.TopExp_Explorer(element['geometry'], OCC.TopAbs.TopAbs_VERTEX) + points = [] + while exp.More(): + point = OCC.BRep.BRep_Tool.Pnt(topods.Vertex(exp.Current())) + points.append((point.X() * self.scale, -point.Y() * self.scale)) + exp.Next() + self.svg.add(self.svg.line(start=points[0], end=points[1], + stroke=svgwrite.rgb(0, 0, 0, '%'), stroke_width=0.2)) + + def draw_polygon(self, polygon, fill, stroke_width): + exp = OCC.BRepTools.BRepTools_WireExplorer(polygon['geometry']) + points = [] + while exp.More(): + point = OCC.BRep.BRep_Tool.Pnt(exp.CurrentVertex()) + points.append((point.X() * self.scale, -point.Y() * self.scale)) + exp.Next() + self.svg.add(self.svg.polygon( + points=points, fill=fill, + stroke=svgwrite.rgb(0, 0, 0, '%'), + stroke_width=stroke_width, + stroke_linejoin='round' + )) + +print('# Starting process') +start_time = time.time() +ifc_cutter = IfcCutter() +#ifc_cutter = IfcCutterDebug() +svg_writer = SvgWriter(ifc_cutter) +ifc_cutter.cut() +svg_writer.write() +print('# Process finished in {:.2f} seconds'.format(time.time() - start_time))