diff --git a/server/map/__pycache__/get_relation.cpython-311.pyc b/server/map/__pycache__/get_relation.cpython-311.pyc
index 96f80953ff49aec41eb8a8c82da8bd1bc35362b4..a7b87b9c6a75148f670c316d70ce731e9b58227a 100644
Binary files a/server/map/__pycache__/get_relation.cpython-311.pyc and b/server/map/__pycache__/get_relation.cpython-311.pyc differ
diff --git a/server/map/get_relation.py b/server/map/get_relation.py
index 40bccf0d40ff030185c86162388bc48d09c2a354..efcf3c2c99fcb4499b17d4c27c54ab06dfe839b2 100644
--- a/server/map/get_relation.py
+++ b/server/map/get_relation.py
@@ -19,7 +19,7 @@ def get_relation(self, body_of_water: str):  # NB: implement body_of_water
     # Divide all polygons into sections
     for polygon in polygons:
 
-        cell_size = 1.5
+        cell_size = 0.4
         # Divide the length and with of polygon into a grid of equally sized parts
         grid_lines = create_grid_coords(polygon, cell_size)
 
@@ -28,6 +28,7 @@ def get_relation(self, body_of_water: str):  # NB: implement body_of_water
 
         # Cut polygon into horizontal sections, bottom to top
         for hrz_line in hrz_lines:
+
             # Split shape into upper and lower section as hrz_line as divider
             divided_poly = cut_polygon_by_points(polygon, hrz_line, cell_size)
             horizontal_section = divided_poly[0]  # Save upper horizontal section
@@ -43,7 +44,7 @@ def get_relation(self, body_of_water: str):  # NB: implement body_of_water
                     break  # Break from loop im remaining section has no coordinates
 
                 # Split the horizontal section into two vertical parts
-                vertical_parts = cut_polygon_by_points(horizontal_section, vrt_line, cell_size)
+                vertical_parts = cut_polygon_by_points(horizontal_section, vrt_line, -0.1)
 
                 divided_map.append(vertical_parts[0])  # Append split vertical sections to final list of shapes
 
@@ -84,8 +85,89 @@ def get_relation(self, body_of_water: str):  # NB: implement body_of_water
     self.wfile.write(tiles_json.encode('utf-8'))
 
 
+'''
 # Takes a polygon and divides its coordinates into two shapes, where divisor is a
 # coordinate that defines the point of division.
+def cut_polygon_by_points(polygon: Polygon, divisor: float, cell_size: float):
+    # Extract polygon exterior coordinates
+    exterior_coords = list(polygon.exterior.coords)
+
+    # Initialize lists to store coordinates of new shapes after split
+    split_shape = []
+    remaining_shape = []
+
+    # Loop through points and check which side of the division line they are
+    if cell_size > 0:  # Horizontal split
+        prev_point = exterior_coords[0]
+        for point in exterior_coords:
+            point = Point(point)  # Convert coordinates to Shapely Point object
+
+            if point.y < divisor < prev_point.y:
+                remaining_shape.append(point)
+                remaining_shape.append((point.x, divisor))
+
+                split_shape.append((prev_point.x, divisor))
+                split_shape.append(prev_point)
+            elif point.y > divisor > prev_point.y:
+                split_shape.append(point)
+                split_shape.append((point.x, divisor))
+
+                remaining_shape.append((prev_point.x, divisor))
+                remaining_shape.append(prev_point)
+            elif point.y < divisor:  # Check if point is over or below divisor
+                split_shape.append(point)  # Append to appropriate shape
+            elif point.y > divisor:
+                remaining_shape.append(point)
+            prev_point = point
+
+        # Populate newly created edges with points to facilitate vertical cutting
+        populated_shapes = populate_new_edge(split_shape, remaining_shape, cell_size, divisor)
+        if populated_shapes is not None:
+            split_shape = populated_shapes[0]
+            remaining_shape = populated_shapes[1]
+
+    else:  # Vertical split
+        prev_point = exterior_coords[1]
+        for point in exterior_coords:
+            point = Point(point)  # Convert coordinates to Shapely Point object
+
+            if point.x < divisor < prev_point:
+                remaining_shape.append(point)
+                remaining_shape.append((divisor, point.y))
+
+                split_shape.append((divisor, prev_point.y))
+                split_shape.append(prev_point)
+            elif point.x > divisor > prev_point:
+                split_shape.append(point)
+                split_shape.append((divisor, point.y))
+
+                remaining_shape.append((divisor, prev_point.y))
+                remaining_shape.append(prev_point)
+            elif point.x < divisor:  # Check if point is over or below divisor
+                split_shape.append(point)  # Append to appropriate shape
+            elif point.x > divisor:
+                remaining_shape.append(point)
+            prev_point = point
+
+    # Check if the split_shape has enough coordinates to create a polygon
+    if len(split_shape) < 3:
+        # print("Not enough coordinates to create valid polygon: Split shape: ", len(split_shape))
+        split_shape = None
+    else:
+        split_shape.append(split_shape[0])  # Append first coord to create closed loop
+
+    # Check if the remaining_shape has enough coordinates to create a polygon
+    if len(remaining_shape) < 3:
+        # print("Not enough coordinates to create valid polygon: Remaining shape: ", len(remaining_shape))
+        remaining_shape = None
+    else:
+        remaining_shape.append(remaining_shape[0])  # Append first coord to create closed loop
+
+    # Return split polygons as Shapely Polygon objects
+    return Polygon(split_shape), Polygon(remaining_shape)
+'''
+
+
 def cut_polygon_by_points(polygon: Polygon, divisor: float, cell_size: float):
     # Extract polygon exterior coordinates
     exterior_coords = list(polygon.exterior.coords)
