upadte: approach 3: split into horizontal, then vertical

server/map/get_relation.py

 import geopandas as gpd
-from shapely.geometry import Polygon, Point, LineString
+from shapely.geometry import Polygon, Point, LineString, MultiPolygon
 import matplotlib.pyplot as plt
 from shapely.ops import linemerge, unary_union, polygonize
 import matplotlib.ticker as ticker
@@ -21,20 +21,41 @@ def get_relation(self, body_of_water: str):
         print("Failed to convert to polygons")
         return
 
-    print("Performing div call")
-    test_line = LineString([(10.4600, 60.7451), (11.2000, 60.7451)])
+    print("Creating grid and cutting polygons")
+    divided_map = []
 
-    new_polygons = []
+    # Divide all polygons into sections
     for polygon in polygons:
-        new_polygon = cut_polygon_by_line(polygon, test_line)
-        new_polygons.extend(new_polygon)
+        # Initialise grid
+        grid_lines = create_grid(polygon, cell_size=0.1)
 
-    tiles = gpd.GeoDataFrame(geometry=new_polygons)
+        hrz_lines = grid_lines[0]  # Horizontal lines
+        vrt_lines = grid_lines[1]  # Vertical lines
 
+        # Cut polygon into horizontal sections
+        for hrz_line in hrz_lines:
+            # Splits shape into upper and lower section by hrz_line
+            cut_poly_1 = cut_polygon_in_two(polygon, hrz_line)
+
+            polygon_part = cut_poly_1[0]  # Save upper horizontal section
+
+            # Cut each horizontal section into vertical sections
+            for vrt_line in vrt_lines:
+                cut_poly_2 = cut_polygon_in_two(polygon_part, vrt_line)
+                divided_map.extend(cut_poly_2[0])  # Append part to final list of shapes
+
+                # Set polygon_part to the remaining, un-split, horizontal section for next iteration
+                polygon_part = cut_poly_2[1]
+
+            polygon = cut_poly_1[1]  # Set polygon to the remaining, un-split shape for next iteration
+
+    tiles = gpd.GeoDataFrame(geometry=divided_map)
+
+    # NB: test plots
     fig, ax = plt.subplots()
     ax.set_aspect(1.5)
     ax.xaxis.set_major_locator(ticker.MultipleLocator(0.2))
-    tiles.plot(ax=ax, color='blue', edgecolor='blue', alpha=0.5, label='Original Polygon')
+    tiles.plot(ax=ax, color='blue', edgecolor='orange', alpha=0.5, label='Original Polygon')
     plt.show()
 
     # Convert GeoDataFrame to GeoJSON
@@ -49,7 +70,30 @@ def get_relation(self, body_of_water: str):
     self.wfile.write(tiles_json.encode('utf-8'))
 
 
+def cut_polygon_in_two(polygon, line):
+    points = list(polygon.exterior.coords)  # Extract polygon coordinates
+
+    shape_1 = []
+    shape_2 = []
+
+    # Split points into separate shapes based on their location relative to line
+    for point in points:
+        point = Point(point)  # Convert coordinates to Shapely Point object
+        if line.distance(point) < 1e-10:  # Floating point error tolerance
+            if line.contains(point):  # Determine if point goes into shape 1 or shape 2
+                shape_1.append(point)
+            else:
+                shape_2.append(point)
+
+    # Convert the lists of points to Polygon objects
+    poly_1 = Polygon(shape_1)
+    poly_2 = Polygon(shape_2)
+
+    return poly_1, poly_2
+
+
 # From https://stackoverflow.com/questions/39338550/cut-a-polygon-with-two-lines-in-shapely
+# Splits a polygon into two parts based on a dividing line
 def cut_polygon_by_line(polygon, line):
     if polygon.geom_type == 'Polygon':
         polygon = [polygon]
@@ -58,3 +102,21 @@ def cut_polygon_by_line(polygon, line):
     polygons = polygonize(borders)
     return list(polygons)
 
+
+# Create a grid of lines spaced by a cell size
+def create_grid(polygon, cell_size):
+    min_x, min_y, max_x, max_y = polygon.bounds
+    x_coords = np.arange(min_x, max_x, cell_size)
+    y_coords = np.arange(min_y, max_y, cell_size)
+    horizontal_lines = []
+    vertical_lines = []
+
+    # Vertical lines from top to bottom
+    for x in x_coords:
+        vertical_lines.append(LineString([(x, min_y), (x, max_y)]))
+
+    # Horizontal lines from top to bottom
+    for y in y_coords:
+        horizontal_lines.append(LineString([(min_x, y), (max_x, y)]))
+
+    return horizontal_lines, vertical_lines