From ccf9baca863ad56c042fa551a9adece437af2901 Mon Sep 17 00:00:00 2001 From: Justin Bronn Date: Thu, 14 Feb 2008 16:33:16 +0000 Subject: [PATCH] gis: geos: Added GeoJSON input/output support for `GEOSGeometry` (requires GDAL 1.5+); Polygons may now be instantiated with tuples of coordinates instead of only `LinearRing` instances; added `extent` property; added `coords` property alias for `tuple` (and added `tuple` for `GeometryCollection`); made small optimizations to KML and coordinate array generators. Small GDAL protochange that didn't make last commit: `get_envelope` is now generated w/the previously unused `env_func`. git-svn-id: http://code.djangoproject.com/svn/django/branches/gis@7114 bcc190cf-cafb-0310-a4f2-bffc1f526a37 --- django/contrib/gis/gdal/prototypes/geom.py | 6 +- django/contrib/gis/geos/base.py | 69 +++++++++++++----- django/contrib/gis/geos/collections.py | 10 ++- django/contrib/gis/geos/coordseq.py | 8 +-- django/contrib/gis/geos/geometries.py | 83 +++++++++++++--------- django/contrib/gis/tests/test_geos.py | 42 +++++++++-- 6 files changed, 150 insertions(+), 68 deletions(-) diff --git a/django/contrib/gis/gdal/prototypes/geom.py b/django/contrib/gis/gdal/prototypes/geom.py index 8a9b9ca82d..7757f557ad 100644 --- a/django/contrib/gis/gdal/prototypes/geom.py +++ b/django/contrib/gis/gdal/prototypes/geom.py @@ -105,7 +105,5 @@ geom_transform = void_output(lgdal.OGR_G_Transform, [c_void_p, c_void_p]) geom_transform_to = void_output(lgdal.OGR_G_TransformTo, [c_void_p, c_void_p]) # For retrieving the envelope of the geometry. -get_envelope = lgdal.OGR_G_GetEnvelope -get_envelope.restype = None -get_envelope.argtypes = [c_void_p, POINTER(OGREnvelope)] -get_envelope.errcheck = check_envelope +get_envelope = env_func(lgdal.OGR_G_GetEnvelope, [c_void_p, POINTER(OGREnvelope)]) + diff --git a/django/contrib/gis/geos/base.py b/django/contrib/gis/geos/base.py index 3c7ddc9a97..fdb04d0b78 100644 --- a/django/contrib/gis/geos/base.py +++ b/django/contrib/gis/geos/base.py @@ -5,7 +5,7 @@ # Python, ctypes and types dependencies. import re from ctypes import addressof, byref, c_double, c_size_t -from types import StringType, UnicodeType, IntType, FloatType, BufferType +from types import UnicodeType # GEOS-related dependencies. from django.contrib.gis.geos.coordseq import GEOSCoordSeq @@ -20,16 +20,17 @@ from django.contrib.gis.geos.prototypes import * # Trying to import GDAL libraries, if available. Have to place in # try/except since this package may be used outside GeoDjango. try: - from django.contrib.gis.gdal import OGRGeometry, SpatialReference - HAS_GDAL=True + from django.contrib.gis.gdal import OGRGeometry, SpatialReference, GEOJSON + HAS_GDAL = True except: - HAS_GDAL=False + HAS_GDAL, GEOJSON = False, False # Regular expression for recognizing HEXEWKB and WKT. A prophylactic measure # to prevent potentially malicious input from reaching the underlying C # library. Not a substitute for good web security programming practices. hex_regex = re.compile(r'^[0-9A-F]+$', re.I) wkt_regex = re.compile(r'^(SRID=(?P\d+);)?(?P(POINT|LINESTRING|LINEARRING|POLYGON|MULTIPOINT|MULTILINESTRING|MULTIPOLYGON|GEOMETRYCOLLECTION)[ACEGIMLONPSRUTY\d,\.\-\(\) ]+)$', re.I) +json_regex = re.compile(r'^\{.+\}$') class GEOSGeometry(object): "A class that, generally, encapsulates a GEOS geometry." @@ -50,24 +51,29 @@ class GEOSGeometry(object): The `srid` keyword is used to specify the Source Reference Identifier (SRID) number for this Geometry. If not set, the SRID will be None. """ - if isinstance(geo_input, UnicodeType): - # Encoding to ASCII, WKT or HEXEWKB doesn't need any more. - geo_input = geo_input.encode('ascii') - if isinstance(geo_input, StringType): - if hex_regex.match(geo_input): - # If the regex matches, the geometry is in HEX form. + if isinstance(geo_input, basestring): + if isinstance(geo_input, UnicodeType): + # Encoding to ASCII, WKT or HEXEWKB doesn't need any more. + geo_input = geo_input.encode('ascii') + + wkt_m = wkt_regex.match(geo_input) + if wkt_m: + # Handling WKT input. + if wkt_m.group('srid'): srid = int(wkt_m.group('srid')) + g = from_wkt(wkt_m.group('wkt')) + elif hex_regex.match(geo_input): + # Handling HEXEWKB input. g = from_hex(geo_input, len(geo_input)) + elif GEOJSON and json_regex.match(geo_input): + # Handling GeoJSON input. + wkb_input = str(OGRGeometry(geo_input).wkb) + g = from_wkb(wkb_input, len(wkb_input)) else: - m = wkt_regex.match(geo_input) - if m: - if m.group('srid'): srid = int(m.group('srid')) - g = from_wkt(m.group('wkt')) - else: - raise ValueError('String or unicode input unrecognized as WKT EWKT, and HEXEWKB.') + raise ValueError('String or unicode input unrecognized as WKT EWKT, and HEXEWKB.') elif isinstance(geo_input, GEOM_PTR): # When the input is a pointer to a geomtry (GEOM_PTR). g = geo_input - elif isinstance(geo_input, BufferType): + elif isinstance(geo_input, buffer): # When the input is a buffer (WKB). wkb_input = str(geo_input) g = from_wkb(wkb_input, len(wkb_input)) @@ -191,7 +197,8 @@ class GEOSGeometry(object): @property def coord_seq(self): "Returns a clone of the coordinate sequence for this Geometry." - return self._cs.clone() + if self.has_cs: + return self._cs.clone() #### Geometry Info #### @property @@ -307,7 +314,7 @@ class GEOSGeometry(object): Returns true if the elements in the DE-9IM intersection matrix for the two Geometries match the elements in pattern. """ - if not isinstance(pattern, StringType) or len(pattern) > 9: + if not isinstance(pattern, str) or len(pattern) > 9: raise GEOSException('invalid intersection matrix pattern') return geos_relatepattern(self.ptr, other.ptr, pattern) @@ -360,6 +367,15 @@ class GEOSGeometry(object): # str(self.wkb).encode('hex') return to_hex(self.ptr, byref(c_size_t())) + @property + def json(self): + """ + Returns GeoJSON representation of this Geometry if GDAL 1.5+ + is installed. + """ + if GEOJSON: return self.ogr.json + geojson = json + @property def wkb(self): "Returns the WKB of the Geometry as a buffer." @@ -521,6 +537,21 @@ class GEOSGeometry(object): raise TypeError('distance() works only on other GEOS Geometries.') return geos_distance(self.ptr, other.ptr, byref(c_double())) + @property + def extent(self): + """ + Returns the extent of this geometry as a 4-tuple, consisting of + (xmin, ymin, xmax, ymax). + """ + env = self.envelope + if isinstance(env, Point): + xmin, ymin = env.tuple + xmax, ymax = xmin, ymin + else: + xmin, ymin = env[0][0] + xmax, ymax = env[0][2] + return (xmin, ymin, xmax, ymax) + @property def length(self): """ diff --git a/django/contrib/gis/geos/collections.py b/django/contrib/gis/geos/collections.py index 7f27bf57fc..a69b2e7c84 100644 --- a/django/contrib/gis/geos/collections.py +++ b/django/contrib/gis/geos/collections.py @@ -85,9 +85,13 @@ class GeometryCollection(GEOSGeometry): @property def kml(self): "Returns the KML for this Geometry Collection." - kml = '' - for g in self: kml += g.kml - return kml + '' + return '%s' % ''.join([g.kml for g in self]) + + @property + def tuple(self): + "Returns a tuple of all the coordinates in this Geometry Collection" + return tuple([g.tuple for g in self]) + coords = tuple # MultiPoint, MultiLineString, and MultiPolygon class definitions. class MultiPoint(GeometryCollection): diff --git a/django/contrib/gis/geos/coordseq.py b/django/contrib/gis/geos/coordseq.py index 5d47985192..bc0c4794d4 100644 --- a/django/contrib/gis/geos/coordseq.py +++ b/django/contrib/gis/geos/coordseq.py @@ -153,14 +153,12 @@ class GEOSCoordSeq(object): # a Z dimension. if self.hasz: substr = '%s,%s,%s ' else: substr = '%s,%s,0 ' - kml = '' - for i in xrange(len(self)): - kml += substr % self[i] - return kml.strip() + '' + return '%s' % \ + ''.join([substr % self[i] for i in xrange(len(self))]).strip() @property def tuple(self): "Returns a tuple version of this coordinate sequence." n = self.size if n == 1: return self[0] - else: return tuple(self[i] for i in xrange(n)) + else: return tuple([self[i] for i in xrange(n)]) diff --git a/django/contrib/gis/geos/geometries.py b/django/contrib/gis/geos/geometries.py index 5a7ef5eabd..c5420e93af 100644 --- a/django/contrib/gis/geos/geometries.py +++ b/django/contrib/gis/geos/geometries.py @@ -4,7 +4,6 @@ GEOSGeometry. """ from ctypes import c_uint, byref -from types import FloatType, IntType, ListType, TupleType from django.contrib.gis.geos.base import GEOSGeometry from django.contrib.gis.geos.coordseq import GEOSCoordSeq from django.contrib.gis.geos.error import GEOSException, GEOSIndexError @@ -24,15 +23,15 @@ class Point(GEOSGeometry): >>> p = Point(5, 23, 8) # 3D point, passed in with individual parameters """ - if isinstance(x, (TupleType, ListType)): + if isinstance(x, (tuple, list)): # Here a tuple or list was passed in under the `x` parameter. ndim = len(x) if ndim < 2 or ndim > 3: raise TypeError('Invalid sequence parameter: %s' % str(x)) coords = x - elif isinstance(x, (IntType, FloatType)) and isinstance(y, (IntType, FloatType)): + elif isinstance(x, (int, float, long)) and isinstance(y, (int, float, long)): # Here X, Y, and (optionally) Z were passed in individually, as parameters. - if isinstance(z, (IntType, FloatType)): + if isinstance(z, (int, float, long)): ndim = 3 coords = [x, y, z] else: @@ -103,7 +102,7 @@ class Point(GEOSGeometry): # The tuple and coords properties tuple = property(get_coords, set_coords) - coords = property(get_coords, set_coords) + coords = tuple class LineString(GEOSGeometry): @@ -124,7 +123,7 @@ class LineString(GEOSGeometry): if len(args) == 1: coords = args[0] else: coords = args - if isinstance(coords, (TupleType, ListType)): + if isinstance(coords, (tuple, list)): # Getting the number of coords and the number of dimensions -- which # must stay the same, e.g., no LineString((1, 2), (1, 2, 3)). ncoords = len(coords) @@ -133,7 +132,7 @@ class LineString(GEOSGeometry): self._checkdim(ndim) # Incrementing through each of the coordinates and verifying for i in xrange(1, ncoords): - if not isinstance(coords[i], (TupleType, ListType, Point)): + if not isinstance(coords[i], (tuple, list, Point)): raise TypeError('each coordinate should be a sequence (list or tuple)') if len(coords[i]) != ndim: raise TypeError('Dimension mismatch.') numpy_coords = False @@ -193,6 +192,7 @@ class LineString(GEOSGeometry): def tuple(self): "Returns a tuple version of the geometry from the coordinate sequence." return self._cs.tuple + coords = tuple def _listarr(self, func): """ @@ -236,12 +236,17 @@ class Polygon(GEOSGeometry): def __init__(self, *args, **kwargs): """ Initializes on an exterior ring and a sequence of holes (both - instances of LinearRings. + instances may be either LinearRing instances, or a tuple/list + that may be constructed into a LinearRing). - Below are some examples of initialization, where shell, hole1, and - hole2 are valid LinearRing geometries: + Examples of initialization, where shell, hole1, and hole2 are + valid LinearRing geometries: >>> poly = Polygon(shell, hole1, hole2) >>> poly = Polygon(shell, (hole1, hole2)) + + Example where a tuple parameters are used: + >>> poly = Polygon(((0, 0), (0, 10), (10, 10), (0, 10), (0, 0)), + ((4, 4), (4, 6), (6, 6), (6, 4), (4, 4))) """ if not args: raise TypeError('Must provide at list one LinearRing instance to initialize Polygon.') @@ -249,32 +254,38 @@ class Polygon(GEOSGeometry): # Getting the ext_ring and init_holes parameters from the argument list ext_ring = args[0] init_holes = args[1:] - if len(init_holes) == 1 and isinstance(init_holes[0], (TupleType, ListType)): + n_holes = len(init_holes) + + # If initialized as Polygon(shell, (LinearRing, LinearRing)) [for backward-compatibility] + if n_holes == 1 and isinstance(init_holes[0], (tuple, list)) and \ + (len(init_holes[0]) == 0 or isinstance(init_holes[0][0], LinearRing)): init_holes = init_holes[0] + n_holes = len(init_holes) - # Ensuring the exterior ring parameter is a LinearRing object - if not isinstance(ext_ring, LinearRing): - raise TypeError('First argument for Polygon initialization must be a LinearRing.') + # Ensuring the exterior ring and holes parameters are LinearRing objects + # or may be instantiated into LinearRings. + ext_ring = self._construct_ring(ext_ring, 'Exterior parameter must be a LinearRing or an object that can initialize a LinearRing.') + holes_list = [] # Create new list, cause init_holes is a tuple. + for i in xrange(n_holes): + holes_list.append(self._construct_ring(init_holes[i], 'Holes parameter must be a sequence of LinearRings or objects that can initialize to LinearRings')) - # Making sure all of the holes are LinearRing objects - if False in [isinstance(hole, LinearRing) for hole in init_holes]: - raise TypeError('Holes parameter must be a sequence of LinearRings.') - - # Getting the holes array. - nholes = len(init_holes) - holes = get_pointer_arr(nholes) - for i in xrange(nholes): holes[i] = geom_clone(init_holes[i].ptr) + # Why another loop? Because if a TypeError is raised, cloned pointers will + # be around that can't be cleaned up. + holes = get_pointer_arr(n_holes) + for i in xrange(n_holes): holes[i] = geom_clone(holes_list[i].ptr) - # Getting the shell pointer address, + # Getting the shell pointer address. shell = geom_clone(ext_ring.ptr) # Calling with the GEOS createPolygon factory. - super(Polygon, self).__init__(create_polygon(shell, byref(holes), c_uint(nholes)), **kwargs) + super(Polygon, self).__init__(create_polygon(shell, byref(holes), c_uint(n_holes)), **kwargs) def __getitem__(self, index): """ - Returns the ring at the specified index. The first index, 0, will always - return the exterior ring. Indices > 0 will return the interior ring. + Returns the ring at the specified index. The first index, 0, will + always return the exterior ring. Indices > 0 will return the + interior ring at the given index (e.g., poly[1] and poly[2] would + return the first and second interior ring, respectively). """ if index == 0: return self.exterior_ring @@ -330,6 +341,15 @@ class Polygon(GEOSGeometry): if index < 0 or index >= len(self): raise GEOSIndexError('invalid Polygon ring index: %s' % index) + def _construct_ring(self, param, msg=''): + "Helper routine for trying to construct a ring from the given parameter." + if isinstance(param, LinearRing): return param + try: + ring = LinearRing(param) + return ring + except TypeError: + raise TypeError(msg) + def get_interior_ring(self, ring_i): """ Gets the interior ring at the specified index, 0 is for the first @@ -355,18 +375,17 @@ class Polygon(GEOSGeometry): # properties for the exterior ring/shell exterior_ring = property(get_ext_ring, set_ext_ring) - shell = property(get_ext_ring, set_ext_ring) + shell = exterior_ring @property def tuple(self): "Gets the tuple for each ring in this Polygon." - return tuple(self[i].tuple for i in xrange(len(self))) + return tuple([self[i].tuple for i in xrange(len(self))]) + coords = tuple @property def kml(self): "Returns the KML representation of this Polygon." - inner_kml = '' - if self.num_interior_rings > 0: - for i in xrange(self.num_interior_rings): - inner_kml += "%s" % self[i+1].kml + inner_kml = ''.join(["%s" % self[i+1].kml + for i in xrange(self.num_interior_rings)]) return "%s%s" % (self[0].kml, inner_kml) diff --git a/django/contrib/gis/tests/test_geos.py b/django/contrib/gis/tests/test_geos.py index 9b0d2b9ce1..19d28cef88 100644 --- a/django/contrib/gis/tests/test_geos.py +++ b/django/contrib/gis/tests/test_geos.py @@ -5,7 +5,7 @@ from django.contrib.gis.geos.base import HAS_GDAL from django.contrib.gis.tests.geometries import * if HAS_NUMPY: from numpy import array -if HAS_GDAL: from django.contrib.gis.gdal import OGRGeometry, SpatialReference, CoordTransform +if HAS_GDAL: from django.contrib.gis.gdal import OGRGeometry, SpatialReference, CoordTransform, GEOJSON class GEOSTest(unittest.TestCase): @@ -85,7 +85,16 @@ class GEOSTest(unittest.TestCase): self.assertEqual(srid, poly.shell.srid) self.assertEqual(srid, fromstr(poly.ewkt).srid) # Checking export - def test01i_eq(self): + def test01i_json(self): + "Testing GeoJSON input/output (via GDAL)." + if not HAS_GDAL or not GEOJSON: return + for g in json_geoms: + geom = GEOSGeometry(g.wkt) + self.assertEqual(g.json, geom.json) + self.assertEqual(g.json, geom.geojson) + self.assertEqual(GEOSGeometry(g.wkt), GEOSGeometry(geom.json)) + + def test01j_eq(self): "Testing equivalence with WKT." p = fromstr('POINT(5 23)') self.assertEqual(p, p.wkt) @@ -268,7 +277,7 @@ class GEOSTest(unittest.TestCase): # Testing __getitem__ and __setitem__ on invalid indices self.assertRaises(GEOSIndexError, poly.__getitem__, len(poly)) - #self.assertRaises(GEOSIndexError, poly.__setitem__, len(poly), False) + self.assertRaises(GEOSIndexError, poly.__setitem__, len(poly), False) self.assertRaises(GEOSIndexError, poly.__getitem__, -1) # Testing __iter__ @@ -279,8 +288,16 @@ class GEOSTest(unittest.TestCase): # Testing polygon construction. self.assertRaises(TypeError, Polygon.__init__, 0, [1, 2, 3]) self.assertRaises(TypeError, Polygon.__init__, 'foo') + + # Polygon(shell, (hole1, ... holeN)) rings = tuple(r for r in poly) self.assertEqual(poly, Polygon(rings[0], rings[1:])) + + # Polygon(shell_tuple, hole_tuple1, ... , hole_tupleN) + ring_tuples = tuple(r.tuple for r in poly) + self.assertEqual(poly, Polygon(*ring_tuples)) + + # Constructing with tuples of LinearRings. self.assertEqual(poly.wkt, Polygon(*tuple(r for r in poly)).wkt) self.assertEqual(poly.wkt, Polygon(*tuple(LinearRing(r.tuple) for r in poly)).wkt) @@ -686,12 +703,27 @@ class GEOSTest(unittest.TestCase): t2.transform(SpatialReference('EPSG:2774')) ct = CoordTransform(SpatialReference('WGS84'), SpatialReference(2774)) t3.transform(ct) - + prec = 3 for p in (t1, t2, t3): - prec = 3 self.assertAlmostEqual(trans.x, p.x, prec) self.assertAlmostEqual(trans.y, p.y, prec) + def test24_extent(self): + "Testing `extent` method." + # The xmin, ymin, xmax, ymax of the MultiPoint should be returned. + mp = MultiPoint(Point(5, 23), Point(0, 0), Point(10, 50)) + self.assertEqual((0.0, 0.0, 10.0, 50.0), mp.extent) + pnt = Point(5.23, 17.8) + # Extent of points is just the point itself repeated. + self.assertEqual((5.23, 17.8, 5.23, 17.8), pnt.extent) + # Testing on the 'real world' Polygon. + poly = fromstr(polygons[3].wkt) + ring = poly.shell + x, y = ring.x, ring.y + xmin, ymin = min(x), min(y) + xmax, ymax = max(x), max(y) + self.assertEqual((xmin, ymin, xmax, ymax), poly.extent) + def suite(): s = unittest.TestSuite() s.addTest(unittest.makeSuite(GEOSTest))