From e03fc6e1b374702d7318b9a181a8c852a9d7d213 Mon Sep 17 00:00:00 2001
From: guille <guillermo.gutierrezmorote@concordia.ca>
Date: Fri, 17 Mar 2023 12:45:53 -0400
Subject: [PATCH] export to obj completed

---
 .../geometry/helpers/geometry_helper.py       | 61 +++++++++++++++++++
 hub/unittests/test_geometry_factory.py        | 21 ++++---
 2 files changed, 72 insertions(+), 10 deletions(-)

diff --git a/hub/imports/geometry/helpers/geometry_helper.py b/hub/imports/geometry/helpers/geometry_helper.py
index 9c72d855..54aa0a5e 100644
--- a/hub/imports/geometry/helpers/geometry_helper.py
+++ b/hub/imports/geometry/helpers/geometry_helper.py
@@ -4,6 +4,8 @@ SPDX - License - Identifier: LGPL - 3.0 - or -later
 Copyright © 2022 Concordia CERC group
 Project Coder Pilar Monsalvete Alvarez de Uribarri pilar.monsalvete@concordia.ca
 """
+import math
+import sys
 
 import numpy as np
 
@@ -52,3 +54,62 @@ class GeometryHelper:
     for point in points:
       res.insert(0,point)
     return res
+
+  @staticmethod
+  def ground_area(points):
+    """
+    Get ground surface area in square meters
+    :return: float
+    """
+    # New method to calculate area
+
+    if len(points) < 3:
+      sys.stderr.write('Warning: the area of a line or point cannot be calculated 1. Area = 0\n')
+      return 0
+    alpha = 0
+    vec_1 = points[1] - points[0]
+    for i in range(2, len(points)):
+      vec_2 = points[i] - points[0]
+      alpha += GeometryHelper.angle_between_vectors(vec_1, vec_2)
+    if alpha == 0:
+      sys.stderr.write('Warning: the area of a line or point cannot be calculated 2. Area = 0\n')
+      return 0
+    #
+    horizontal_points = points
+    area = 0
+    for i in range(0, len(horizontal_points) - 1):
+      point = horizontal_points[i]
+      next_point = horizontal_points[i + 1]
+      area += (next_point[1] + point[1]) / 2 * (next_point[0] - point[0])
+    next_point = horizontal_points[0]
+    point = horizontal_points[len(horizontal_points) - 1]
+    area += (next_point[1] + point[1]) / 2 * (next_point[0] - point[0])
+    _area = abs(area)
+    return _area
+
+  @staticmethod
+  def angle_between_vectors(vec_1, vec_2):
+    """
+    angle between vectors in radians
+    :param vec_1: vector
+    :param vec_2: vector
+    :return: float
+    """
+    if np.linalg.norm(vec_1) == 0 or np.linalg.norm(vec_2) == 0:
+      sys.stderr.write("Warning: impossible to calculate angle between planes' normal. Return 0\n")
+      return 0
+    cosine = np.dot(vec_1, vec_2) / np.linalg.norm(vec_1) / np.linalg.norm(vec_2)
+    if cosine > 1 and cosine - 1 < 1e-5:
+      cosine = 1
+    elif cosine < -1 and cosine + 1 > -1e-5:
+      cosine = -1
+    alpha = math.acos(cosine)
+    return alpha
+
+
+  @staticmethod
+  def invert_points(points):
+    res = []
+    for point in points:
+      res.insert(0,point)
+    return res
diff --git a/hub/unittests/test_geometry_factory.py b/hub/unittests/test_geometry_factory.py
index 35e439dd..2449d674 100644
--- a/hub/unittests/test_geometry_factory.py
+++ b/hub/unittests/test_geometry_factory.py
@@ -4,11 +4,13 @@ SPDX - License - Identifier: LGPL - 3.0 - or -later
 Copyright © 2022 Concordia CERC group
 Project Coder Pilar Monsalvete Alvarez de Uribarri pilar.monsalvete@concordia.ca
 """
-
+import datetime
 from pathlib import Path
 from unittest import TestCase
 from hub.helpers.geometry_helper import GeometryHelper
 
+from numpy import inf
+
 import hub.exports.exports_factory
 from hub.imports.construction_factory import ConstructionFactory
 from hub.imports.geometry_factory import GeometryFactory
@@ -132,13 +134,18 @@ class TestGeometryFactory(TestCase):
     """
     Test geojson import
     """
-    file = 'concordia.geojson'
+    file = '10000_buildings.geojson'
+    start = datetime.datetime.now()
     city = self._get_city(file, 'geojson',
-                          height_field='citygml_me',
+                          height_field='building_height',
                           year_of_construction_field='ANNEE_CONS',
                           function_field='CODE_UTILI')
-
+    end = datetime.datetime.now()
+    print(f'geometry load in {end-start} s')
+    start = datetime.datetime.now()
     hub.exports.exports_factory.ExportsFactory('obj', city, self._output_path).export()
+    end = datetime.datetime.now()
+    print(f'geometry export in {end - start} s')
     self.assertEqual(195, len(city.buildings), 'wrong number of buildings')
     self._check_buildings(city)
 
@@ -147,11 +154,6 @@ class TestGeometryFactory(TestCase):
     Test neighbours map creation
     """
     file = 'neighbours.geojson'
-    city = self._get_city(file, 'geojson',
-                          year_of_construction_field='ANNEE_CONS',
-                          function_field='LIBELLE_UT')
-    info_lod0 = GeometryHelper.city_mapping(city, plot=False)
-
     city = self._get_city(file, 'geojson',
                           height_field='citygml_me',
                           year_of_construction_field='ANNEE_CONS',
@@ -161,7 +163,6 @@ class TestGeometryFactory(TestCase):
     self.assertEqual(info_lod0, info_lod1)
     for building in city.buildings:
       self.assertEqual(2, len(building.neighbours))
-      print(building.volume)
 
     self.assertEqual('2_part_0_zone_0',city.city_object('1_part_0_zone_0').neighbours[0].name)
     self.assertEqual('3_part_0_zone_0',city.city_object('1_part_0_zone_0').neighbours[1].name)