[ckan-dev] ckanext-spatial, asShape, Precision

[hotz]

Hello,

we have following problem with ckanext-spatial:

During harvesting the current master branch of Ckanext-spatial raises
the Validation error:

    raise p.toolkit.ValidationError(error_dict,
error_summary=package_error_summary(error_dict))
ValidationError: {'spatial': [u'Error: (InternalError) GEOSEquals:
TopologyException: side location conflict at 9.6314640334309907
53.394984999999998\nCONTEXT:  SQL function "st_equals" statement 1\n
\'SELECT ST_Equals(ST_GeomFromText(%(ST_GeomFromText_1)s,
%(ST_GeomFromText_2)s), ST_GeomFromEWKB(%(ST_GeomFromEWKB_1)s)) AS
"ST_Equals_1"\' {\'ST_GeomFromText_1\': \'MULTIPOLYGON
(((10.3951693282136 53.3936072531247, 10.3951693282136 53.7466021480132,
9.631464033430991 53.7466021480132, 9.631464033430991 53.3936072531247,
10.3951693282136 53.3936072531247)), ((10.326304 53.394985, 10.326304
53.964153, 8.420551 53.964153, 8.420551 53.394985, 10.326304
53.394985)))\', \'ST_GeomFromEWKB_1\': <read-only buffer for
0x7efe16241b10, size 199, offset 0 at 0x7efe1586aab0>,
\'ST_GeomFromText_2\': 4326}']}

Have a look at: 9.6314640334309907 vs. 9.631464033430991 and
53.394984999999998 vs. 53.394985


The process is as follows:

1) Harvesting of

http://www.metaver.de/csw?REQUEST=GetRecordById&SERVICE=CSW&VERSION=2.0.2&id=FDE79B82-201C-4A1A-A1E5-5A10E01E3701&iplug=/ingrid-group:ige-iplug-HH&elementSetName=full

provides the json geometry: {u'type': u'multipolygon', u'coordinates':
[[[[10.4255093212651, 53.3936072531247], [10.4255093212651,
53.7462455729101], [9.63146403343099, 53.7462455729101],
[9.63146403343099, 53.3936072531247], [10.4255093212651,
53.3936072531247]]], [[[10.326304, 53.394985], [10.326304, 53.964153],
[8.420551, 53.964153], [8.420551, 53.394985], [10.326304, 53.394985]]]]}

2)
This code in ckanext-spatial/ckanext/spatial/lib/__init__.py:
        shape = asShape(geometry)

creates:

shape: MULTIPOLYGON (((10.4255093212651 53.3936072531247,
10.4255093212651 53.7462455729101, 9.631464033430991 53.7462455729101,
9.631464033430991 53.3936072531247, 10.4255093212651 53.3936072531247)),
((10.326304 53.394985, 10.326304 53.964153, 8.420551 53.964153, 8.420551
53.394985, 10.326304 53.394985)))

which changes:
9.63146403343099

to:
9.631464033430991

(see also [1] for batch check)

3)

Current master of ckanext-spatial at github:
Test for equality in
/usr/lib/ckan/default/src/ckanext-spatial/ckanext/spatial/lib/__init__.py:

            if not compare_geometry_fields(package_extent.the_geom,
existing_package_extent.the_geom):

And compare_geometry_fields:
   
/usr/lib/ckan/default/src/ckanext-spatial/ckanext/spatial/geoalchemy_common.py
    return Session.scalar(ST_Equals(geom_field1, geom_field2))

And ST_Equals raises the error (see above).

4)

An earlier version of ckanext-spatial (Version 1.0.4) has no problems here:

Test for unequality withoug ST_Equals in
/usr/lib/ckan/default/src/ckanext-spatial/ckanext/spatial/lib/__init__.py :
Line 64:
     if Session.scalar(package_extent.the_geom.wkt) <>
Session.scalar(existing_package_extent.the_geom.wkt):

5)

Proposal (current quick fix):

In
/usr/lib/ckan/default/src/ckanext-spatial/ckanext/spatial/lib/__init__.py :

def save_package_extent...

+    try:
+               res = compare_geometry_fields(package_extent.the_geom,
existing_package_extent.the_geom)
+    except Exception,e:
+                #just unchange the extend by saying new extend and
existing extend from DB are equal!
+                log.debug('IGNORE NEW GEOMETRY: Found error in compare.
handle it as equal because new geometry seems to be erroneous. Package
%s' % package_id)
+                # because of erroness SQLAlchemy rollback so that next
updates will fit
+                Session.rollback()
+                res = True
-    #if not compare_geometry_fields(package_extent.the_geom,
existing_package_extent.the_geom):
+    if not res:
+       # Update extent
+       existing_package_extent.the_geom = package_extent.the_geom

6) asShape seems to have unsolved precision problems:

https://gis.stackexchange.com/questions/10201/arcpy-geometry-geo-interface-and-asshape-function-loss-of-precision-and-h

How do we cope with this?

Best wishes,

Lothar


[1]

from shapely.geometry import asShape
>>> geometry={u'type': u'multipolygon', u'coordinates':
[[[[10.4255093212651, 53.3936072531247], [10.4255093212651,
53.7462455729101], [9.63146403343099, 53.7462455729101],
[9.63146403343099, 53.3936072531247], [10.4255093212651,
53.3936072531247]]], [[[10.326304, 53.394985], [10.326304, 53.964153],
[8.420551, 53.964153], [8.420551, 53.394985], [10.326304, 53.394985]]]]}
>>> shape = asShape(geometry)
>>> str(shape)
'MULTIPOLYGON (((10.4255093212651 53.3936072531247, 10.4255093212651
53.7462455729101, 9.631464033430991 53.7462455729101, 9.631464033430991
53.3936072531247, 10.4255093212651 53.3936072531247)), ((10.326304
53.394985, 10.326304 53.964153, 8.420551 53.964153, 8.420551 53.394985,
10.326304 53.394985)))'
>>>

-- 
Hamburger Informatik Technologie-Center e.V. (HITeC)
c/o FB Informatik, Universität Hamburg
Dr. Lothar Hotz
Vogt-Kölln-Str. 30 , D-22527 Hamburg
Tel: 040/42883-2605/ Fax: -2572
http://www.hitec-hamburg.de
Private Seite: kogs-www.informatik.uni-hamburg.de/~hotz