Rounding all coordinates in shapely?

Question

Following on from How can I add points to a LineString in shapely?, I'm hitting a problem where floating point representations mean objects aren't precisely where I expect them be after they are moved using affinity.translate.

This would be fixed if there was a way of rounding all the coordinates to the nearest cm (or other insignificantly-small value).

Does this exist?

For what it's worth, I'm doing all my intersection and union operations first before the translation and that seems to be working.

gene · Accepted Answer · 2022-02-01T09:17:12.680

Shapely and GEOS cannot reduce precision (floating-point precision problem) but you can use other functions as numpy.round() or round(),with the GeoJSON format.

For polygons

from shapely.geometry import shape, mapping
import numpy as np
# one polygon
print poly.wkt
POLYGON ((101.2200527190607 -114.6493019170767, 146.6225142079163 -114.4488495484725, 185.0301801801801 -114.2792792792793, 184.6581081081081 -217.7153153153153, 14.99324324324321 -218.4594594594595, 16.4815315315315 -115.0234234234234, 101.2200527190607 -114.6493019170767))
# convert to GeoJSON format
geojson = mapping(poly)
print geojson
{'type': 'Polygon', 'coordinates': (((101.2200527190607, -114.6493019170767), (146.6225142079163, -114.4488495484725), (185.0301801801801, -114.2792792792793), (184.6581081081081, -217.7153153153153), (14.99324324324321, -218.4594594594595), (16.4815315315315, -115.0234234234234), (101.2200527190607, -114.6493019170767)),)}
 geojson['coordinates'] = np.round(np.array(geojson['coordinates']),2)
 print geojson
 {'type': 'Polygon', 'coordinates': array([[[ 101.22, -114.65],
    [ 146.62, -114.45],
    [ 185.03, -114.28],
    [ 184.66, -217.72],
    [  14.99, -218.46],
    [  16.48, -115.02],
    [ 101.22, -114.65]]])}
 print shape(geojson)
 POLYGON ((101.22 -114.65, 146.62 -114.45, 185.03 -114.28, 184.66 -217.72, 14.99 -218.46, 16.48 -115.02, 101.22 -114.65))

If you don't want to use Numpy, you can adapt the function def _set_precision(coords, precision) of the geojson-precision module

 def set_precision(coords, precision):
    result = []
    try:
        return round(coords, int(precision))
    except TypeError:
        for coord in coords:
            result.append(set_precision(coord, precision))
    return result
 geojson = mapping(poly)
 geojson['coordinates']= set_precision(geojson['coordinates'], 2)
 print geojson
 {'type': 'Polygon', 'coordinates': [[[101.22, -114.65], [146.62, -114.45], [185.03, -114.28], [184.66, -217.72], [14.99, -218.46], [16.48, -115.02], [101.22, -114.65]]]}
 print shape(geojson)
 POLYGON ((101.22 -114.65, 146.62 -114.45, 185.03 -114.28, 184.66 -217.72, 14.99 -218.46, 16.48 -115.02, 101.22 -114.65))

isn't there a typo in your alternate code: result.append(_set_precision(coord, precision)) should be result.append(set_precision(coord, precision)) ??? — Riccardo, Jan 31 '22 at 21:48

score 16 · Answer 2 · edited Sep 24 '19 at 16:59

16

There are a few instances where @gene's answer does not work.

For example, the using the overprecise value -73.92391000000001

geojson = {'type': 'Polygon', 'coordinates': [[[-73.92391, 41.31064], [-73.92391, 41.31069], [-73.92388, 41.31069], [-73.92388, 41.31064], [-73.92391, 41.31064]]]}
print(str(shape(geojson)))
POLYGON ((-73.92391000000001 41.31064, -73.92391000000001 41.31069, -73.92388 41.31069, -73.92388 41.31064, -73.92391000000001 41.31064))

EDIT (MAY 2019): I recommend a roundtrip using shapely.wkt.dumps method instead, where rounding_precision=n is the n'th decimal point to which you want the coordinates rounded:

