city_retrofit/city_model_structure/surface.py

from __future__ import annotations
import numpy as np
import pyny3d.geoms as pn
from helpers.geometry import Geometry


class Surface:
  def __init__(self, coordinates, surface_type=None, name=None, swr='0.2', remove_last=True, is_projected=False):
    self._coordinates = coordinates
    self._type = surface_type
    self._name = name
    self._swr = swr
    self._remove_last = remove_last
    self._is_projected = is_projected
    self._geometry = Geometry()
    self._polygon = None
    self._area = None
    self._points = None
    self._points_list = None
    self._normal = None
    self._azimuth = None
    self._inclination = None
    self._area_above_ground = None
    self._area_below_ground = None
    self._parent = None
    self._shapely = None
    self._projected_surface = None
    self._shared_surfaces = None
    self._global_irradiance_hour = np.zeros(8760)
    self._global_irradiance_month = np.zeros(12)

  def parent(self, parent, surface_id):
    self._parent = parent
    self._name = str(surface_id)

  @property
  def name(self):
    if self._name is None:
      raise Exception('surface has no name')
    return self._name

  @property
  def swr(self):
    return self._swr

  @swr.setter
  def swr(self, value):
    self._swr = value

  @property
  def points(self):
    if self._points is None:
      self._points = np.fromstring(self._coordinates, dtype=float, sep=' ')
      self._points = Geometry.to_points_matrix(self._points, self._remove_last)
    return self._points

  @property
  def points_list(self):
    if self._points_list is None:
      s = self.points
      self._points_list = np.reshape(s, len(s) * 3)
    return self._points_list

  @property
  def polygon(self):
    if self._polygon is None:
      try:
        self._polygon = pn.Polygon(self.points)
      except ValueError:
        # is not really a polygon but a line so just return none
        self._polygon = None
    return self._polygon

  @property
  def area(self):
    if self._area is None:
      self._area = self.polygon.get_area()
    return self._area

  def _is_almost_same_terrain(self, terrain_points, ground_points):
    equal = 0
    for t in terrain_points:
      for g in ground_points:
        if self._geometry.almost_equal(t, g):
          equal += 1
    return equal == len(terrain_points)

  @property
  def _is_terrain(self):
    for t_points in self._parent.terrains:
      if len(t_points) == len(self.points) and self._is_almost_same_terrain(t_points, self.points):
        return True
    return False

  @property
  def area_above_ground(self):
    if self._area_above_ground is None:
      self._area_above_ground = self.area - self.area_below_ground
    return self._area_above_ground

  @property
  def area_below_ground(self):
    if self._area_below_ground is None:
      self._area_below_ground = 0.0
      if self._is_terrain:
        self._area_below_ground = self.area
    return self._area_below_ground

  @property
  def normal(self):
    if self._normal is None:
      points = self.points
      n = np.cross(points[1] - points[0], points[2] - points[0])
      self._normal = n / np.linalg.norm(n)
    return self._normal

  @property
  def azimuth(self):
    if self._azimuth is None:
      normal = self.normal
      self._azimuth = np.arctan2(normal[1], normal[0])
    return self._azimuth

  @property
  def inclination(self):
    if self._inclination is None:
      self._inclination = np.arccos(self.normal[2])
    return self._inclination

  @property
  def type(self):
    if self._type is None:
      grad = np.rad2deg(self.inclination)
      if 170 <= grad:
        self._type = 'Ground'
      elif 80 <= grad <= 100:
        self._type = 'Wall'
      else:
        self._type = 'Roof'
    return self._type

  def shared(self, surface):
    if self.type is not 'Wall':
      return
    if self._geometry.is_almost_same_surface(self, surface):
      self._shared_surfaces.append((100, surface))
      surface.shared(self)

  @property
  def global_irradiance_hour(self):
    return self._global_irradiance_hour

  @global_irradiance_hour.setter
  def global_irradiance_hour(self, value):
    self._global_irradiance_hour = value

  @property
  def global_irradiance_month(self):
    return self._global_irradiance_month

  @global_irradiance_month.setter
  def global_irradiance_month(self, value):
    self._global_irradiance_month = value

  @property
  def shapely(self):
    if self.polygon is None:
      return None
    if self._shapely is None:
      self._shapely = self.polygon.get_shapely()
    return self._shapely

  @property
  def projection(self) -> Surface:
    if self._is_projected:
      return self
    if self._projected_surface is None:
      coordinates = ''
      for coordinate in self.shapely.exterior.coords:
        if coordinates != '':
          coordinates = coordinates + ' '
        coordinates = coordinates + str(coordinate[0]) + ' ' + str(coordinate[1]) + ' 0.0'
      self._projected_surface = Surface(coordinates)
    return self._projected_surface