Skip to content

Commit

Permalink
Merge branch '920-simplification-du-check-input' into 'master'
Browse files Browse the repository at this point in the history
Resolve "simplification du check input"

Closes #920

See merge request 3d/cars-park/cars!764
  • Loading branch information
dyoussef committed Nov 25, 2024
2 parents 8f74306 + 8b4fc0e commit 956d423
Show file tree
Hide file tree
Showing 2 changed files with 66 additions and 103 deletions.
137 changes: 60 additions & 77 deletions cars/core/projection.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,100 +45,83 @@


def compute_dem_intersection_with_poly( # noqa: C901
srtm_dir: str, ref_poly: Polygon, ref_epsg: int
srtm_file: str, ref_poly: Polygon, ref_epsg: int
) -> Polygon:
"""
Compute the intersection polygon between the defined dem regions
and the reference polygon in input
:raise Exception: when the input dem doesn't intersect the reference polygon
:param srtm_dir: srtm directory
:param srtm_file: srtm file
:param ref_poly: reference polygon
:param ref_epsg: reference epsg code
:return: The intersection polygon between the defined dem regions
and the reference polygon in input
"""
dem_poly = None

srtm_files_dir = []
if os.path.isfile(os.path.abspath(srtm_dir)):
srtm_files_dir.append([srtm_dir])
else:
for _, _, srtm_files in os.walk(srtm_dir):
srtm_files_dir.append(srtm_files)

for srtm_files in srtm_files_dir:
logging.info(
"Browsing all files of the srtm dir. "
"Some files might be unreadable by rasterio (non blocking matter)."
)
for file in srtm_files:
unsupported_formats = [".omd"]
_, ext = os.path.splitext(file)
if ext not in unsupported_formats:
if inputs.rasterio_can_open(os.path.join(srtm_dir, file)):
with rio.open(os.path.join(srtm_dir, file)) as data:
xmin = min(data.bounds.left, data.bounds.right)
ymin = min(data.bounds.bottom, data.bounds.top)
xmax = max(data.bounds.left, data.bounds.right)
ymax = max(data.bounds.bottom, data.bounds.top)

try:
file_epsg = data.crs.to_epsg()
file_bb = Polygon(
[
(xmin, ymin),
(xmin, ymax),
(xmax, ymax),
(xmax, ymin),
(xmin, ymin),
]
)

# transform polygon if needed
if os.path.isdir(os.path.abspath(srtm_file)):
raise RuntimeError("srtm input should be a file and not a directory")

if inputs.rasterio_can_open(srtm_file):
with rio.open(srtm_file) as data:
xmin = min(data.bounds.left, data.bounds.right)
ymin = min(data.bounds.bottom, data.bounds.top)
xmax = max(data.bounds.left, data.bounds.right)
ymax = max(data.bounds.bottom, data.bounds.top)

try:
file_epsg = data.crs.to_epsg()
file_bb = Polygon(
[
(xmin, ymin),
(xmin, ymax),
(xmax, ymax),
(xmax, ymin),
(xmin, ymin),
]
)

# transform polygon if needed
if ref_epsg != file_epsg:
file_bb = polygon_projection(file_bb, file_epsg, ref_epsg)

# if the srtm tile intersects the reference polygon
if file_bb.intersects(ref_poly):
local_dem_poly = None

# retrieve valid polygons
for poly, val in shapes(
data.dataset_mask(),
transform=data.transform,
):
if val != 0:
poly = shape(poly)
poly = poly.buffer(0)
if ref_epsg != file_epsg:
file_bb = polygon_projection(
file_bb, file_epsg, ref_epsg
poly = polygon_projection(
poly, file_epsg, ref_epsg
)

# if the srtm tile intersects the reference polygon
if file_bb.intersects(ref_poly):
local_dem_poly = None

# retrieve valid polygons
for poly, val in shapes(
data.dataset_mask(),
transform=data.transform,
):
if val != 0:
poly = shape(poly)
poly = poly.buffer(0)
if ref_epsg != file_epsg:
poly = polygon_projection(
poly, file_epsg, ref_epsg
)

# combine valid polygons
if local_dem_poly is None:
local_dem_poly = poly
else:
local_dem_poly = poly.union(
local_dem_poly
)

# combine the tile valid polygon to the other
# tiles' ones
if dem_poly is None:
dem_poly = local_dem_poly
else:
dem_poly = dem_poly.union(local_dem_poly)

except AttributeError as attribute_error:
logging.warning(
"Impossible to read the SRTM"
"tile epsg code: {}".format(attribute_error)
)
# combine valid polygons
if local_dem_poly is None:
local_dem_poly = poly
else:
local_dem_poly = poly.union(local_dem_poly)

# combine the tile valid polygon to the other
# tiles' ones
if dem_poly is None:
dem_poly = local_dem_poly
else:
dem_poly = dem_poly.union(local_dem_poly)

except AttributeError as attribute_error:
logging.warning(
"Impossible to read the SRTM"
"tile epsg code: {}".format(attribute_error)
)

# compute dem coverage polygon over the reference polygon
if dem_poly is None or not dem_poly.intersects(ref_poly):
Expand Down
32 changes: 6 additions & 26 deletions tests/core/test_projection.py
Original file line number Diff line number Diff line change
Expand Up @@ -90,34 +90,12 @@ def test_compute_dem_intersection_with_poly():
)

dem_inter_poly, cover = projection.compute_dem_intersection_with_poly(
absolute_data_path("input/phr_ventoux/srtm"), inter_poly, inter_epsg
)
assert dem_inter_poly.equals(inter_poly)
assert cover == 100.0

# test partial coverage over with several srtm tiles with no data holes
inter_poly = Polygon(
[(4.8, 44.2), (4.8, 44.3), (6.2, 44.3), (6.2, 44.2), (4.8, 44.2)]
)
dem_inter_poly, cover = projection.compute_dem_intersection_with_poly(
absolute_data_path("input/utils_input/srtm_with_hole"),
absolute_data_path("input/phr_ventoux/srtm/N44E005.hgt"),
inter_poly,
inter_epsg,
)

ref_dem_inter_poly = Polygon(
[
(4.999583333333334, 44.3),
(6.2, 44.3),
(6.2, 44.2),
(4.999583333333334, 44.2),
(4.999583333333334, 44.3),
]
)

assert ref_dem_inter_poly.exterior.equals(dem_inter_poly.exterior)
assert len(list(dem_inter_poly.interiors)) == 6
assert cover == 85.72172619047616
assert dem_inter_poly.equals(inter_poly)
assert cover == 100.0

# test no coverage
inter_poly = Polygon(
Expand All @@ -126,7 +104,9 @@ def test_compute_dem_intersection_with_poly():

with pytest.raises(Exception) as intersect_error:
dem_inter_poly, cover = projection.compute_dem_intersection_with_poly(
absolute_data_path("input/phr_ventoux/srtm"), inter_poly, inter_epsg
absolute_data_path("input/phr_ventoux/srtm/N44E005.hgt"),
inter_poly,
inter_epsg,
)
assert (
str(intersect_error.value) == "The input DEM does not intersect "
Expand Down

0 comments on commit 956d423

Please sign in to comment.