system_assignation/imports/geometry/citygml.py

"""
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

from city_model_structure.city import City
from city_model_structure.building import Building
from helpers.geometry_helper import GeometryHelper
from imports.geometry.citygml_lod2 import CityGmlLod2
from imports.geometry.citygml_lod1 import CityGmlLod1


class CityGml:
  """
  CityGml class
  """
  def __init__(self, path):
    self._city = None
    self._lod1_tags = ['lod1Solid', 'lod1MultiSurface']
    self._lod2_tags = ['lod2Solid', 'lod2MultiSurface', 'lod2MultiCurve']
    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

  def _create_building(self, city_object):
    name = city_object['@id']
    function = None
    year_of_construction = None
    if 'yearOfConstruction' in city_object:
      year_of_construction = city_object['yearOfConstruction']
    if 'function' in city_object:
      function = city_object['function']
    if any(key in city_object for key in self._lod1_tags):
      lod = 1
      surfaces = CityGmlLod1(city_object).surfaces
    elif any(key in city_object for key in self._lod2_tags):
      lod = 2
      surfaces = CityGmlLod2(city_object).surfaces
    else:
      raise NotImplementedError("Not supported level of detail")
    return Building(name, lod, surfaces, year_of_construction, function, self._lower_corner, [])

  def _create_parts_consisting_building(self):
    raise NotImplementedError

  @property
  def city(self) -> City:
    """
    City model structure enriched with the geometry information
    :return: City
    """

    if self._city is None:
      self._city = City(self._lower_corner, self._upper_corner, self._srs_name)

      for o in self._gml['CityModel']['cityObjectMember']:
        city_object = o['Building']
        if 'consistsOfBuildingPart' in city_object:
          raise NotImplementedError("Building cluster")
        '''
        if 'lod1Solid' in o['Building']:
          lod += 1
          surfaces = CityGmlLod1.lod1_solid(o)
        elif 'lod1MultiSurface' in o['Building']:
          lod += 1
          surfaces = CityGmlLod1.lod1_multi_surface(o)
        
        elif 'lod2Solid' in o['Building'] :
          lod += 1
          surfaces = CityGmlLod2.lod2_solid(o)
        elif 'lod2MultiSurface' in o['Building']:
          # todo: check if this is a real case or a miss-formed citygml
          lod = 2
          surfaces = surfaces + CityGmlLod2.lod2_solid_multi_surface(o)
        else:
          self._city.add_city_object(self._create_building(city_object))
    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 _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