city_retrofit/imports/geometry/gpandas.py

"""
gpandas module parses geopandas input table and import the geometry into the city model structure
SPDX - License - Identifier: LGPL - 3.0 - or -later
Copyright © 2022 Concordia CERC group
Project Coder:  Milad Aghamohamadnia --- milad.aghamohamadnia@concordia.ca
"""

import trimesh.exchange.load
from trimesh import Scene
import trimesh.geometry
from shapely.geometry import Polygon as ShapelyPoly
from shapely.geometry import Point
import trimesh
from pyproj import CRS
from city_model_structure.city import City
from city_model_structure.building import Building
from city_model_structure.building_demand.surface import Surface
from city_model_structure.attributes.polygon import Polygon
import numpy as np


class GPandas:
    """
    GeoPandas class
    """
    def __init__(self, gdataframe, srs_name='EPSG:26911'):
        """_summary_

        Arguments:
            gdataframe {Geopandas.Dataframe} -- input geometry data in geopandas table

        Keyword Arguments:
            srs_name {str} -- coordinate system of coordinate system (default: {'EPSG:26911'})
        """
        self._srs_name = srs_name
        self._city = None
        self._scene = gdataframe
        self._scene = self._scene.to_crs(self._srs_name)
        minx, miny, maxx, maxy = self._scene.total_bounds 
        self._lower_corner = [minx, miny, 0]
        self._upper_corner = [maxx, maxy, 0]
        
    @property
    def scene(self) -> Scene:
        """
        Get GeoPandas scene
        """
        return self._scene
        
    @property
    def city(self) -> City:
        """
        Get city out of a GeoPandas Table
        """
        if self._city is None:            
            self._city = City(self._lower_corner, self._upper_corner, self._srs_name)
            for ix, bldg in self._scene.iterrows():
                polygon = bldg['geometry']
                geom = bldg.geom
                polygon = ShapelyPoly(geom['coordinates'][0])
                height = float(bldg['height_mean'])
                building_mesh = trimesh.creation.extrude_polygon(polygon, height)
                trimesh.repair.fill_holes(building_mesh)
                trimesh.repair.fix_winding(building_mesh)
                year_of_construction = int(bldg['year_built'])
                name = str(ix)
                lod = 1
                if year_of_construction>2000:
                    function = 'residential'
                else:
                    function = 'industry'
                average_storey_height = 3
                
                surfaces = []
                face_normal = building_mesh.face_normals
                for ix, face in enumerate(building_mesh.faces):
                    points = []
                    for vertex_index in face:
                        points.append(building_mesh.vertices[vertex_index])
                    solid_polygon = Polygon(points)
                    perimeter_polygon = solid_polygon
                    s_type = 'Ground' if face_normal[ix][2]==-1 else ('Roof' if face_normal[ix][2]==1 else 'Wall')
                    surface = Surface(solid_polygon, perimeter_polygon, surface_type=s_type)
                    surfaces.append(surface)
                building = Building(name, lod, surfaces, year_of_construction, function, self._lower_corner, terrains=None)
                self._city.add_city_object(building)
        return self._city
    
    @staticmethod
    def resize_polygon(poly, factor=0.10, expand=False):
        """
            returns the shapely polygon which is smaller or bigger by passed factor.
            
        Arguments:
            poly {shapely.geometry.Polygon} -- an input geometry in shapely polygon format

        Keyword Arguments:
            factor {float} -- factor of expansion (default: {0.10})
            expand {bool} -- If expand = True , then it returns bigger polygon, else smaller (default: {False})

        Returns:
            {shapely.geometry.Polygon} -- output geometry in shapely polygon format
        """
        xs = list(poly.exterior.coords.xy[0])
        ys = list(poly.exterior.coords.xy[1])
        x_center = 0.5 * min(xs) + 0.5 * max(xs)
        y_center = 0.5 * min(ys) + 0.5 * max(ys)
        min_corner = Point(min(xs), min(ys))
        max_corner = Point(max(xs), max(ys))
        center = Point(x_center, y_center)
        shrink_distance = center.distance(min_corner)*factor

        if expand:
            poly_resized = poly.buffer(shrink_distance) #expand
        else:
            poly_resized = poly.buffer(-shrink_distance) #shrink
        return poly_resized