Skip to content

Commit

Permalink
[Fixes #11171] Fix calculation of wkt for regions crossing the dateli…
Browse files Browse the repository at this point in the history
…ne (#11172)

* Fix region wkt calculation when crossing dateline

* simplification

* fix PEP8

* fix failing test

* fixed wrong local path for fixtures

* restore test

* Keep WKT for regions not across IDL return Polygon
  • Loading branch information
giohappy authored Jun 20, 2023
1 parent 5ae4e18 commit 2c1e886
Show file tree
Hide file tree
Showing 3 changed files with 121 additions and 4 deletions.
17 changes: 17 additions & 0 deletions geonode/base/fixtures/regions.json
Original file line number Diff line number Diff line change
Expand Up @@ -4334,6 +4334,23 @@
"bbox_y1": 27.454559
}
},
{
"pk": 256,
"model": "base.region",
"fields": {
"rght": 515,
"code": "PAC",
"name": "Pacific",
"parent": null,
"level": 1,
"lft": 464,
"tree_id": 90,
"bbox_x0": 112.921112,
"bbox_x1": -108.87291,
"bbox_y0": -54.640301,
"bbox_y1": 20.41712
}
},
{
"pk": 257,
"model": "base.region",
Expand Down
69 changes: 67 additions & 2 deletions geonode/tests/test_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@
from unittest.mock import patch
from datetime import datetime, timedelta

from django.contrib.gis.geos import Polygon
from django.contrib.gis.geos import GEOSGeometry, Polygon
from django.contrib.auth import get_user_model
from django.core.management import call_command

Expand All @@ -37,7 +37,7 @@
from geonode.geoserver.helpers import set_attributes
from geonode.tests.base import GeoNodeBaseTestSupport
from geonode.br.management.commands.utils.utils import ignore_time
from geonode.utils import copy_tree, fixup_shp_columnnames, get_supported_datasets_file_types, unzip_file
from geonode.utils import copy_tree, fixup_shp_columnnames, get_supported_datasets_file_types, unzip_file, bbox_to_wkt
from geonode import settings


Expand Down Expand Up @@ -286,3 +286,68 @@ def test_should_replace_the_type_id_if_already_exists(self):
self.assertEqual(len(supported_keys), prev_count)
shp_type = [t for t in supported_types if t["id"] == "shp"][0]
self.assertEqual(shp_type["label"], "Replaced type")


class TestRegionsCrossingDateLine(TestCase):
def setUp(self):
self.test_region_coords = [
[112.921111999999994, -108.872910000000005, -54.640301000000001, 20.417120000000001], # Pacific
[174.866196000000002, -178.203156000000007, -21.017119999999998, -12.466220000000000], # Fiji
[158.418335000000013, -150.208359000000002, -11.437030000000000, 4.719560000000000], # Kiribati
[165.883803999999998, -175.987198000000006, -52.618591000000002, -29.209969999999998], # New Zeland
]

self.region_across_idl = [
112.921111999999994,
-108.872910000000005,
-54.640301000000001,
20.417120000000001,
] # Pacific
self.region_not_across_idl = [
-31.266000999999999,
39.869301000000000,
27.636310999999999,
81.008797000000001,
] # Europe

self.dataset_crossing = GEOSGeometry(
"POLYGON ((-179.30100799101168718 -17.31310259828852693, -170.41740336774466869 -9.63481116511300328, -164.30155495876181249 -19.7237289784715415, \
-177.91712988386967709 -30.43762400150689018, -179.30100799101168718 -17.31310259828852693))",
srid=4326,
)
self.dataset_not_crossing = GEOSGeometry(
"POLYGON ((-41.86461347546176626 5.43160371103088835, 34.20404118809119609 4.3602142087273279, 9.56208263510924894 -48.8521310723496498, \
-42.22174330956295307 -32.42415870369499942, -41.86461347546176626 5.43160371103088835))",
srid=4326,
)

def test_region_across_dateline_do_not_intersect_areas_outside(self):
for i, region_coords in enumerate(self.test_region_coords):
geographic_bounding_box = bbox_to_wkt(
region_coords[0], region_coords[1], region_coords[2], region_coords[3]
)
_, wkt = geographic_bounding_box.split(";")
poly = GEOSGeometry(wkt, srid=4326)

self.assertFalse(
poly.intersection(self.dataset_crossing).empty, f"True intersection not detected for region {i}"
)
self.assertTrue(poly.intersection(self.dataset_not_crossing).empty, "False intersection detected")

def test_region_wkt_multipolygon_if_across_idl(self):
bbox_across_idl = bbox_to_wkt(
self.region_not_across_idl[0],
self.region_not_across_idl[1],
self.region_not_across_idl[2],
self.region_not_across_idl[3],
)
_, wkt = bbox_across_idl.split(";")
poly = GEOSGeometry(wkt, srid=4326)
self.assertEqual(poly.geom_type, "Polygon", f"Expexted 'Polygon' type but received {poly.geom_type}")

bbox_across_idl = bbox_to_wkt(
self.region_across_idl[0], self.region_across_idl[1], self.region_across_idl[2], self.region_across_idl[3]
)
_, wkt = bbox_across_idl.split(";")
poly = GEOSGeometry(wkt, srid=4326)
self.assertEqual(poly.geom_type, "MultiPolygon", f"Expexted 'MultiPolygon' type but received {poly.geom_type}")
39 changes: 37 additions & 2 deletions geonode/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -496,9 +496,44 @@ def bbox_to_wkt(x0, x1, y0, y1, srid="4326", include_srid=True):
if srid and str(srid).startswith("EPSG:"):
srid = srid[5:]
if None not in {x0, x1, y0, y1}:
wkt = "POLYGON(({:f} {:f},{:f} {:f},{:f} {:f},{:f} {:f},{:f} {:f}))".format(
float(x0), float(y0), float(x0), float(y1), float(x1), float(y1), float(x1), float(y0), float(x0), float(y0)
polys = []

# We assume that if x1 is smaller then x0 we're crossing the date line
crossing_idl = x1 < x0
if crossing_idl:
polys.append(
[
(float(x0), float(y0)),
(float(x0), float(y1)),
(180.0, float(y1)),
(180.0, float(y0)),
(float(x0), float(y0)),
]
)
polys.append(
[
(-180.0, float(y0)),
(-180.0, float(y1)),
(float(x1), float(y1)),
(float(x1), float(y0)),
(-180.0, float(y0)),
]
)
else:
polys.append(
[
(float(x0), float(y0)),
(float(x0), float(y1)),
(float(x1), float(y1)),
(float(x1), float(y0)),
(float(x0), float(y0)),
]
)

poly_wkts = ",".join(
["(({}))".format(",".join(["{:f} {:f}".format(coords[0], coords[1]) for coords in poly])) for poly in polys]
)
wkt = f"MULTIPOLYGON({poly_wkts})" if len(polys) > 1 else f"POLYGON{poly_wkts}"
if include_srid:
wkt = f"SRID={srid};{wkt}"
else:
Expand Down

0 comments on commit 2c1e886

Please sign in to comment.