""" CityGml module parses citygml files and import the geometry into the city model structure SPDX - License - Identifier: LGPL - 3.0 - or -later Copyright © 2020 Project Author Guille Gutierrez guillermo.gutierrezmorote@concordia.ca """ import numpy as np import xmltodict import time from city_model_structure.city import City from city_model_structure.building import Building from city_model_structure.attributes.surface import Surface from helpers.geometry_helper import GeometryHelper from city_model_structure.attributes.polygon import Polygon class CityGml: """ CityGml class """ def __init__(self, path): self._city = None with open(path) as gml: # Clean the namespaces is an important task to prevent wrong ns:field due poor citygml implementations self._gml = xmltodict.parse(gml.read(), process_namespaces=True, xml_attribs=True, namespaces={ 'http://www.opengis.net/gml': None, 'http://www.opengis.net/citygml/1.0': None, 'http://www.opengis.net/citygml/building/1.0': None, 'http://schemas.opengis.net/citygml/building/1.0/building.xsd': None, 'http://www.opengis.net/citygml/appearance/1.0': None, 'http://schemas.opengis.net/citygml/appearance/1.0/appearance.xsd': None, 'http://www.opengis.net/citygml/relief/1.0': None, 'http://schemas.opengis.net/citygml/relief/1.0/relief.xsd': None, 'http://www.opengis.net/citygml/generics/1.0': None, 'http://www.w3.org/2001/XMLSchema-instance': None, 'urn:oasis:names:tc:ciq:xsdschema:xAL:2.0': None, 'http://www.w3.org/1999/xlink': None, 'http://www.opengis.net/citygml/relief/2.0': None, 'http://www.opengis.net/citygml/building/2.0': None, 'http://www.opengis.net/citygml/building/2.0 http://schemas.opengis.net/citygml/building/2.0/building.xsd ' 'http://www.opengis.net/citygml/relief/2.0 http://schemas.opengis.net/citygml/relief/2.0/relief.xsd" ' 'xmlns="http://www.opengis.net/citygml/2.0': None, 'http://www.opengis.net/citygml/2.0': None }, force_list=('cityObjectMember', 'curveMember', 'boundedBy', 'surfaceMember')) self._city_objects = None self._geometry = GeometryHelper() for bound in self._gml['CityModel']['boundedBy']: envelope = bound['Envelope'] if '#text' in envelope['lowerCorner']: self._lower_corner = np.fromstring(envelope['lowerCorner']['#text'], dtype=float, sep=' ') self._upper_corner = np.fromstring(envelope['upperCorner']['#text'], dtype=float, sep=' ') else: self._lower_corner = np.fromstring(envelope['lowerCorner'], dtype=float, sep=' ') self._upper_corner = np.fromstring(envelope['upperCorner'], dtype=float, sep=' ') if '@srsName' in envelope: self._srs_name = envelope['@srsName'] @property def content(self): """ CityGml raw content :return: str """ return self._gml @property def city(self) -> City: """ City model structure enriched with the geometry information :return: City """ init = time.process_time_ns() if self._city is None: # todo: refactor this method to clearly choose the gml type self._city = City(self._lower_corner, self._upper_corner, self._srs_name) i = 0 building_part = None for o in self._gml['CityModel']['cityObjectMember']: i += 1 lod = 0 surfaces = [] if 'lod1Solid' in o['Building']: lod += 1 surfaces = CityGml._lod1_solid(o) elif 'lod1MultiSurface' in o['Building']: lod += 1 surfaces = CityGml._lod1_multi_surface(o) elif 'lod2MultiSurface' in o['Building']: # todo: check if this is a real case or a miss-formed citygml lod = 2 surfaces = surfaces + CityGml._lod2_solid_multi_surface(o) else: for bound in o['Building']['boundedBy']: surface_type = next(iter(bound)) if 'lod2MultiSurface' in bound[surface_type]: lod = 2 surfaces = surfaces + CityGml._lod2(bound) if 'lod3Solid' in o['Building']: lod += 4 if 'lod4Solid' in o['Building']: lod += 8 lod_terrain_str = 'lod' + str(lod) + 'TerrainIntersection' terrains = [] if lod_terrain_str in o['Building']: terrains = self._terrains(o, lod_terrain_str) function = None year_of_construction = None if 'consistsOfBuildingPart' in o['Building']: if 'BuildingPart' in o['Building']['consistsOfBuildingPart']: name = o['Building']['consistsOfBuildingPart']['BuildingPart']['name'] if 'yearOfConstruction' in o['Building']['consistsOfBuildingPart']['BuildingPart']: year_of_construction = o['Building']['consistsOfBuildingPart']['BuildingPart']['yearOfConstruction'] if 'function' in o['Building']['consistsOfBuildingPart']['BuildingPart']: function = o['Building']['consistsOfBuildingPart']['BuildingPart']['function'] else: name = o['Building']['@id'] if 'yearOfConstruction' in o['Building']: year_of_construction = o['Building']['yearOfConstruction'] if 'function' in o['Building']: function = o['Building']['function'] self._city.add_city_object(Building(name, lod, surfaces, year_of_construction, function, self._lower_corner, terrains)) return self._city def _terrains(self, city_object, lod_terrain_str): try: curves = [c['LineString']['posList']['#text'] for c in city_object['Building'][lod_terrain_str]['MultiCurve']['curveMember']] except TypeError: curves = [c['LineString']['posList'] for c in city_object['Building'][lod_terrain_str]['MultiCurve']['curveMember']] terrains = [] for curve in curves: curve_points = np.fromstring(curve, dtype=float, sep=' ') curve_points = self._geometry.to_points_matrix(curve_points) terrains.append(curve_points) return terrains @staticmethod def _lod1_solid(o): try: solid_points = [CityGml._solid_points(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList']['#text'])) for s in o['Building']['lod1Solid']['Solid']['exterior']['CompositeSurface']['surfaceMember']] except TypeError: solid_points = [CityGml._solid_points(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'])) for s in o['Building']['lod1Solid']['Solid']['exterior']['CompositeSurface']['surfaceMember']] return [Surface(Polygon(sp),Polygon(sp)) for sp in solid_points] @staticmethod def _lod1_multi_surface(o): solid_points = [CityGml._solid_points(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'])) for s in o['Building']['lod1MultiSurface']['MultiSurface']['surfaceMember']] return [Surface(Polygon(sp),Polygon(sp)) for sp in solid_points] @staticmethod def _lod2_solid_multi_surface(o): if 'boundedBy' in o['Building']['consistsOfBuildingPart']['BuildingPart']: if 'RoofSurface' in o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']: if o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']['RoofSurface']['lod2MultiSurface'] != 'None': polygons = [Polygon(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'])) for s in o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']['RoofSurface']['lod2MultiSurface']['MultiSurface']['surfaceMember']] elif 'WallSurface' in o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']: if o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']['WallSurface']['lod2MultiSurface'] != 'None': polygons = [Polygon(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'])) for s in o['Building']['consistsOfBuildingPart']['BuildingPart']['boundedBy']['WallSurface']['lod2MultiSurface']['MultiSurface']['surfaceMember']] else: polygons = [Polygon(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'])) for s in o['Building']['lod2MultiSurface']['MultiSurface']['surfaceMember']] return [Surface(p,p) for p in polygons] @staticmethod def _lod2_composite_surface(s): solid_points = [CityGml._solid_points((CityGml._remove_last_point(sm['Polygon']['exterior']['LinearRing']['posList']))) for sm in s['CompositeSurface']['surfaceMember']] return [Surface(Polygon(sp),Polygon(sp)) for sp in solid_points] @staticmethod def _lod2_multi_surface(s, surface_type): # todo: this need to be changed into surface bounded? try: solid_points = [CityGml._solid_points(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing']['posList'] ['#text']))] except TypeError: solid_points = [CityGml._solid_points(CityGml._remove_last_point(s['Polygon']['exterior']['LinearRing'] ['posList']))] return [Surface(Polygon(sp),Polygon(sp), surface_type=surface_type) for sp in solid_points] @staticmethod def _lod2(bound): surfaces = [] for surface_type in iter(bound): for s in bound[surface_type]['lod2MultiSurface']['MultiSurface']['surfaceMember']: if 'CompositeSurface' in s: surfaces = surfaces + CityGml._lod2_composite_surface(s) else: surfaces = surfaces + CityGml._lod2_multi_surface(s, surface_type) return surfaces @staticmethod def _remove_last_point(points): array = points.split(' ') res = " " return res.join(array[0:len(array) - 3]) @staticmethod def _solid_points(coordinates) -> np.ndarray: """ Solid surface point matrix [[x, y, z],[x, y, z],...] :parameter coordinates: string from file :return: np.ndarray """ solid_points = np.fromstring(coordinates, dtype=float, sep=' ') solid_points = GeometryHelper.to_points_matrix(solid_points) return solid_points @staticmethod def _holes_points(holes_coordinates) -> [np.ndarray]: """ Holes surfaces point matrices [[[x, y, z],[x, y, z],...]] :parameter holes_coordinates: strings from file :return: [np.ndarray] """ holes_points = [] for hole_coordinates in holes_coordinates: hole_points = np.fromstring(hole_coordinates, dtype=float, sep=' ') hole_points = GeometryHelper.to_points_matrix(hole_points) holes_points.append(hole_points) return holes_points @staticmethod def _perimeter_points(coordinates, solid_points, holes_points) -> np.ndarray: """ Matrix of points of the perimeter in the same order as in coordinates [[x, y, z],[x, y, z],...] :parameter coordinates: string from file :parameter solid_points: np.ndarray points that define the solid part of a surface :parameter holes_points: [np.ndarray] points that define each the holes in a surface :return: np.ndarray """ if holes_points is None: perimeter_points = solid_points else: _perimeter_coordinates = coordinates for hole_points in holes_points: _hole = np.append(hole_points, hole_points[0]) _closed_hole = ' '.join(str(e) for e in [*_hole[:]]) # add a mark 'M' to ensure that the recombination of points does not provoke errors in finding holes _perimeter_coordinates = _perimeter_coordinates.replace(_closed_hole, 'M') _perimeter_coordinates = _perimeter_coordinates.replace('M', '') perimeter_points = np.fromstring(_perimeter_coordinates, dtype=float, sep=' ') perimeter_points = GeometryHelper.to_points_matrix(perimeter_points) # remove duplicated points pv = np.array([perimeter_points[0]]) for point in perimeter_points: duplicated_point = False for p in pv: if GeometryHelper().almost_equal(0.0, p, point): duplicated_point = True if not duplicated_point: pv = np.append(pv, [point], axis=0) perimeter_points = pv return perimeter_points @staticmethod def _holes_polygons(holes_points) -> [Polygon]: """ hole surfaces, a list of hole polygons found in a surface :parameter holes_points: [np.ndarray] that define each of the holes :return: None, [] or [Polygon] None -> not known whether holes exist in reality or not due to low level of detail of input data [] -> no holes in the surface [Polygon] -> one or more holes in the surface """ if holes_points is None: holes_polygons = None else: holes_polygons = [] for hole_points in holes_points: holes_polygons.append(Polygon(hole_points)) return holes_polygons