From 59a0e44e0798f58b2bb1e86daaab0459b7806720 Mon Sep 17 00:00:00 2001
From: Sara <sarasdj@stud.ntnu.no>
Date: Tue, 12 Mar 2024 10:44:41 +0100
Subject: [PATCH] update: approach 4, grid of coordinates

---
 .../__pycache__/get_relation.cpython-311.pyc  | Bin 5828 -> 5315 bytes
 server/map/get_relation.py                    |  98 ++++++++----------
 2 files changed, 45 insertions(+), 53 deletions(-)

diff --git a/server/map/__pycache__/get_relation.cpython-311.pyc b/server/map/__pycache__/get_relation.cpython-311.pyc
index 812af27f4b3e7b4696370419efaf670e47b9b2ff..76c8ff54db3931e76b092b850b48ecad05d4aa01 100644
GIT binary patch
delta 2532
zcmb7FO>7&-6`tk(lFOeJ$>onInv`uRW)jhHq}Gv)Ku#*Vu2ZOJY!rqI!(cTlnF`IN
zmun|hm!T3Q0MY6oW!zW=5`YKC5fC4Q_TrNfpgr`sODtf5fB~V654j07a>*g++oe{K
zl(uLecHjKH`R2{c_ul+<>@S7qJ`aTg2;?6}wv{iVP&y4yxSh(g->p5ja7>{wuf(yS
z$nzZbf52Q$U=I2MN_?KdBD9m3x_}D7wR^fKz)TwcXFq1h5$bZ$pckE{D_TVNdv$~`
zbDu)!3;6r;K80D*q(++it+<0$3#h>w%nZ8naf8!EZ#F^(hxrG>eIH=L-<Jk&@YtXB
zny4W_JFutiHv~YcH3E9bGio|eMgWWV!SOG4oo|SSxDP651T|rQU!7VRHT(cL-3WaT
z2zWswgo6k0!bW&69?Tsw*fx)d7L3y{MAm%SU+?ubxN;hn+uOsf+*rA*@17t+i_Y=8
zVz3WF`%kVB(PQ2ySO&-%lOQ=-q(pKf828~Yd5w;VLC>Tea=KG`*6V;d5wDIT<oEP2
zC4q@t7{C%9z)>7~5WmmDk-^^=3hi5Q09H)Xk#^JpCK)D1#hjXe<-VijZRRk86F5n3
zz?i&a6_zMupv8m-+J@nvOh?-($H0eE;K?BlBOqX!_}S;5UpjWSRIVsEud8|G9bGHU
z>G_&kxxJuP^SGp!^7EQ{JwFd7WSspw6(<Q!^<6wmKIHyDrAda*jwsh_`t4reSEb5L
z<&37OS}!`MDJ5OGx}cTu)j3twaQ)Bc`9D#CbKsW4>G~zY`HPt?K4$ajZA7tWs3+l+
znf{3tp0UF-rZ7V${P`)0opQThvBHyfc+wOmx1#aoH?}zxAKFHT*rTp$@!xPLacqNs
zIBv?*mOO3C(`3z`Y07>u0M$AN27@nd{B-luL)}Wx*y)+p`EF*UUHatJ*3Y+MIq)Eq
zZZB;P-)p$lu9Si`vZLF`9~yJjs_-Oxd^7m)+TU0!J8NfWR|8wg;r8$yXTm_X8z0(|
zGLNOh9qDj;YNK@TRZBYYQ6R^?0Ly<<BqM?EhO`0Lc#p^BTje^OQQ(@Ym6z12UaHK|
z<VGZ(zXWl37wzRIT;wiW7ky-nPRurF2x!I)Y8oPE8G@Sk?%IjO10?1QpLWvV78%dR
zpx$z;TEbx3tiz4ufjJafnSjW0$Hsm*b{^2Ye24xgw6{Do`G;ul=@F8X27`FVt$GYP
zei<Bm$?MPtbm!c51pD^yk3yKffV`d13}zq=Waz&437^{>^yx8Y4W|I)bMrpnK0Euh
zld6`E&#HP}sj4>@uH{{#ma9PAx)Th24(D%`DrM|NM7`jfoz<LtSqP9vp>&hF-Q!BN
zno`9*>YW~4QK~(*rj@I@=D-(t_tLeWL0vl!suzI}Q?#<G?Sw^7tk)`K{c8PMsiyQI
znsU8Vu9m9{I~{CI1+2Pr!}-2w@b(~WWJfh2jv<<6fF;Z5X<~5wq@8)e&b(|T#_h!T
z@-KkJF?sooZX)H@qr=-sWQ(q9u}>s<Rapz(3AWg-C<D#?#cm|N^7eYN6Um#Ae0P`8
zKeD8PEfq}n+mdpRrO}Qw+P-Nur7>F?Gu?0N|D3yLDukyXwoX}m(&m%v!yP_j@|jQT
z?P)7FZs*1~OP$<FGk5Zd5HiElmT=k@PMgB%ZY14oEx{I9=`oY&Fp1UQwolx7*J4I(
zX4GUxyG-mcGuUAU*QYm5*qJelDcVfYWQtG4=;|p;%-Le@v*bq38lJR=CrxqEWG1z^
zbFwn1x!qS08cvROF7(+~pzb;e7D^EE&w+*Je+xjMv5Nz|-W#U;ej0mL2VVRGWPhK6
zhn9EgOZQ2*l$urZNhiw!K^>yx7nG`^0fTV?;x{IXZ-SHB5l{tzd>or1li?xm7BuU1
zawvZ4^{0HG^~!sV4;t%KhfkY)`qK;T)E}>Wa%H2`$sRYe$KCWhZiOf8@PsK$Y>Bbe
zOQx7LnXG$Kg%pX(Z&JsIB99lx0Jq0hZ(gt6)-F4Bx~^-Qlk$7)d_^tkJ?1y6TIrMb
zD0w7j0~Hu_GF<&z@>FguouQ}`?YV8_l>bIw^yR3jHouLa*c@?-FE6^qL3NfAs6&n=
O6q{$<;y)^PBmV-Dc}v^?

delta 3112
zcmZuzU2NOd6~2^6iV`L2$NKS471K!^L2;cd$<n$m>@Kxqw@KivGtv%OQcpG+S&Ach
zm-H4J6r=%t&<X>w%MCGx7SPLFqe%g;1qK5{feqc>#eje{6~>Gf*~4DyEi166b?1_n
zV!1isJ-qjxpL_1*ch2SK$-nd`KJ@uK2xy@HYdPMA&{t&8Ajf{^$=%gNI{FT}9+q6d
zBN7KZDp8Xt>0bSV%DX_%3GK<3PW?P}-67HXXH;)7s1m_GNz6Kh&;w{TQAk3QD2)n?
zMvtPKr2?zE8#076R$?Ez)**lgjUZ6LrMV=oqp4KjVC+6Hc557Xx?b?8JW+bn*3?+d
zbI1{|@hy(Lni~wy7kswizMWU|9kLTNro!s0&M>VB5~p7fI^TxQ!4SO5FM#iZGwkwD
zmVc=>w~bJ?c_<q5(9?P-1-}}!S0_B{*59C`-WzZNg+L*w1uGs{Jkrg#8c|3KX006)
zql%XV*sgYwflXPV(UZLTznndkpQxPFf|5@XB)=4R7+h!Iq@meFp_VWL5JndjBu4tT
zw4mRk87f2~35E_wqOZ~&?lviGADsSs`Y2em@sWn-N|bn^wNNES9Oxg>e*If#fF)Q+
zM)f>%QvaO!64qL0^Lotn8>(IZ*cA`W<4jJi^D|gZsq$1BXQWTRc6~y5#=$Q^jS|rx
zxKHW-af?F^#=pb)szSfPoie#oRqj-c6UyvvC|tU{$D;6adnn2rw^-TrI3iYJn<Lf8
zfDsumBLn(`r)xxTgBV!(wZ|QO8{?ZlHadpQj^Xm9-LBqB>iwbe74Q>3v4<F6zr|L#
zo%rcZ@7AmlA2;LUE1sIZ?fD)X2Or;vx{v=~ai$4k2$qIsWy`+-75C9WY-mV5pxKH(
z9BI*{RN`l!S~V4aOX5-M$wQgfl4&U0{O(y+e}jpBzu<thbV8#w2joSzz^MMlN|2g|
z5*ZShs&E!KKJ^8Cj@y<h{sGLo$!zqn0gh@;3(4X3!y*bTU(;w?IPnbOZOw@t8hvPk
z&ci5F$?@)!wlCFtf>Av!1d32bSj~0I0f3_?md?MPn_o~B@y>#xiqo<<Ih&f7=N82u
z=I0jE`J9-_N#d1EPM%Z&f@v{7ot0-)CH=?Gw*U5b?TewM{;N4z%wusnug;43TcSEE
zi&$0`=2W7t5O0c@!m=penO7H+e7)0?Q8GD2P32}}Oyb1k`0-I-`WHe+k&`mEu+Gio
zb7~57>uhR%Ud~BL4wLFor{^G`x?jGl$~cq9Q#1KImXx~7T5ak~qjR>-%;!nou=B7P
zn77|&6%L$;`-DSe_2UgTUWA_F1-z8dV?MHy{@@n|pD_8v9-^3W>f@!Y<A3!Sr!SkQ
zFWWO#u3a>^F_Rmsa%11HDAu+1qS@7FcAYh%=gjE2(p5+of3P%O69T1i=!bgt5X(UR
z04p;)fyj!y>ixC1%+z@Qz2eZjLu(O(?>71FO79~+QRNdgf7`vnyM?uDhX0u9KemUQ
zz7w@T?7r}xu(n_Xdd)y@)p~ZnGx=Xw&)Y^qePzws!z)(6Lk2f&a>G?_7y=d|D{mRT
zxao_R#&(Ln@YB&QA6$8%JY1zaa2U?z*;)JvNbVyIW#Mg;MBTQs7&|qTrV8jF93{s)
zpgW4zLG?ky0BKa3T7uC!b>-k4S#YS$+=Rq)wBo$rJecMUoHu+HXw|LX6grE5bLWJG
zw3!39kX3=;zDG}<19l%SVh!DN7buN-x{p=IvS$)ua(g7HIyaqP$Vn-@m~`L(Ox0-u
zD~0S-6m1w0e@ls(1$ByuAbqB$7nv!7R(+KX7w15%bOI}(YHZMcHjb~+@B5lV>ztDe
z?g1@KE=o+mib;*fmw<tYqyWAOtj^w%aaxw@u4%x(j1@cxvkEye3*1&r^|od>MJpaF
z9DD<mGDQ2Qz)Gm*5=>Y7+Ut*8y;WClt*vt<S$6ICqIO)H33wZ^K{GbEC2hTF#6}F?
zsOcL8xCk9<Ck&y-6naYIJN}NfgyBDG`j3`g+u{8y+^>e04}d2MOErk>%DMw?Q2JK^
zH4Md;gXY#)D7YsJ9?e~GlH%FfQmQ?g^Oi&NXv{4~+5tED5n;Jcb6FO=#y!izq48ir
z6K-1+!GePt>V_KH{AWN&qFS*9{}i0&wu+=<Vkt3^7iaT0vqTEg94U7{m$90efj;@T
zifCWVQj)E+7D|-5BR3Ca3#W2v8Bf4w8&!^mZ|K#1p}~oI;HgXNqbPNu;dXya-K|rL
ziWkhS54z86tQ{jx2~$Y>r=AJDT&%@9s&u#(iIpZG?}aGI`%%i~Wva$Sx4HNuF8)Dw
z<IIM<adPu#TfLjvYVv|{<W=*?s|Gh>awApiX_}8&=3^8wZ{Oy+A93Bk^K4w(ys&X;
zEAl65%Tw(eF~m_*95uMXizaumYCXGr7g_eD%CPmaiq4A$H(+uDRc-*JYj0QizAD{k
z13sqzBlIIGp}*92?!@!3pE{FUxHG?qr%2zas2Gz!9lJ~u=;CAgZ`<M?vUQTv%C`RJ
ow(?sqQ)enOdk9{KnMTE_&OOv*qZCWVO%@JO&+nn{GAp|O1JA<Gr2qf`

diff --git a/server/map/get_relation.py b/server/map/get_relation.py
index b17ab008..15d0590e 100644
--- a/server/map/get_relation.py
+++ b/server/map/get_relation.py
@@ -13,7 +13,7 @@ def get_relation(self, body_of_water: str):
     # Filter only polygons, exclude points and other feature types to reduce response size
     polygon_data = geo_data[geo_data['geometry'].geom_type == 'Polygon']
 
-    if not polygon_data:
+    if polygon_data.empty:
         raise ValueError("Failed to extract polygon data from file")
 
     # Extract coordinates from polygons and create polygon objects
@@ -26,22 +26,30 @@ def get_relation(self, body_of_water: str):
 
     # Divide all polygons into sections
     for polygon in polygons:
-        # Initialise grid
-        grid_lines = create_grid(polygon, cell_size=0.1)
 
-        hrz_lines = grid_lines[0]  # Horizontal lines
-        vrt_lines = grid_lines[1]  # Vertical lines
+        # Divide the length and with of polygon into a grid
+        grid_lines = create_grid_coords(polygon, cell_size=0.1)
+
+        hrz_lines = grid_lines[0]  # Horizontal coordinates
+        vrt_lines = grid_lines[1]  # Vertical coordinates
 
         # 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)
 
+            # Split shape into upper and lower section as hrz_line as divider
+            cut_poly_1 = cut_polygon_in_two(polygon, hrz_line, True)
             polygon_part = cut_poly_1[0]  # Save upper horizontal section
 
+            if not polygon_part or not cut_poly_1[0]:
+                continue
+
             # Cut each horizontal section into vertical sections
             for vrt_line in vrt_lines:
-                cut_poly_2 = cut_polygon_in_two(polygon_part, vrt_line)
+                cut_poly_2 = cut_polygon_in_two(polygon_part, vrt_line, False)
+
+                if not cut_poly_2[0]:
+                    continue
+
                 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
@@ -70,66 +78,50 @@ def get_relation(self, body_of_water: str):
     self.wfile.write(tiles_json.encode('utf-8'))
 
 
-def cut_polygon_in_two(polygon, line):
-    # Check for invalid parameters
-    if not isinstance(polygon, Polygon) or not isinstance(line, LineString):
-        print("Inputs must be Shapely Polygon and LineString obects")
-
+# Takes a polygon and divides its coordinates into two shapes, where divisor is a
+# coordinate that defines the point of division
+def cut_polygon_in_two(polygon: Polygon, divisor: float, horizontal: bool):
     # Extract polygon exterior coordinates
     exterior_coords = list(polygon.exterior.coords)
 
-    # Initialize lists to store points for the two shapes
-    shape_1 = []
-    shape_2 = []
+    # Initialize lists to store coordinates of new shapes after split
+    split_shape = []
+    remaining_shape = []
 
-    # Split points into separate shapes based on their location relative to the line
+    # Loop through points and check which side of the division line they are
     for point in exterior_coords:
         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 belongs to shape 1 or shape 2
-                shape_1.append(point)
+        if horizontal:  # Horizontal split
+            if point.y < divisor:
+                split_shape.append(point)
             else:
-                shape_2.append(point)
-
-    # Check if shapes are empty
-    if not shape_1 or not shape_2:
-        raise ValueError("One or both shapes are empty after cutting")
-
-    # Create Polygon objects from the lists of points
-    poly_1 = Polygon(shape_1)
-    poly_2 = Polygon(shape_2)
+                remaining_shape.append(point)
+        else:  # Vertical split
+            if point.x < divisor:
+                split_shape.append(point)
+            else:
+                remaining_shape.append(point)
 
-    return poly_1, poly_2
+    # Check if polygons have enough coordinates
+    if len(split_shape) < 3 or len(remaining_shape) < 3:
+        print("Not enough coordinates to create valid polygons")
+        return None, None
 
+    # Append the first coordinate of the shapes to create closed loop
+    split_shape.append(split_shape[0]) # NB: may not be necessary?
+    remaining_shape.append(remaining_shape[0])
 
-# 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]
-    merged = linemerge([poly.boundary for poly in polygon] + [line])
-    borders = unary_union(merged)
-    polygons = polygonize(borders)
-    return list(polygons)
+    return Polygon(split_shape), Polygon(remaining_shape)  # Return split polygons as Shapely Polygon objects
 
 
-# Create a grid of lines spaced by a cell size
-def create_grid(polygon, cell_size):
+# Generate grid of equally spaced x and y coordinates where the grid size is determined by cell_size
+def create_grid_coords(polygon: Polygon, cell_size: float):
     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)]))
+    if len(x_coords) == 0 or len(y_coords) == 0:
+        raise ValueError("No grid points generated")
 
-    if not horizontal_lines or not vertical_lines:
-        raise ValueError("List of horizontal or vertical lines is empty")
+    return x_coords, y_coords
 
-    return horizontal_lines, vertical_lines
-- 
GitLab