@@ -174,8 +256,80 @@ def circle_polygon():
     return Polygon(circle_points)
 
 
+def populate_new_edge(polygon, cell_size: float, divisor: float):
+    if polygon is not None and polygon_min_x is not None:
+
+        sub_division = []
+        final_shapes = []
+        last_sub_division = []
+        left_most_point = None
+
+        for i in range(len(polygon)):
+            x, y = polygon[i]
+            if x > divisor:
+                sub_division.append((x, y))
+            elif x < divisor:
+                last_sub_division.append((x, y))
+
+        if len(sub_division) < 1 or len(last_sub_division) < 1:
+            return None
+
+        # Find most left point in the polygon
+        for point in last_sub_division:
+            if not left_most_point:
+                left_most_point = point[0]
+            elif point[0] < left_most_point:
+                left_most_point = point[0]
+
+        # Sort sub_division from bottom to top
+        if sub_division[0][1] > sub_division[-1][1]:
+            sub_division = sorted(sub_division, key=lambda point: sub_division[1])
+
+        # Append corners to edge piece
+        sub_division.append((sub_division[-1][0] - cell_size, sub_division[-1][1]))
+        sub_division.append((sub_division[-1][0] - cell_size, sub_division[0][1]))
+
+        # Append first section to list as Shapely Polygon object
+        final_shapes.append(Polygon(sub_division))
+
+        # Save last point of previous tile (upper left corner) as the starting coordinate for a new tile
+        next_point = sub_division[-1]
+
+        # Clear subdivision list for next tile and add point as a reference
+        sub_division = [next_point]
+
+        # Continue making tiles on the horizontal section while the left edge is not reached
+        while sub_division[-1][0] > left_most_point:
+            # Append three more corners to new tile
+            sub_division.append((sub_division[0][0], sub_division[0][1] - cell_size))
+            sub_division.append((sub_division[0][0] - cell_size, sub_division[-1][1]))
+            sub_division.append((sub_division[-1][0] - cell_size, sub_division[0][0]))
+            # NB maybe have to close loop???
+
+            # Append new shape, save next starting coordinate, and clear sub_div for next tile
+            final_shapes.append(Polygon(sub_division))
+            next_point = sub_division[-1]
+            sub_division = [next_point]
+
+            print("sub_div: ", sub_division[-1], " left_most_point: ", left_most_point)
+
+        # Sort last coordinates of last tile from bottom to top
+        if last_sub_division[0][1] > last_sub_division[-1][1]:
+            last_sub_division = sorted(last_sub_division, key=lambda point: last_sub_division[1])
+
+        # Append two corners to complete the last tile
+        last_sub_division.append((sub_division[-1][0] - cell_size, sub_division[-1][1]))
+        last_sub_division.append((sub_division[-1][0] - cell_size, sub_division[-1][1] - cell_size))
+        # NB maybe have to close loop???
+
+        final_shapes.append(Polygon(last_sub_division))
+
+        return final_shapes  # Return list of all tiles for current horizontal section
+
+
+'''
 # Adds equally spaced points along an edge that is created after a polygon is cut.
-def populate_new_edge(split_shape, cell_size: float, divisor: float):
+def populate_new_edge(split_shape, remaining_shape, cell_size: float, divisor: float):
     # Prepend new points onto the newly created edge to facilitate vertical splitting
     if split_shape is not None and polygon_min_x is not None:
         # Define starting point with an x-value that will be common for all polygons
@@ -184,23 +338,31 @@ def populate_new_edge(split_shape, cell_size: float, divisor: float):
         # Create list of corners
         corners = []
 
+        divisor_range = ((divisor - 0.1),(divisor + 0.1)) # Define tolerance
+
         # Find all corners of split shape. Corner = point that intersects divisor
         for point in split_shape:
-            if point[1] == divisor:
+            if divisor_range[0] < point.y < divisor_range[1]:
                 corners.append(point)
 
         if not corners or len(corners) < 2:
             return None
 
         # Sort corners in ascending order (left to right) based on x coordinate
-        sorted_corners = sorted(corners, key=lambda p: p[0])
+        sorted_corners = sorted(corners, key=lambda p: p.x)
 
-        while starting_point < sorted_corners[0][0]:  # Increment starting point until it is withing the polygons bounds
+        while starting_point < sorted_corners[0].x:  # Increment starting point until it is withing the polygons bounds
             starting_point += cell_size
 
+        # if starting_point < left_corner.x: # NB: optimised substitute for previous while loop, requires testing
+        #    starting_point += cell_size * math.floor(starting_point - left_corner.x)
+
         # Insert new points with cell_size spacing while starting_point is within bounds
-        while starting_point < sorted_corners[-1][0]:
+        while starting_point < sorted_corners[-1].x:
             split_shape.insert(0, (starting_point, divisor))  # NB may have to add/subtract small offset of 0.00001
+            remaining_shape.insert(-1, (starting_point, divisor))  # Prepend new point to shape
+
             starting_point += cell_size
 
-        return split_shape
+        return split_shape, remaining_shape
+'''