Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
37 changes: 26 additions & 11 deletions dem_handler/utils/spatial.py
Original file line number Diff line number Diff line change
Expand Up @@ -593,21 +593,36 @@ def check_dem_type_in_bounds(
logger.info(f"Searching {dem_index_path}")

gdf = gpd.read_file(dem_index_path, layer=layer)
bounding_box = shapely.geometry.box(*bounds)
intersecting_tiles = 0

if gdf.crs is not None:
# ensure same crs
bounding_box = (
gpd.GeoSeries([bounding_box], crs="EPSG:4326").to_crs(gdf.crs).iloc[0]
# handle antimeridian crossing
if bounds[0] > bounds[2]:
logger.warning(
"bounds cross the antimeridian, searching for intersections either side."
)
bounds_east, bounds_west = split_bounds_at_antimeridian(bounds)
search_shapes = [
shapely.geometry.box(*bounds_east),
shapely.geometry.box(*bounds_west),
]
else:
logger.info('No crs found for index file. Assuming EPSG:4326"')
bounding_box = gpd.GeoSeries([bounding_box], crs="EPSG:4326")
# Find rows that intersect with the bounding box
intersecting_tiles = gdf[gdf.intersects(bounding_box)]
if len(intersecting_tiles) == 0:
search_shapes = [shapely.geometry.box(*bounds)]

for bounding_box in search_shapes:
if gdf.crs is not None:
# ensure same crs
bounding_box = (
gpd.GeoSeries([bounding_box], crs="EPSG:4326").to_crs(gdf.crs).iloc[0]
)
else:
logger.info('No crs found for index file. Assuming EPSG:4326"')
bounding_box = gpd.GeoSeries([bounding_box], crs="EPSG:4326")
# Count rows that intersect with the bounding box
intersecting_tiles += len(gdf[gdf.intersects(bounding_box)])

if intersecting_tiles == 0:
logger.info(f"No intersecting tiles found")
return False
else:
logger.info(f"{len(intersecting_tiles)} intersecting tiles found")
logger.info(f"{intersecting_tiles} intersecting tiles found")
return True
21 changes: 21 additions & 0 deletions tests/utils/test_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -214,6 +214,25 @@ class BoundsDEMCheckCase:
is_error=False,
)

# crosses AM over the coast, no intersect
test_antimeridian_no_intersect_with_rema = BoundsDEMCheckCase(
dem_type="rema",
resolution=32,
bounds=(173.430893, -71.618423, -178.032867, -68.765106),
in_bounds=False,
is_error=False,
)

# crosses AM over the coast, no intersect
test_antimeridian_no_intersect_with_cop30 = BoundsDEMCheckCase(
dem_type="cop30",
resolution=30,
bounds=(173.430893, -71.618423, -178.032867, -68.765106),
in_bounds=False,
is_error=False,
)


test_cases = [
test_invalid_dem_type,
test_invalid_dem_resolution,
Expand All @@ -224,6 +243,8 @@ class BoundsDEMCheckCase:
test_heard_island_bounds_with_cop30,
test_antimeridian_with_cop30,
test_antimeridian_with_rema,
test_antimeridian_no_intersect_with_rema,
test_antimeridian_no_intersect_with_cop30,
]


Expand Down