shapely.wkt.loads(shapely.wkt.dumps(geom, rounding_precision=n))

edited Sep 24 '19 at 16:59

Jamie Bull

349
1
4
14

answered Mar 20 '18 at 17:45

koshmaster

636
4
10

1

And to remove sequential duplicates created by rounding: shapely.wkt.loads(...).simplify(0). – Ian Apr 30 '20 at 15:31
I like the simplicity of this answer. However what about performances ? – Astariul Mar 24 '21 at 02:48

score 2 · Answer 3 · edited Apr 30 '23 at 12:02

while there isn't a built in method or function in shapely you can use shapely.ops.transform to write a utility function to round coordinates for any geometry, including multipart geometries.

from shapley.ops import transform
def round_coordinates(geom, ndigits=2):
def _round_coords(x, y, z=None):
      x = round(x, ndigits)
      y = round(y, ndigits)
  if z is not None:
      z = round(z, ndigits) # corrected from x to z

  return [c for c in (x, y, z) if c is not None]


return transform(_round_coords, geom)

then running it like so

>>> _p = Polygon(((0.01, 1.234), (5.4321, 2.9348343943), (0.00000001, 9.23923), (0.01, 1.234)))
>>> round_coordinates(_p, 2).wkt
'POLYGON ((0.01 1.23, 5.43 2.93, 0 9.24, 0.01 1.23))'

imo this is the best answer since it uses the built in shapely transform method instead of converting to some other format. But it throws a weird error for me complaining about a generator. I fixed it by dropping the one line for loop. see here https://gis.stackexchange.com/a/445283 — Shawn, Nov 16 '22 at 14:12

Istopopoki · Answer 4 · 2019-12-12T14:16:53.867

In my case, the accepted answer did not work as my geojson was a multi polygon. Doing np.round(np.array(geosjon['coordinates']),2) throws a Python exception.

Instead, I wrote a short little recursive function that rounds all floats in a nested list. Meaning, that it also works on a list of list of list of list... for any depth.

def round_nested_list(l, precision):
    def round_element(e):
        if isinstance(e, list):
            return round_nested_list(e, precision)
        else:
            return round(e, precision)
    return [round_element(e) for e in l]

# test
l = [[[1.1123, 1.1123],[[1.1123, 1.1123],[[1.1123, 1.1123],[1.1123, 1.1123]]]], [[1.1123, 1.1123],[[1.1231, 1.1123],[[1.1231, 1.1231],[1.1231, 1.1231]]]]]

round_nested_list(l, 2)

# result is [[[1.11, 1.11], [[1.11, 1.11], [[1.11, 1.11], [1.11, 1.11]]]],
 [[1.11, 1.11], [[1.12, 1.11], [[1.12, 1.12], [1.12, 1.12]]]]]

I tried on a multipolygon with 20k+ points, and it took a few seconds. Good for me.

Edit : And if the coordinates are in tuples (at the last level), you can use :

def round_nested_list(l, precision):
    def round_element(e):
        if isinstance(e, list):
            return round_nested_list(e, precision)
        elif isinstance(e, tuple):
            return tuple(round(i, precision) for i in e)
        else:
            return round(e, precision)
    return [round_element(e) for e in l]

score 1 · Answer 5 · answered Apr 25 '23 at 13:03

As of shapely 2.0 (see release notes), a precision model is now available: the set_precision function and the grid_size parameter for overlay functions.

However, it's not necessarily behaving as you might expect, see some unexpected results with contains here.

score 0 · Answer 6 · answered Apr 08 '16 at 15:26

I don't see any way to reduce precision with Shapely even though GEOS, the base geometry engine, has precision reducer classes. You might try ogr2ogr. The geoJson driver has a "Coordinate Precision" option. You can export to geojson then create a shapefile from geojson. Just an idea.

Rounding all coordinates in shapely?

6 Answers6

Linked