diff --git a/CHANGELOG.md b/CHANGELOG.md index 0d976dc..4df7e59 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,6 +6,17 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## Unreleased +## 0.4.2 - 2020-09-18 +## Added: +- faster implementation of spatial selection + +## Changed: +- all elements in WKT and shapefiles containing multiple polygons are considered (not only the first element) +- the validity conditions for multipolygons have been relaxed: valid adjacent polygons are now accepted + +## Fixed: +- bug in copy_point_cloud with masks and zero-filled feature-arrays + ## 0.4.1 - 2020-05-20 ## Added: - select_equal filter accepts list of values to compare to the points' attributes diff --git a/CITATION.cff b/CITATION.cff index 6dc4632..2bfe962 100644 --- a/CITATION.cff +++ b/CITATION.cff @@ -52,7 +52,7 @@ authors: family-names: Koma given-names: Zsófia cff-version: "1.0.3" -date-released: 2020-05-20 +date-released: 2020-09-18 doi: "10.5281/zenodo.1219422" keywords: - "airborne laser scanning" @@ -62,5 +62,5 @@ keywords: license: "Apache-2.0" message: "If you use this software, please cite it using these metadata." title: "Laserchicken: toolkit for ALS point clouds" -version: "0.4.1" +version: "0.4.2" ... diff --git a/laserchicken/_version.txt b/laserchicken/_version.txt index 267577d..2b7c5ae 100644 --- a/laserchicken/_version.txt +++ b/laserchicken/_version.txt @@ -1 +1 @@ -0.4.1 +0.4.2 diff --git a/laserchicken/filter.py b/laserchicken/filter.py index 21b9d1d..2e1e804 100644 --- a/laserchicken/filter.py +++ b/laserchicken/filter.py @@ -9,6 +9,7 @@ from shapely.errors import WKTReadingError from shapely.wkt import loads from shapely.geometry import box +from shapely.vectorized import contains import numpy as np from laserchicken.keys import point @@ -115,9 +116,9 @@ def select_polygon(point_cloud, polygon_string, read_from_file=False, return_mas else: polygon = _load_polygon(polygon_string) - if isinstance(polygon, shapely.geometry.polygon.Polygon) and polygon.is_valid: + if isinstance(polygon, shapely.geometry.polygon.Polygon): points_in = _contains(point_cloud, polygon) - elif isinstance(polygon,shapely.geometry.multipolygon.MultiPolygon) and polygon.is_valid: + elif isinstance(polygon,shapely.geometry.multipolygon.MultiPolygon): points_in = [] count=1 for poly in polygon: @@ -145,17 +146,17 @@ def _read_wkt_file(path): with open(path) as f: content = f.readlines() - content = [x.strip() for x in content] - return _load_polygon(content[0]) + content = [_load_polygon(x.strip()) for x in content] + geom = shapely.geometry.MultiPolygon(content) if len(content) > 1 else content[0] + return geom def _read_shp_file(path): shape = shapefile.Reader(path) - # first feature of the shapefile - feature = shape.shapeRecords()[0] - first = feature.shape.__geo_interface__ - # or shp_geom = shape(first) with PyShp) - shp_geom = shapely.geometry.shape(first) + features = shape.shapeRecords() + shp_geoms = [shapely.geometry.shape(feature.shape.__geo_interface__) + for feature in features] + shp_geom = shapely.geometry.MultiPolygon(shp_geoms) if len(shp_geoms) > 1 else shp_geoms[0] return shp_geom @@ -192,6 +193,9 @@ def _contains(pc, polygon): y = pc[point]['y']['data'] points_in = [] + if not polygon.is_valid: + raise ValueError('Invalid polygon in input') + mbr = polygon.envelope point_box = box(np.min(x), np.min(y), np.max(x), np.max(y)) @@ -204,10 +208,8 @@ def _contains(pc, polygon): tree = kd_tree.get_kdtree_for_pc(pc) indices = np.sort(tree.query_ball_point(x=p, r=rad)) - point_id = 0 - for i in indices: - if polygon.contains(Point(x[i], y[i])): - points_in.append(i) - point_id += 1 + if len(indices) > 0: + mask = contains(polygon, x[indices], y[indices]) + points_in.extend(indices[mask]) return points_in diff --git a/laserchicken/test_filter.py b/laserchicken/test_filter.py index de6fea6..b2b7e2b 100644 --- a/laserchicken/test_filter.py +++ b/laserchicken/test_filter.py @@ -220,13 +220,6 @@ def test_points_in_polygon_wkt_MultiLine(): with pytest.raises(ValueError): select_polygon(pc_in, "MULTILINESTRING((3 4,10 50,20 25),(-5 -8,-10 -8,-15 -4))") - @staticmethod - def test_points_in_polygon_wkt_MultiPolygon(): - pc_in = load("testdata/AHN2.las") - with pytest.raises(ValueError): - select_polygon(pc_in, - "MULTIPOLYGON(((1 1,5 1,5 5,1 5,1 1),(2 2, 3 2, 3 3, 2 3,2 2)),((3 3,6 2,6 4,3 3)))") - @staticmethod def test_points_in_polygon_wkt_Collection(): pc_in = load("testdata/AHN2.las") @@ -330,12 +323,6 @@ def test_points_in_polygon_wkt_MultiLine(): with pytest.raises(ValueError): select_polygon(pc_in, "testdata/ahn2_geometries_wkt/multiline.wkt", read_from_file=True) - @staticmethod - def test_points_in_polygon_wkt_MultiPolygon(): - pc_in = load("testdata/AHN2.las") - with pytest.raises(ValueError): - select_polygon(pc_in, "testdata/ahn2_geometries_wkt/multipolygon.wkt", read_from_file=True) - @staticmethod def test_points_in_polygon_wkt_Collection(): pc_in = load("testdata/AHN2.las") diff --git a/laserchicken/utils.py b/laserchicken/utils.py index de538c5..ec4ebe1 100644 --- a/laserchicken/utils.py +++ b/laserchicken/utils.py @@ -122,7 +122,8 @@ def copy_point_cloud(source_point_cloud, array_mask=None): new_value = copy_point_cloud(value, array_mask) elif isinstance(value, np.ndarray): if array_mask is not None: - new_value = value[array_mask] if any(value) else np.copy(value) + new_value = (value[array_mask] + if value.size > 0 else np.copy(value)) else: new_value = np.copy(value) else: diff --git a/testdata/ahn2_geometries_wkt/multipolygon.wkt b/testdata/ahn2_geometries_wkt/multipolygon.wkt deleted file mode 100644 index fad63cf..0000000 --- a/testdata/ahn2_geometries_wkt/multipolygon.wkt +++ /dev/null @@ -1 +0,0 @@ -MULTIPOLYGON(((1 1,5 1,5 5,1 5,1 1),(2 2, 3 2, 3 3, 2 3,2 2)),((3 3,6 2,6 4,3 3))) \ No newline at end of file