Fixes handling of points in Near.to

Author: bethard

src/normit/geo/ops.py

@@ -101,11 +101,15 @@ def to(geometry: shapely.geometry.base.BaseGeometry,
         # project to UTM where we can measure distance in meters
         proj = utm_proj(geometry)
         geometry = shapely.ops.transform(proj, geometry)
+        # if radius is specified, convert it to meters
+        if radius is not None:
+            radius_m = radius.to(UNITS.meter).magnitude
         # if no radius is specified, infer it from the geometry's area
-        if radius is None:
+        elif geometry.boundary is not None:
             radius_m = _radius_by_area_in_meters(geometry)
+        # if no radius can be inferred (e.g., a linestring), assume 1km
         else:
-            radius_m = radius.to(UNITS.meter).magnitude
+            radius_m = 1000
         # if no distance is specified, buffer by one radius
         if distance is None:
             result = geometry.buffer(radius_m)

test/data/openstreetmap/37/3739700410

@@ -0,0 +1,6 @@
+{
+"type": "FeatureCollection",
+"features": [
+{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ 12.5653736, 55.6724527 ] } }
+]
+}

test/data/openstreetmap/50/5096749954

@@ -0,0 +1,6 @@
+{
+"type": "FeatureCollection",
+"features": [
+{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ 12.564607, 55.6721858 ] } }
+]
+}

test/data/openstreetmap/91/91260821

@@ -0,0 +1,6 @@
+{
+"type": "FeatureCollection",
+"features": [
+{ "type": "Feature", "properties": { }, "geometry": { "type": "MultiPolygon", "coordinates": [ [ [ [ 12.563476, 55.6727045 ], [ 12.5638715, 55.6723517 ], [ 12.5640886, 55.6721559 ], [ 12.5641398, 55.6721707 ], [ 12.5641976, 55.6721899 ], [ 12.5642582, 55.6722106 ], [ 12.5642206, 55.6722424 ], [ 12.5642634, 55.6722569 ], [ 12.5643273, 55.6722787 ], [ 12.5643912, 55.6723004 ], [ 12.5645411, 55.6723514 ], [ 12.5645983, 55.6723709 ], [ 12.5646555, 55.6723903 ], [ 12.5648512, 55.6724569 ], [ 12.5651171, 55.6725474 ], [ 12.5653832, 55.6726382 ], [ 12.5656343, 55.6727233 ], [ 12.5657414, 55.6727597 ], [ 12.5658558, 55.6726523 ], [ 12.5666975, 55.6718619 ], [ 12.5667251, 55.6718716 ], [ 12.5669246, 55.6719328 ], [ 12.5668533, 55.6719908 ], [ 12.5666413, 55.6721917 ], [ 12.5666192, 55.6722104 ], [ 12.5664561, 55.672365 ], [ 12.5664423, 55.6723781 ], [ 12.5664476, 55.6723948 ], [ 12.566335, 55.6725093 ], [ 12.5663053, 55.6724993 ], [ 12.566051, 55.6727403 ], [ 12.5659762, 55.6728137 ], [ 12.5660255, 55.6728286 ], [ 12.5658149, 55.6730204 ], [ 12.5655494, 55.6732678 ], [ 12.5655138, 55.6732585 ], [ 12.5654346, 55.6733334 ], [ 12.5651457, 55.6736069 ], [ 12.5650944, 55.6736555 ], [ 12.5648684, 55.6735671 ], [ 12.5650796, 55.6733727 ], [ 12.5651829, 55.6732736 ], [ 12.5646179, 55.6730912 ], [ 12.5645641, 55.6731429 ], [ 12.5645198, 55.6731818 ], [ 12.5645025, 55.6731988 ], [ 12.5639659, 55.6730158 ], [ 12.5639835, 55.6729991 ], [ 12.5639912, 55.6729922 ], [ 12.5640817, 55.6729078 ], [ 12.5635226, 55.6727201 ], [ 12.563476, 55.6727045 ] ] ] ] } }
+]
+}

test/test_geo_ops.py

@@ -116,6 +116,10 @@ def test_outer_regions(georeader: GeoJsonDirReader):
     albuquerque = georeader.read(171262)
     assert Near.to(colorado_river).intersection(albuquerque).area == 0
 
+    # near to Copenhagen central station contains the station
+    copenhagen_central_station = georeader.read(3739700410, 5096749954, 91260821)
+    assert Near.to(copenhagen_central_station).intersection(copenhagen_central_station).area > 0
+
 
 def test_distance(georeader: GeoJsonDirReader):
     de = georeader.read(51477)