1
0
mirror of https://github.com/django/django.git synced 2025-07-04 01:39:20 +00:00

gis: geos: Re-factored GEOS interface. Improvements include:

(1) All interactions with the GEOS library take place through 
     pre-defined ctypes prototypes, abstracting away return-value
     error-checking.
 (2) Mutability is now safe. Because GEOS geometry pointers are
     marked as `const` (Thanks Hobu), previous modification techniques 
     caused unpredictable instability when geometries were constructed 
     from the pointers of other geometries (e.g., a Polygon with rings 
     from another Polygon).  Geometry components are cloned first before 
     creation of the modified geometry.
 (3) Fixed memory leaks by freeing pointers from strings allocated
     in GEOS because ctypes only frees pointers allocated in Python.

Backwards-Incompatibility Notice:
 All children geometries (e.g., rings from a Polygon, geometries from a 
 collection) are now clones -- modifications will not propagate up to
 the parent geometry as before.


git-svn-id: http://code.djangoproject.com/svn/django/branches/gis@6653 bcc190cf-cafb-0310-a4f2-bffc1f526a37
This commit is contained in:
Justin Bronn 2007-11-06 00:39:36 +00:00
parent 2694fffcc0
commit 1c7ad200c7
16 changed files with 817 additions and 1050 deletions

View File

@ -1,34 +1,33 @@
"""
The goal of this module is to be a ctypes wrapper around the GEOS library
that will work on both *NIX and Windows systems. Specifically, this uses
the GEOS C api.
The goal of this module is to be a ctypes wrapper around the GEOS library
that will work on both *NIX and Windows systems. Specifically, this uses
the GEOS C api.
I have several motivations for doing this:
(1) The GEOS SWIG wrapper is no longer maintained, and requires the
installation of SWIG.
(2) The PCL implementation is over 2K+ lines of C and would make
PCL a requisite package for the GeoDjango application stack.
(3) Windows and Mac compatibility becomes substantially easier, and does not
require the additional compilation of PCL or GEOS and SWIG -- all that
is needed is a Win32 or Mac compiled GEOS C library (dll or dylib)
in a location that Python can read (e.g. 'C:\Python25').
I have several motivations for doing this:
(1) The GEOS SWIG wrapper is no longer maintained, and requires the
installation of SWIG.
(2) The PCL implementation is over 2K+ lines of C and would make
PCL a requisite package for the GeoDjango application stack.
(3) Windows and Mac compatibility becomes substantially easier, and does not
require the additional compilation of PCL or GEOS and SWIG -- all that
is needed is a Win32 or Mac compiled GEOS C library (dll or dylib)
in a location that Python can read (e.g. 'C:\Python25').
In summary, I wanted to wrap GEOS in a more maintainable and portable way using
only Python and the excellent ctypes library (now standard in Python 2.5).
In summary, I wanted to wrap GEOS in a more maintainable and portable way using
only Python and the excellent ctypes library (now standard in Python 2.5).
In the spirit of loose coupling, this library does not require Django or
GeoDjango. Only the GEOS C library and ctypes are needed for the platform
of your choice.
In the spirit of loose coupling, this library does not require Django or
GeoDjango. Only the GEOS C library and ctypes are needed for the platform
of your choice.
For more information about GEOS:
http://geos.refractions.net
For more information about GEOS:
http://geos.refractions.net
For more info about PCL and the discontinuation of the Python GEOS
library see Sean Gillies' writeup (and subsequent update) at:
http://zcologia.com/news/150/geometries-for-python/
http://zcologia.com/news/429/geometries-for-python-update/
For more info about PCL and the discontinuation of the Python GEOS
library see Sean Gillies' writeup (and subsequent update) at:
http://zcologia.com/news/150/geometries-for-python/
http://zcologia.com/news/429/geometries-for-python-update/
"""
from django.contrib.gis.geos.base import GEOSGeometry, wkt_regex, hex_regex
from django.contrib.gis.geos.geometries import Point, LineString, LinearRing, Polygon, HAS_NUMPY
from django.contrib.gis.geos.collections import GeometryCollection, MultiPoint, MultiLineString, MultiPolygon

View File

@ -1,19 +1,21 @@
"""
This module contains the 'base' GEOSGeometry object -- all GEOS geometries
inherit from this object.
This module contains the 'base' GEOSGeometry object -- all GEOS Geometries
inherit from this object.
"""
# ctypes and types dependencies.
from ctypes import \
byref, string_at, create_string_buffer, pointer, \
c_char_p, c_double, c_int, c_size_t
# 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
# Python and GEOS-related dependencies.
import re
from django.contrib.gis.geos.coordseq import GEOSCoordSeq, create_cs
# GEOS-related dependencies.
from django.contrib.gis.geos.coordseq import GEOSCoordSeq
from django.contrib.gis.geos.error import GEOSException, GEOSGeometryIndexError
from django.contrib.gis.geos.libgeos import lgeos, HAS_NUMPY
from django.contrib.gis.geos.pointer import GEOSPointer, NULL_GEOM
from django.contrib.gis.geos.libgeos import GEOM_PTR
# All other functions in this module come from the ctypes
# prototypes module -- which handles all interaction with
# the underlying GEOS library.
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.
@ -24,29 +26,29 @@ except:
HAS_GDAL=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.
# 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'^(POINT|LINESTRING|LINEARRING|POLYGON|MULTIPOINT|MULTILINESTRING|MULTIPOLYGON|GEOMETRYCOLLECTION)[ACEGIMLONPSRUTY\d,\.\-\(\) ]+$', re.I)
class GEOSGeometry(object):
"A class that, generally, encapsulates a GEOS geometry."
# Initially, all geometries use a NULL pointer.
_ptr = NULL_GEOM
# Initially, the geometry pointer is NULL
_ptr = None
#### Python 'magic' routines ####
def __init__(self, geo_input, srid=None):
"""
The base constructor for GEOS geometry objects, and may take the
following inputs:
following inputs:
* string: WKT
* string: HEXEWKB (a PostGIS-specific canonical form)
* buffer: WKB
The `srid` keyword is used to specify the Source Reference Identifier
(SRID) number for this Geometry. If not set, the SRID will be None.
(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.
@ -54,61 +56,53 @@ class GEOSGeometry(object):
if isinstance(geo_input, StringType):
if hex_regex.match(geo_input):
# If the regex matches, the geometry is in HEX form.
sz = c_size_t(len(geo_input))
buf = create_string_buffer(geo_input)
g = lgeos.GEOSGeomFromHEX_buf(buf, sz)
g = from_hex(geo_input, len(geo_input))
elif wkt_regex.match(geo_input):
# Otherwise, the geometry is in WKT form.
g = lgeos.GEOSGeomFromWKT(c_char_p(geo_input))
g = from_wkt(geo_input)
else:
raise GEOSException('String or unicode input unrecognized as WKT or HEXEWKB.')
elif isinstance(geo_input, (IntType, GEOSPointer)):
# When the input is either a memory address (an integer), or a
# GEOSPointer object.
raise ValueError('String or unicode input unrecognized as WKT or 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):
# When the input is a buffer (WKB).
wkb_input = str(geo_input)
sz = c_size_t(len(wkb_input))
g = lgeos.GEOSGeomFromWKB_buf(c_char_p(wkb_input), sz)
g = from_wkb(wkb_input, len(wkb_input))
else:
# Invalid geometry type.
raise TypeError, 'Improper geometry input type: %s' % str(type(geo_input))
raise TypeError('Improper geometry input type: %s' % str(type(geo_input)))
if bool(g):
# Setting the pointer object with a valid pointer.
self._ptr = GEOSPointer(g)
self._ptr = g
else:
raise GEOSException, 'Could not initialize GEOS Geometry with given input.'
raise GEOSException('Could not initialize GEOS Geometry with given input.')
# Setting the SRID, if given.
if srid and isinstance(srid, int): self.srid = srid
# Setting the class type (e.g., 'Point', 'Polygon', etc.)
self.__class__ = GEOS_CLASSES[self.geom_type]
# Setting the class type (e.g., Point, Polygon, etc.)
self.__class__ = GEOS_CLASSES[self.geom_typeid]
# Setting the coordinate sequence for the geometry (will be None on
# geometries that do not have coordinate sequences)
self._set_cs()
# _populate() needs to be called for parent Geometries.
if isinstance(self, (Polygon, GeometryCollection)): self._populate()
def __del__(self):
"""
Destroys this Geometry; in other words, frees the memory used by the
GEOS C++ object -- but only if the pointer is not a child Geometry
(e.g., don't delete the LinearRings spawned from a Polygon).
GEOS C++ object.
"""
#print 'base: Deleting %s (parent=%s, valid=%s)' % (self.__class__.__name__, self._ptr.parent, self._ptr.valid)
if not self._ptr.child: self._ptr.destroy()
if self._ptr: destroy_geom(self._ptr)
def __str__(self):
"WKT is used for the string representation."
return self.wkt
def __repr__(self):
return '<%s object>' % self.geom_type
"Short-hand representation because WKT may be very large."
return '<%s object at %s>' % (self.geom_type, hex(addressof(self._ptr)))
# Comparison operators
def __eq__(self, other):
@ -165,47 +159,6 @@ class GEOSGeometry(object):
"""
return self.sym_difference(other)
#### Internal GEOSPointer-related routines. ####
def _nullify(self):
"""
Returns the address of this Geometry, and nullifies any related pointers.
This function is called if this Geometry is used in the initialization
of another Geometry.
"""
# First getting the memory address of the geometry.
address = self._ptr()
# If the geometry is a child geometry, then the parent geometry pointer is
# nullified.
if self._ptr.child:
p = self._ptr.parent
# If we have a grandchild (a LinearRing from a MultiPolygon or
# GeometryCollection), then nullify the collection as well.
if p.child: p.parent.nullify()
p.nullify()
# Nullifying the geometry pointer
self._ptr.nullify()
return address
def _reassign(self, new_geom):
"Reassigns the internal pointer to that of the new Geometry."
# Only can re-assign when given a pointer or a geometry.
if not isinstance(new_geom, (GEOSPointer, GEOSGeometry)):
raise TypeError, 'cannot reassign geometry on given type: %s' % type(new_geom)
gtype = new_geom.geom_type
# Re-assigning the internal GEOSPointer to the new geometry, nullifying
# the new Geometry in the process.
if isinstance(new_geom, GEOSPointer): self._ptr = new_geom
else: self._ptr = GEOSPointer(new_geom._nullify())
# The new geometry class may be different from the original, so setting
# the __class__ and populating the internal geometry or ring dictionary.
self.__class__ = GEOS_CLASSES[gtype]
if isinstance(self, (Polygon, GeometryCollection)): self._populate()
#### Coordinate Sequence Routines ####
@property
def has_cs(self):
@ -219,41 +172,35 @@ class GEOSGeometry(object):
def _set_cs(self):
"Sets the coordinate sequence for this Geometry."
if self.has_cs:
if not self._ptr.coordseq_valid:
self._ptr.set_coordseq(lgeos.GEOSGeom_getCoordSeq(self._ptr()))
self._cs = GEOSCoordSeq(self._ptr, self.hasz)
self._cs = GEOSCoordSeq(get_cs(self._ptr), self.hasz)
else:
self._cs = None
@property
def coord_seq(self):
"Returns the coordinate sequence for this Geometry."
return self._cs
"Returns a clone of the coordinate sequence for this Geometry."
return self._cs.clone()
#### Geometry Info ####
@property
def geom_type(self):
"Returns a string representing the Geometry type, e.g. 'Polygon'"
return string_at(lgeos.GEOSGeomType(self._ptr()))
return geos_type(self._ptr)
@property
def geom_typeid(self):
"Returns an integer representing the Geometry type."
return lgeos.GEOSGeomTypeId(self._ptr())
return geos_typeid(self._ptr)
@property
def num_geom(self):
"Returns the number of geometries in the Geometry."
n = lgeos.GEOSGetNumGeometries(self._ptr())
if n == -1: raise GEOSException, 'Error getting number of geometries.'
else: return n
return get_num_geoms(self._ptr)
@property
def num_coords(self):
"Returns the number of coordinates in the Geometry."
n = lgeos.GEOSGetNumCoordinates(self._ptr())
if n == -1: raise GEOSException, 'Error getting number of coordinates.'
else: return n
return get_num_coords(self._ptr)
@property
def num_points(self):
@ -263,151 +210,126 @@ class GEOSGeometry(object):
@property
def dims(self):
"Returns the dimension of this Geometry (0=point, 1=line, 2=surface)."
return lgeos.GEOSGeom_getDimensions(self._ptr())
return get_dims(self._ptr)
def normalize(self):
"Converts this Geometry to normal form (or canonical form)."
status = lgeos.GEOSNormalize(self._ptr())
if status == -1: raise GEOSException, 'failed to normalize geometry'
## Internal for GEOS unary & binary predicate functions ##
def _unary_predicate(self, func):
"""
Returns the result, or raises an exception for the given unary predicate
function.
"""
val = func(self._ptr())
if val == 0: return False
elif val == 1: return True
else: raise GEOSException, '%s: exception occurred.' % func.__name__
def _binary_predicate(self, func, other, *args):
"""
Returns the result, or raises an exception for the given binary
predicate function.
"""
if not isinstance(other, GEOSGeometry):
raise TypeError, 'Binary predicate operation ("%s") requires another GEOSGeometry instance.' % func.__name__
val = func(self._ptr(), other._ptr(), *args)
if val == 0: return False
elif val == 1: return True
else: raise GEOSException, '%s: exception occurred.' % func.__name__
return geos_normalize(self._ptr)
#### Unary predicates ####
@property
def empty(self):
"""
Returns a boolean indicating whether the set of points in this Geometry
are empty.
are empty.
"""
return self._unary_predicate(lgeos.GEOSisEmpty)
@property
def valid(self):
"This property tests the validity of this Geometry."
return self._unary_predicate(lgeos.GEOSisValid)
@property
def simple(self):
"Returns false if the Geometry not simple."
return self._unary_predicate(lgeos.GEOSisSimple)
@property
def ring(self):
"Returns whether or not the geometry is a ring."
return self._unary_predicate(lgeos.GEOSisRing)
return geos_isempty(self._ptr)
@property
def hasz(self):
"Returns whether the geometry has a 3D dimension."
return self._unary_predicate(lgeos.GEOSHasZ)
return geos_hasz(self._ptr)
@property
def ring(self):
"Returns whether or not the geometry is a ring."
return geos_isring(self._ptr)
@property
def simple(self):
"Returns false if the Geometry not simple."
return geos_issimple(self._ptr)
@property
def valid(self):
"This property tests the validity of this Geometry."
return geos_isvalid(self._ptr)
#### Binary predicates. ####
def relate_pattern(self, other, pattern):
"""
Returns true if the elements in the DE-9IM intersection matrix for the
two Geometries match the elements in pattern.
"""
if len(pattern) > 9:
raise GEOSException, 'invalid intersection matrix pattern'
return self._binary_predicate(lgeos.GEOSRelatePattern, other, c_char_p(pattern))
def disjoint(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is FF*FF****.
"""
return self._binary_predicate(lgeos.GEOSDisjoint, other)
def touches(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is FT*******, F**T***** or F***T****.
"""
return self._binary_predicate(lgeos.GEOSTouches, other)
def intersects(self, other):
"Returns true if disjoint returns false."
return self._binary_predicate(lgeos.GEOSIntersects, other)
def contains(self, other):
"Returns true if other.within(this) returns true."
return geos_contains(self._ptr, other._ptr)
def crosses(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*T****** (for a point and a curve,a point and an area or a line and
an area) 0******** (for two curves).
is T*T****** (for a point and a curve,a point and an area or a line and
an area) 0******** (for two curves).
"""
return self._binary_predicate(lgeos.GEOSCrosses, other)
return geos_crosses(self._ptr, other._ptr)
def within(self, other):
def disjoint(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*F**F***.
is FF*FF****.
"""
return self._binary_predicate(lgeos.GEOSWithin, other)
def contains(self, other):
"Returns true if other.within(this) returns true."
return self._binary_predicate(lgeos.GEOSContains, other)
def overlaps(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*T***T** (for two points or two surfaces) 1*T***T** (for two curves).
"""
return self._binary_predicate(lgeos.GEOSOverlaps, other)
return geos_disjoint(self._ptr, other._ptr)
def equals(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*F**FFF*.
is T*F**FFF*.
"""
return self._binary_predicate(lgeos.GEOSEquals, other)
return geos_equals(self._ptr, other._ptr)
def equals_exact(self, other, tolerance=0):
"""
Returns true if the two Geometries are exactly equal, up to a
specified tolerance.
specified tolerance.
"""
return self._binary_predicate(lgeos.GEOSEqualsExact, other,
c_double(tolerance))
return geos_equalsexact(self._ptr, other._ptr, float(tolerance))
def intersects(self, other):
"Returns true if disjoint returns false."
return geos_intersects(self._ptr, other._ptr)
def overlaps(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*T***T** (for two points or two surfaces) 1*T***T** (for two curves).
"""
return geos_overlaps(self._ptr, other._ptr)
def relate_pattern(self, other, pattern):
"""
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:
raise GEOSException('invalid intersection matrix pattern')
return geos_relatepattern(self._ptr, other._ptr, pattern)
def touches(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is FT*******, F**T***** or F***T****.
"""
return geos_touches(self._ptr, other._ptr)
def within(self, other):
"""
Returns true if the DE-9IM intersection matrix for the two Geometries
is T*F**F***.
"""
return geos_within(self._ptr, other._ptr)
#### SRID Routines ####
def get_srid(self):
"Gets the SRID for the geometry, returns None if no SRID is set."
s = lgeos.GEOSGetSRID(self._ptr())
s = geos_get_srid(self._ptr)
if s == 0: return None
else: return s
def set_srid(self, srid):
"Sets the SRID for the geometry."
lgeos.GEOSSetSRID(self._ptr(), c_int(srid))
geos_set_srid(self._ptr, srid)
srid = property(get_srid, set_srid)
#### Output Routines ####
@property
def wkt(self):
"Returns the WKT (Well-Known Text) of the Geometry."
return string_at(lgeos.GEOSGeomToWKT(self._ptr()))
return to_wkt(self._ptr)
@property
def hex(self):
@ -418,16 +340,13 @@ class GEOSGeometry(object):
"""
# A possible faster, all-python, implementation:
# str(self.wkb).encode('hex')
sz = c_size_t()
h = lgeos.GEOSGeomToHEX_buf(self._ptr(), byref(sz))
return string_at(h, sz.value)
return to_hex(self._ptr, byref(c_size_t()))
@property
def wkb(self):
"Returns the WKB of the Geometry as a buffer."
sz = c_size_t()
h = lgeos.GEOSGeomToWKB_buf(self._ptr(), byref(sz))
return buffer(string_at(h, sz.value))
bin = to_wkb(self._ptr, byref(c_size_t()))
return buffer(bin)
@property
def kml(self):
@ -455,154 +374,138 @@ class GEOSGeometry(object):
else:
return None
#### Topology Routines ####
def _unary_topology(self, func, *args):
"""
Returns a GEOSGeometry for the given unary (takes only one Geomtery
as a paramter) topological operation.
"""
return GEOSGeometry(func(self._ptr(), *args), srid=self.srid)
@property
def crs(self):
"Alias for `srs` property."
return self.srs
def _binary_topology(self, func, other, *args):
"""
Returns a GEOSGeometry for the given binary (takes two Geometries
as parameters) topological operation.
"""
return GEOSGeometry(func(self._ptr(), other._ptr(), *args), srid=self.srid)
#### Topology Routines ####
def _topology(self, gptr):
"Helper routine to return Geometry from the given pointer."
return GEOSGeometry(gptr, srid=self.srid)
@property
def boundary(self):
"Returns the boundary as a newly allocated Geometry object."
return self._topology(geos_boundary(self._ptr))
def buffer(self, width, quadsegs=8):
"""
Returns a geometry that represents all points whose distance from this
Geometry is less than or equal to distance. Calculations are in the
Spatial Reference System of this Geometry. The optional third parameter sets
the number of segment used to approximate a quarter circle (defaults to 8).
(Text from PostGIS documentation at ch. 6.1.3)
Geometry is less than or equal to distance. Calculations are in the
Spatial Reference System of this Geometry. The optional third parameter sets
the number of segment used to approximate a quarter circle (defaults to 8).
(Text from PostGIS documentation at ch. 6.1.3)
"""
if not isinstance(width, (FloatType, IntType)):
raise TypeError, 'width parameter must be a float'
if not isinstance(quadsegs, IntType):
raise TypeError, 'quadsegs parameter must be an integer'
return self._unary_topology(lgeos.GEOSBuffer, c_double(width), c_int(quadsegs))
@property
def envelope(self):
"Return the envelope for this geometry (a polygon)."
return self._unary_topology(lgeos.GEOSEnvelope)
return self._topology(geos_buffer(self._ptr, width, quadsegs))
@property
def centroid(self):
"""
The centroid is equal to the centroid of the set of component Geometries
of highest dimension (since the lower-dimension geometries contribute zero
"weight" to the centroid).
of highest dimension (since the lower-dimension geometries contribute zero
"weight" to the centroid).
"""
return self._unary_topology(lgeos.GEOSGetCentroid)
@property
def boundary(self):
"Returns the boundary as a newly allocated Geometry object."
return self._unary_topology(lgeos.GEOSBoundary)
return self._topology(geos_centroid(self._ptr))
@property
def convex_hull(self):
"""
Returns the smallest convex Polygon that contains all the points
in the Geometry.
in the Geometry.
"""
return self._unary_topology(lgeos.GEOSConvexHull)
return self._topology(geos_convexhull(self._ptr))
def difference(self, other):
"""
Returns a Geometry representing the points making up this Geometry
that do not make up other.
"""
return self._topology(geos_difference(self._ptr, other._ptr))
@property
def envelope(self):
"Return the envelope for this geometry (a polygon)."
return self._topology(geos_envelope(self._ptr))
def intersection(self, other):
"Returns a Geometry representing the points shared by this Geometry and other."
return self._topology(geos_intersection(self._ptr, other._ptr))
@property
def point_on_surface(self):
"Computes an interior point of this Geometry."
return self._unary_topology(lgeos.GEOSPointOnSurface)
return self._topology(geos_pointonsurface(self._ptr))
def relate(self, other):
"Returns the DE-9IM intersection matrix for this Geometry and the other."
return geos_relate(self._ptr, other._ptr)
def simplify(self, tolerance=0.0, preserve_topology=False):
"""
Returns the Geometry, simplified using the Douglas-Peucker algorithm
to the specified tolerance (higher tolerance => less points). If no
tolerance provided, defaults to 0.
to the specified tolerance (higher tolerance => less points). If no
tolerance provided, defaults to 0.
By default, this function does not preserve topology - e.g. polygons can
be split, collapse to lines or disappear holes can be created or
disappear, and lines can cross. By specifying preserve_topology=True,
the result will have the same dimension and number of components as the
input. This is significantly slower.
be split, collapse to lines or disappear holes can be created or
disappear, and lines can cross. By specifying preserve_topology=True,
the result will have the same dimension and number of components as the
input. This is significantly slower.
"""
if preserve_topology:
return self._unary_topology(lgeos.GEOSTopologyPreserveSimplify, c_double(tolerance))
return self._topology(geos_preservesimplify(self._ptr, tolerance))
else:
return self._unary_topology(lgeos.GEOSSimplify, c_double(tolerance))
def relate(self, other):
"Returns the DE-9IM intersection matrix for this Geometry and the other."
return string_at(lgeos.GEOSRelate(self._ptr(), other._ptr()))
def difference(self, other):
"""Returns a Geometry representing the points making up this Geometry
that do not make up other."""
return self._binary_topology(lgeos.GEOSDifference, other)
return self._topology(geos_simplify(self._ptr, tolerance))
def sym_difference(self, other):
"""
Returns a set combining the points in this Geometry not in other,
and the points in other not in this Geometry.
and the points in other not in this Geometry.
"""
return self._binary_topology(lgeos.GEOSSymDifference, other)
def intersection(self, other):
"Returns a Geometry representing the points shared by this Geometry and other."
return self._binary_topology(lgeos.GEOSIntersection, other)
return self._topology(geos_symdifference(self._ptr, other._ptr))
def union(self, other):
"Returns a Geometry representing all the points in this Geometry and other."
return self._binary_topology(lgeos.GEOSUnion, other)
return self._topology(geos_union(self._ptr, other._ptr))
#### Other Routines ####
@property
def area(self):
"Returns the area of the Geometry."
a = c_double()
status = lgeos.GEOSArea(self._ptr(), byref(a))
if status != 1: return None
else: return a.value
return geos_area(self._ptr, byref(c_double()))
def distance(self, other):
"""
Returns the distance between the closest points on this Geometry
and the other. Units will be in those of the coordinate system of
the Geometry.
and the other. Units will be in those of the coordinate system of
the Geometry.
"""
if not isinstance(other, GEOSGeometry):
raise TypeError, 'distance() works only on other GEOS Geometries.'
dist = c_double()
status = lgeos.GEOSDistance(self._ptr(), other._ptr(), byref(dist))
if status != 1: return None
else: return dist.value
raise TypeError('distance() works only on other GEOS Geometries.')
return geos_distance(self._ptr, other._ptr, byref(c_double()))
@property
def length(self):
"""
Returns the length of this Geometry (e.g., 0 for point, or the
circumfrence of a Polygon).
circumfrence of a Polygon).
"""
l = c_double()
status = lgeos.GEOSLength(self._ptr(), byref(l))
if status != 1: return None
else: return l.value
return geos_length(self._ptr, byref(c_double()))
def clone(self):
"Clones this Geometry."
return GEOSGeometry(lgeos.GEOSGeom_clone(self._ptr()), srid=self.srid)
return GEOSGeometry(geom_clone(self._ptr), srid=self.srid)
# Class mapping dictionary
from django.contrib.gis.geos.geometries import Point, Polygon, LineString, LinearRing
from django.contrib.gis.geos.collections import GeometryCollection, MultiPoint, MultiLineString, MultiPolygon
GEOS_CLASSES = {'Point' : Point,
'Polygon' : Polygon,
'LineString' : LineString,
'LinearRing' : LinearRing,
'MultiPoint' : MultiPoint,
'MultiLineString' : MultiLineString,
'MultiPolygon' : MultiPolygon,
'GeometryCollection' : GeometryCollection,
GEOS_CLASSES = {0 : Point,
1 : LineString,
2 : LinearRing,
3 : Polygon,
4 : MultiPoint,
5 : MultiLineString,
6 : MultiPolygon,
7 : GeometryCollection,
}

View File

@ -1,14 +1,14 @@
"""
This module houses the Geometry Collection objects:
GeometryCollection, MultiPoint, MultiLineString, and MultiPolygon
This module houses the Geometry Collection objects:
GeometryCollection, MultiPoint, MultiLineString, and MultiPolygon
"""
from ctypes import c_int, c_uint, byref, cast
from ctypes import c_int, c_uint, byref
from types import TupleType, ListType
from django.contrib.gis.geos.base import GEOSGeometry
from django.contrib.gis.geos.error import GEOSException, GEOSGeometryIndexError
from django.contrib.gis.geos.geometries import Point, LineString, LinearRing, Polygon
from django.contrib.gis.geos.libgeos import lgeos, get_pointer_arr, GEOM_PTR
from django.contrib.gis.geos.pointer import GEOSPointer
from django.contrib.gis.geos.libgeos import get_pointer_arr, GEOM_PTR
from django.contrib.gis.geos.prototypes import create_collection, destroy_geom, geom_clone, geos_typeid, get_cs, get_geomn
class GeometryCollection(GEOSGeometry):
_allowed = (Point, LineString, LinearRing, Polygon)
@ -33,66 +33,40 @@ class GeometryCollection(GEOSGeometry):
# Ensuring that only the permitted geometries are allowed in this collection
if False in [isinstance(geom, self._allowed) for geom in init_geoms]:
raise TypeError, 'Invalid Geometry type encountered in the arguments.'
raise TypeError('Invalid Geometry type encountered in the arguments.')
# Creating the geometry pointer array, and populating each element in
# the array with the address of the Geometry returned by _nullify().
ngeom = len(init_geoms)
geoms = get_pointer_arr(ngeom)
for i in xrange(ngeom):
geoms[i] = cast(init_geoms[i]._nullify(), GEOM_PTR)
# Calling the parent class, using the pointer returned from the
# GEOS createCollection() factory.
addr = lgeos.GEOSGeom_createCollection(c_int(self._typeid),
byref(geoms), c_uint(ngeom))
super(GeometryCollection, self).__init__(addr, **kwargs)
def __del__(self):
"Overloaded deletion method for Geometry Collections."
#print 'collection: Deleting %s (parent=%s, valid=%s)' % (self.__class__.__name__, self._ptr.parent, self._ptr.valid)
# If this geometry is still valid, it hasn't been modified by others.
if self._ptr.valid:
# Nullifying pointers to internal Geometries, preventing any
# attempted future access.
for g in self._ptr: g.nullify()
else:
# Internal memory has become part of other Geometry objects; must
# delete the internal objects which are still valid individually,
# because calling the destructor on the entire geometry will result
# in an attempted deletion of NULL pointers for the missing
# components (which may crash Python).
for g in self._ptr:
if len(g) > 0:
# The collection geometry is a Polygon, destroy any leftover
# LinearRings.
for r in g: r.destroy()
g.destroy()
super(GeometryCollection, self).__del__()
# Creating the geometry pointer array.
ngeoms = len(init_geoms)
geoms = get_pointer_arr(ngeoms)
for i in xrange(ngeoms): geoms[i] = geom_clone(init_geoms[i]._ptr)
super(GeometryCollection, self).__init__(create_collection(c_int(self._typeid), byref(geoms), c_uint(ngeoms)), **kwargs)
def __getitem__(self, index):
"Returns the Geometry from this Collection at the given index (0-based)."
# Checking the index and returning the corresponding GEOS geometry.
self._checkindex(index)
return GEOSGeometry(self._ptr[index], srid=self.srid)
return GEOSGeometry(geom_clone(get_geomn(self._ptr, index)), srid=self.srid)
def __setitem__(self, index, geom):
"Sets the Geometry at the specified index."
self._checkindex(index)
if not isinstance(geom, self._allowed):
raise TypeError, 'Incompatible Geometry for collection.'
raise TypeError('Incompatible Geometry for collection.')
# Constructing the list of geometries that will go in the collection.
new_geoms = []
for i in xrange(len(self)):
if i == index: new_geoms.append(geom)
else: new_geoms.append(self[i])
ngeoms = len(self)
geoms = get_pointer_arr(ngeoms)
for i in xrange(ngeoms):
if i == index:
geoms[i] = geom_clone(geom._ptr)
else:
geoms[i] = geom_clone(get_geomn(self._ptr, i))
# Creating a new geometry collection from the list, and
# re-assigning the pointers.
new_collection = self.__class__(*new_geoms, **{'srid':self.srid})
self._reassign(new_collection)
# Creating a new collection, and destroying the contents of the previous poiner.
prev_ptr = self._ptr
srid = self.srid
self._ptr = create_collection(c_int(self._typeid), byref(geoms), c_uint(ngeoms))
if srid: self.srid = srid
destroy_geom(prev_ptr)
def __iter__(self):
"Iterates over each Geometry in the Collection."
@ -106,28 +80,7 @@ class GeometryCollection(GEOSGeometry):
def _checkindex(self, index):
"Checks the given geometry index."
if index < 0 or index >= self.num_geom:
raise GEOSGeometryIndexError, 'invalid GEOS Geometry index: %s' % str(index)
def _nullify(self):
"Overloaded from base method to nullify geometry references in this Collection."
# Nullifying the references to the internal Geometry objects from this Collection.
for g in self._ptr: g.nullify()
return super(GeometryCollection, self)._nullify()
def _populate(self):
"Internal routine that populates the internal children geometries list."
ptr_list = []
for i in xrange(len(self)):
# Getting the geometry pointer for the geometry at the index.
geom_ptr = lgeos.GEOSGetGeometryN(self._ptr(), c_int(i))
# Adding the coordinate sequence to the list, or using None if the
# collection Geometry doesn't support coordinate sequences.
if lgeos.GEOSGeomTypeId(geom_ptr) in (0, 1, 2):
ptr_list.append((geom_ptr, lgeos.GEOSGeom_getCoordSeq(geom_ptr)))
else:
ptr_list.append((geom_ptr, None))
self._ptr.set_children(ptr_list)
raise GEOSGeometryIndexError('invalid GEOS Geometry index: %s' % str(index))
@property
def kml(self):

View File

@ -1,13 +1,13 @@
"""
This module houses the GEOSCoordSeq object, and is used internally
by GEOSGeometry to house the actual coordinates of the Point,
LineString, and LinearRing geometries.
This module houses the GEOSCoordSeq object, which is used internally
by GEOSGeometry to house the actual coordinates of the Point,
LineString, and LinearRing geometries.
"""
from django.contrib.gis.geos.error import GEOSException, GEOSGeometryIndexError
from django.contrib.gis.geos.libgeos import lgeos, HAS_NUMPY
from django.contrib.gis.geos.pointer import GEOSPointer
from ctypes import c_double, c_int, c_uint, byref
from ctypes import c_double, c_uint, byref
from types import ListType, TupleType
from django.contrib.gis.geos.error import GEOSException, GEOSGeometryIndexError
from django.contrib.gis.geos.libgeos import CS_PTR, HAS_NUMPY
from django.contrib.gis.geos.prototypes import cs_clone, cs_getdims, cs_getordinate, cs_getsize, cs_setordinate
if HAS_NUMPY: from numpy import ndarray
class GEOSCoordSeq(object):
@ -16,15 +16,15 @@ class GEOSCoordSeq(object):
#### Python 'magic' routines ####
def __init__(self, ptr, z=False):
"Initializes from a GEOS pointer."
if not isinstance(ptr, GEOSPointer):
raise TypeError, 'Coordinate sequence should initialize with a GEOSPointer.'
if not isinstance(ptr, CS_PTR):
raise TypeError('Coordinate sequence should initialize with a CS_PTR.')
self._ptr = ptr
self._z = z
def __iter__(self):
"Iterates over each point in the coordinate sequence."
for i in xrange(self.size):
yield self.__getitem__(i)
yield self[i]
def __len__(self):
"Returns the number of points in the coordinate sequence."
@ -49,7 +49,7 @@ class GEOSCoordSeq(object):
elif HAS_NUMPY and isinstance(value, ndarray):
pass
else:
raise TypeError, 'Must set coordinate with a sequence (list, tuple, or numpy array).'
raise TypeError('Must set coordinate with a sequence (list, tuple, or numpy array).')
# Checking the dims of the input
if self.dims == 3 and self._z:
n_args = 3
@ -58,7 +58,7 @@ class GEOSCoordSeq(object):
n_args = 2
set_3d = False
if len(value) != n_args:
raise TypeError, 'Dimension of value does not match.'
raise TypeError('Dimension of value does not match.')
# Setting the X, Y, Z
self.setX(index, value[0])
self.setY(index, value[1])
@ -69,43 +69,25 @@ class GEOSCoordSeq(object):
"Checks the given index."
sz = self.size
if (sz < 1) or (index < 0) or (index >= sz):
raise GEOSGeometryIndexError, 'invalid GEOS Geometry index: %s' % str(index)
raise GEOSGeometryIndexError('invalid GEOS Geometry index: %s' % str(index))
def _checkdim(self, dim):
"Checks the given dimension."
if dim < 0 or dim > 2:
raise GEOSException, 'invalid ordinate dimension "%d"' % dim
raise GEOSException('invalid ordinate dimension "%d"' % dim)
#### Ordinate getting and setting routines ####
def getOrdinate(self, dimension, index):
"Returns the value for the given dimension and index."
self._checkindex(index)
self._checkdim(dimension)
# Wrapping the dimension, index
dim = c_uint(dimension)
idx = c_uint(index)
# 'd' is the value of the ordinate, passed in by reference
d = c_double()
status = lgeos.GEOSCoordSeq_getOrdinate(self._ptr.coordseq(), idx, dim, byref(d))
if status == 0:
raise GEOSException, 'could not retrieve %sth ordinate at index: %s' % (dimension, index)
return d.value
return cs_getordinate(self._ptr, index, dimension, byref(c_double()))
def setOrdinate(self, dimension, index, value):
"Sets the value for the given dimension and index."
self._checkindex(index)
self._checkdim(dimension)
# Wrapping the dimension, index
dim = c_uint(dimension)
idx = c_uint(index)
# Setting the ordinate
status = lgeos.GEOSCoordSeq_setOrdinate(self._ptr.coordseq(), idx, dim, c_double(value))
if status == 0:
raise GEOSException, 'Could not set the ordinate for (dim, index): (%d, %d)' % (dimension, index)
cs_setordinate(self._ptr, index, dimension, value)
def getX(self, index):
"Get the X value at the index."
@ -135,34 +117,25 @@ class GEOSCoordSeq(object):
@property
def size(self):
"Returns the size of this coordinate sequence."
n = c_uint(0)
status = lgeos.GEOSCoordSeq_getSize(self._ptr.coordseq(), byref(n))
if status == 0:
raise GEOSException, 'Could not get CoordSeq size.'
return n.value
return cs_getsize(self._ptr, byref(c_uint()))
@property
def dims(self):
"Returns the dimensions of this coordinate sequence."
n = c_uint(0)
status = lgeos.GEOSCoordSeq_getDimensions(self._ptr.coordseq(), byref(n))
if status == 0:
raise GEOSException, 'Could not get CoordSeq dimensions.'
return n.value
return cs_getdims(self._ptr, byref(c_uint()))
@property
def hasz(self):
"""
Returns whether this coordinate sequence is 3D. This property value is
inherited from the parent Geometry.
inherited from the parent Geometry.
"""
return self._z
### Other Methods ###
@property
def clone(self):
"Clones this coordinate sequence."
return GEOSCoordSeq(GEOSPointer(0, lgeos.GEOSCoordSeq_clone(self._ptr.coordseq())), self.hasz)
return GEOSCoordSeq(cs_clone(self._ptr), self.hasz)
@property
def kml(self):
@ -180,11 +153,5 @@ class GEOSCoordSeq(object):
def tuple(self):
"Returns a tuple version of this coordinate sequence."
n = self.size
if n == 1:
return self.__getitem__(0)
else:
return tuple(self.__getitem__(i) for i in xrange(n))
# ctypes prototype for the Coordinate Sequence creation factory
create_cs = lgeos.GEOSCoordSeq_create
create_cs.argtypes = [c_uint, c_uint]
if n == 1: return self[0]
else: return tuple(self[i] for i in xrange(n))

View File

@ -1,6 +1,6 @@
"""
This module houses the GEOS exceptions, specifically, GEOSException and
GEOSGeometryIndexError.
This module houses the GEOS exceptions, specifically, GEOSException and
GEOSGeometryIndexError.
"""
class GEOSException(Exception):

View File

@ -1,15 +1,15 @@
"""
This module houses the Point, LineString, LinearRing, and Polygon OGC
geometry classes. All geometry classes in this module inherit from
GEOSGeometry.
This module houses the Point, LineString, LinearRing, and Polygon OGC
geometry classes. All geometry classes in this module inherit from
GEOSGeometry.
"""
from ctypes import c_double, c_int, c_uint, byref, cast
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, create_cs
from django.contrib.gis.geos.libgeos import lgeos, get_pointer_arr, GEOM_PTR, HAS_NUMPY
from django.contrib.gis.geos.pointer import GEOSPointer
from django.contrib.gis.geos.coordseq import GEOSCoordSeq
from django.contrib.gis.geos.error import GEOSException, GEOSGeometryIndexError
from django.contrib.gis.geos.libgeos import get_pointer_arr, GEOM_PTR, HAS_NUMPY
from django.contrib.gis.geos.prototypes import *
if HAS_NUMPY: from numpy import ndarray, array
class Point(GEOSGeometry):
@ -17,16 +17,18 @@ class Point(GEOSGeometry):
def __init__(self, x, y=None, z=None, srid=None):
"""
The Point object may be initialized with either a tuple, or individual
parameters. For example:
>>> p = Point((5, 23)) # 2D point, passed in as a tuple
>>> p = Point(5, 23, 8) # 3D point, passed in with individual parameters
parameters.
For Example:
>>> p = Point((5, 23)) # 2D point, passed in as a tuple
>>> p = Point(5, 23, 8) # 3D point, passed in with individual parameters
"""
if isinstance(x, (TupleType, ListType)):
# 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)
raise TypeError('Invalid sequence parameter: %s' % str(x))
coords = x
elif isinstance(x, (IntType, FloatType)) and isinstance(y, (IntType, FloatType)):
# Here X, Y, and (optionally) Z were passed in individually, as parameters.
@ -37,20 +39,17 @@ class Point(GEOSGeometry):
ndim = 2
coords = [x, y]
else:
raise TypeError, 'Invalid parameters given for Point initialization.'
raise TypeError('Invalid parameters given for Point initialization.')
# Creating the coordinate sequence, and setting X, Y, [Z]
cs = create_cs(c_uint(1), c_uint(ndim))
status = lgeos.GEOSCoordSeq_setX(cs, c_uint(0), c_double(coords[0]))
if not status: raise GEOSException, 'Could not set X during Point initialization.'
status = lgeos.GEOSCoordSeq_setY(cs, c_uint(0), c_double(coords[1]))
if not status: raise GEOSException, 'Could not set Y during Point initialization.'
if ndim == 3:
status = lgeos.GEOSCoordSeq_setZ(cs, c_uint(0), c_double(coords[2]))
cs_setx(cs, 0, coords[0])
cs_sety(cs, 0, coords[1])
if ndim == 3: cs_setz(cs, 0, coords[2])
# Initializing using the address returned from the GEOS
# createPoint factory.
super(Point, self).__init__(lgeos.GEOSGeom_createPoint(cs), srid=srid)
super(Point, self).__init__(create_point(cs), srid=srid)
def __len__(self):
"Returns the number of dimensions for this Point (either 0, 2 or 3)."
@ -58,43 +57,35 @@ class Point(GEOSGeometry):
if self.hasz: return 3
else: return 2
def _getOrdinate(self, dim, idx):
"The coordinate sequence getOrdinate() wrapper."
return self._cs.getOrdinate(dim, idx)
def _setOrdinate(self, dim, idx, value):
"The coordinate sequence setOrdinate() wrapper."
self._cs.setOrdinate(dim, idx, value)
def get_x(self):
"Returns the X component of the Point."
return self._getOrdinate(0, 0)
return self._cs.getOrdinate(0, 0)
def set_x(self, value):
"Sets the X component of the Point."
self._setOrdinate(0, 0, value)
self._cs.setOrdinate(0, 0, value)
def get_y(self):
"Returns the Y component of the Point."
return self._getOrdinate(1, 0)
return self._cs.getOrdinate(1, 0)
def set_y(self, value):
"Sets the Y component of the Point."
self._setOrdinate(1, 0, value)
self._cs.setOrdinate(1, 0, value)
def get_z(self):
"Returns the Z component of the Point."
if self.hasz:
return self._getOrdinate(2, 0)
return self._cs.getOrdinate(2, 0)
else:
return None
def set_z(self, value):
"Sets the Z component of the Point."
if self.hasz:
self._setOrdinate(2, 0, value)
self._cs.setOrdinate(2, 0, value)
else:
raise GEOSException, 'Cannot set Z on 2D Point.'
raise GEOSException('Cannot set Z on 2D Point.')
# X, Y, Z properties
x = property(get_x, set_x)
@ -120,14 +111,14 @@ class LineString(GEOSGeometry):
def __init__(self, *args, **kwargs):
"""
Initializes on the given sequence -- may take lists, tuples, NumPy arrays
of X,Y pairs, or Point objects. If Point objects are used, ownership is
_not_ transferred to the LineString object.
of X,Y pairs, or Point objects. If Point objects are used, ownership is
_not_ transferred to the LineString object.
Examples:
ls = LineString((1, 1), (2, 2))
ls = LineString([(1, 1), (2, 2)])
ls = LineString(array([(1, 1), (2, 2)]))
ls = LineString(Point(1, 1), Point(2, 2))
ls = LineString((1, 1), (2, 2))
ls = LineString([(1, 1), (2, 2)])
ls = LineString(array([(1, 1), (2, 2)]))
ls = LineString(Point(1, 1), Point(2, 2))
"""
# If only one argument provided, set the coords array appropriately
if len(args) == 1: coords = args[0]
@ -138,44 +129,44 @@ class LineString(GEOSGeometry):
# must stay the same, e.g., no LineString((1, 2), (1, 2, 3)).
ncoords = len(coords)
if coords: ndim = len(coords[0])
else: raise TypeError, 'Cannot initialize on empty sequence.'
else: raise TypeError('Cannot initialize on empty sequence.')
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)):
raise TypeError, 'each coordinate should be a sequence (list or tuple)'
if len(coords[i]) != ndim: raise TypeError, 'Dimension mismatch.'
raise TypeError('each coordinate should be a sequence (list or tuple)')
if len(coords[i]) != ndim: raise TypeError('Dimension mismatch.')
numpy_coords = False
elif HAS_NUMPY and isinstance(coords, ndarray):
shape = coords.shape # Using numpy's shape.
if len(shape) != 2: raise TypeError, 'Too many dimensions.'
if len(shape) != 2: raise TypeError('Too many dimensions.')
self._checkdim(shape[1])
ncoords = shape[0]
ndim = shape[1]
numpy_coords = True
else:
raise TypeError, 'Invalid initialization input for LineStrings.'
raise TypeError('Invalid initialization input for LineStrings.')
# Creating a coordinate sequence object because it is easier to
# set the points using GEOSCoordSeq.__setitem__().
cs = GEOSCoordSeq(GEOSPointer(0, create_cs(c_uint(ncoords), c_uint(ndim))), z=bool(ndim==3))
cs = GEOSCoordSeq(create_cs(ncoords, ndim), z=bool(ndim==3))
for i in xrange(ncoords):
if numpy_coords: cs[i] = coords[i,:]
elif isinstance(coords[i], Point): cs[i] = coords[i].tuple
else: cs[i] = coords[i]
# Getting the initialization function
# Getting the correct initialization function
if kwargs.get('ring', False):
func = lgeos.GEOSGeom_createLinearRing
func = create_linearring
else:
func = lgeos.GEOSGeom_createLineString
func = create_linestring
# If SRID was passed in with the keyword arguments
srid = kwargs.get('srid', None)
# Calling the base geometry initialization with the returned pointer
# from the function.
super(LineString, self).__init__(func(cs._ptr.coordseq()), srid=srid)
super(LineString, self).__init__(func(cs._ptr), srid=srid)
def __getitem__(self, index):
"Gets the point at the specified index."
@ -187,15 +178,15 @@ class LineString(GEOSGeometry):
def __iter__(self):
"Allows iteration over this LineString."
for i in xrange(self.__len__()):
yield self.__getitem__(i)
for i in xrange(len(self)):
yield self[i]
def __len__(self):
"Returns the number of points in this LineString."
return len(self._cs)
def _checkdim(self, dim):
if dim not in (2, 3): raise TypeError, 'Dimension mismatch.'
if dim not in (2, 3): raise TypeError('Dimension mismatch.')
#### Sequence Properties ####
@property
@ -245,16 +236,15 @@ class Polygon(GEOSGeometry):
def __init__(self, *args, **kwargs):
"""
Initializes on an exterior ring and a sequence of holes (both
instances of LinearRings. All LinearRing instances used for creation
will become owned by this Polygon.
instances of LinearRings.
Below are some examples of initialization, where shell, hole1, and
hole2 are valid LinearRing geometries:
>>> poly = Polygon(shell, hole1, hole2)
>>> poly = Polygon(shell, (hole1, hole2))
hole2 are valid LinearRing geometries:
>>> poly = Polygon(shell, hole1, hole2)
>>> poly = Polygon(shell, (hole1, hole2))
"""
if not args:
raise TypeError, 'Must provide at list one LinearRing instance to initialize Polygon.'
raise TypeError('Must provide at list one LinearRing instance to initialize Polygon.')
# Getting the ext_ring and init_holes parameters from the argument list
ext_ring = args[0]
@ -264,43 +254,22 @@ class Polygon(GEOSGeometry):
# 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.'
raise TypeError('First argument for Polygon initialization must be a LinearRing.')
# 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.'
raise TypeError('Holes parameter must be a sequence of LinearRings.')
# Getting the holes
# Getting the holes array.
nholes = len(init_holes)
holes = get_pointer_arr(nholes)
for i in xrange(nholes):
# Casting to the Geometry Pointer type
holes[i] = cast(init_holes[i]._nullify(), GEOM_PTR)
for i in xrange(nholes): holes[i] = geom_clone(init_holes[i]._ptr)
# Getting the shell pointer address,
shell = ext_ring._nullify()
shell = geom_clone(ext_ring._ptr)
# Calling with the GEOS createPolygon factory.
super(Polygon, self).__init__(lgeos.GEOSGeom_createPolygon(shell, byref(holes), c_uint(nholes)), **kwargs)
def __del__(self):
"Overloaded deletion method for Polygons."
#print 'polygon: Deleting %s (parent=%s, valid=%s)' % (self.__class__.__name__, self._ptr.parent, self._ptr.valid)
# Not performed on children Polygons from MultiPolygon or GeometryCollection objects.
if not self._ptr.child:
# If this geometry is still valid, it hasn't been modified by others.
if self._ptr.valid:
# Nulling the pointers to internal rings, preventing any
# attempted future access.
for r in self._ptr: r.nullify()
else:
# Internal memory has become part of other Geometry objects; must
# delete the internal objects which are still valid individually,
# because calling the destructor on entire geometry will result
# in an attempted deletion of NULL pointers for the missing
# components (which may crash Python).
for r in self._ptr: r.destroy()
super(Polygon, self).__del__()
super(Polygon, self).__init__(create_polygon(shell, byref(holes), c_uint(nholes)), **kwargs)
def __getitem__(self, index):
"""
@ -318,22 +287,39 @@ class Polygon(GEOSGeometry):
# Checking the index and ring parameters.
self._checkindex(index)
if not isinstance(ring, LinearRing):
raise TypeError, 'must set Polygon index with a LinearRing object'
raise TypeError('must set Polygon index with a LinearRing object')
# Constructing the ring parameters
new_rings = []
for i in xrange(len(self)):
if index == i: new_rings.append(ring)
else: new_rings.append(self[i])
# Getting the shell
if index == 0:
shell = geom_clone(ring._ptr)
else:
shell = geom_clone(get_extring(self._ptr))
# Constructing the new Polygon with the ring parameters, and reassigning the internals.
new_poly = Polygon(*new_rings, **{'srid':self.srid})
self._reassign(new_poly)
# Getting the interior rings (holes)
nholes = len(self)-1
if nholes > 0:
holes = get_pointer_arr(nholes)
for i in xrange(nholes):
if i == (index-1):
holes[i] = geom_clone(ring._ptr)
else:
holes[i] = geom_clone(get_intring(self._ptr, i))
holes_param = byref(holes)
else:
holes_param = None
# Getting the current pointer, replacing with the newly constructed
# geometry, and destroying the old geometry.
prev_ptr = self._ptr
srid = self.srid
self._ptr = create_polygon(shell, holes_param, c_uint(nholes))
if srid: self.srid = srid
destroy_geom(prev_ptr)
def __iter__(self):
"Iterates over each ring in the polygon."
for i in xrange(len(self)):
yield self.__getitem__(i)
yield self[i]
def __len__(self):
"Returns the number of rings in this Polygon."
@ -342,51 +328,26 @@ class Polygon(GEOSGeometry):
def _checkindex(self, index):
"Internal routine for checking the given ring index."
if index < 0 or index >= len(self):
raise GEOSGeometryIndexError, 'invalid Polygon ring index: %s' % index
def _nullify(self):
"Overloaded from base method to nullify ring references as well."
# Nullifying the references to the internal rings of this Polygon.
for r in self._ptr: r.nullify()
return super(Polygon, self)._nullify()
def _populate(self):
"Internal routine for populating the internal ring pointers."
# Only populate if there aren't already children pointers.
if len(self._ptr) == 0:
# Getting the exterior ring pointer address.
ring_list = [lgeos.GEOSGetExteriorRing(self._ptr())]
# Getting the interior ring pointer addresses.
ring_list += [lgeos.GEOSGetInteriorRingN(self._ptr(), c_int(i)) for i in xrange(self.num_interior_rings)]
# Getting the coordinate sequence pointer address for each of the rings.
ptr_list = [(ring_ptr, lgeos.GEOSGeom_getCoordSeq(ring_ptr)) for ring_ptr in ring_list]
# Setting the children pointers.
self._ptr.set_children(ptr_list)
raise GEOSGeometryIndexError('invalid Polygon ring index: %s' % index)
def get_interior_ring(self, ring_i):
"""
Gets the interior ring at the specified index, 0 is for the first
interior ring, not the exterior ring.
"""
# Returning the ring from the internal ring dictionary (have to add one
# to index since all internal rings come after the exterior ring)
self._checkindex(ring_i+1)
return GEOSGeometry(self._ptr[ring_i+1], srid=self.srid)
return GEOSGeometry(geom_clone(get_intring(self._ptr, ring_i)), srid=self.srid)
#### Polygon Properties ####
@property
def num_interior_rings(self):
"Returns the number of interior rings."
# Getting the number of rings
n = lgeos.GEOSGetNumInteriorRings(self._ptr())
# -1 indicates an exception occurred
if n == -1: raise GEOSException, 'Error getting the number of interior rings.'
else: return n
return get_nrings(self._ptr)
def get_ext_ring(self):
"Gets the exterior ring of the Polygon."
return GEOSGeometry(self._ptr[0], srid=self.srid)
return GEOSGeometry(geom_clone(get_extring(self._ptr)), srid=self.srid)
def set_ext_ring(self, ring):
"Sets the exterior ring of the Polygon."
@ -399,7 +360,7 @@ class Polygon(GEOSGeometry):
@property
def tuple(self):
"Gets the tuple for each ring in this Polygon."
return tuple(self.__getitem__(i).tuple for i in xrange(self.__len__()))
return tuple(self[i].tuple for i in xrange(len(self)))
@property
def kml(self):

View File

@ -1,15 +1,14 @@
"""
This module houses the ctypes initialization procedures, as well
as the notice and error handler function callbacks (get called
when an error occurs in GEOS).
This module houses the ctypes initialization procedures, as well
as the notice and error handler function callbacks (get called
when an error occurs in GEOS).
This module also houses GEOS Pointer utilities, including
get_pointer_arr(), and GEOM_PTR.
This module also houses GEOS Pointer utilities, including
get_pointer_arr(), and GEOM_PTR.
"""
from django.contrib.gis.geos.error import GEOSException
from ctypes import c_char_p, c_int, string_at, CDLL, CFUNCTYPE, POINTER, Structure
import atexit, os, sys
from ctypes import c_char_p, string_at, Structure, CDLL, CFUNCTYPE, POINTER
from django.contrib.gis.geos.error import GEOSException
# NumPy supported?
try:
@ -25,16 +24,14 @@ if os.name == 'nt':
lib_name = 'libgeos_c-1.dll'
elif os.name == 'posix':
platform = os.uname()[0] # Using os.uname()
if platform in ('Linux', 'SunOS'):
# Linux or Solaris shared library
lib_name = 'libgeos_c.so'
elif platform == 'Darwin':
if platform == 'Darwin':
# Mac OSX Shared Library (Thanks Matt!)
lib_name = 'libgeos_c.dylib'
else:
raise GEOSException, 'Unknown POSIX platform "%s"' % platform
# Attempting to use the .so extension for all other platforms
lib_name = 'libgeos_c.so'
else:
raise GEOSException, 'Unsupported OS "%s"' % os.name
raise GEOSException('Unsupported OS "%s"' % os.name)
# Getting the GEOS C library. The C interface (CDLL) is used for
# both *NIX and Windows.
@ -71,13 +68,13 @@ lgeos.initGEOS(notice_h, error_h)
#### GEOS Geometry C data structures, and utility functions. ####
# Opaque GEOS geometry structure
class GEOSGeom_t(Structure):
"Opaque structure used when arrays of geometries are needed as parameters."
pass
# Opaque GEOS geometry structures, used for GEOM_PTR and CS_PTR
class GEOSGeom_t(Structure): pass
class GEOSCoordSeq_t(Structure): pass
# Pointer to opaque geometry structure
# Pointers to opaque GEOS geometry structures.
GEOM_PTR = POINTER(GEOSGeom_t)
CS_PTR = POINTER(GEOSCoordSeq_t)
# Used specifically by the GEOSGeom_createPolygon and GEOSGeom_createCollection
# GEOS routines

View File

@ -1,265 +0,0 @@
"""
This module houses the GEOSPointer class, used by GEOS Geometry objects
internally for memory management. Do not modify unless you _really_
know what you're doing.
"""
from ctypes import cast, c_int, c_void_p, pointer, POINTER, Structure
from django.contrib.gis.geos.error import GEOSException
from django.contrib.gis.geos.libgeos import lgeos
# This C data structure houses the memory addresses (integers) of the
# pointers returned from the GEOS C routines.
class GEOSStruct(Structure): pass
GEOSStruct._fields_ = [("geom", POINTER(c_int)), # for the geometry
("cs", POINTER(c_int)), # for geometries w/coordinate sequences
("parent", POINTER(GEOSStruct)), # points to the GEOSStruct of the parent
("child", c_void_p), # points to an array of GEOSStructs
("nchild", c_int), # the number of children
]
class GEOSPointer(object):
"""
The GEOSPointer provides a layer of abstraction in accessing the values
returned by GEOS geometry creation routines. Memory addresses (integers)
are kept in a C pointer, which allows parent geometries to be 'nullified'
when a child's memory is used in construction of another geometry.
This object is also used to store pointers for any associated coordinate
sequence and may store the pointers for any children geometries.
"""
#### Python 'magic' routines ####
def __init__(self, address, coordseq=0):
"""
Initializes on an address (an integer), another GEOSPointer, or a
GEOSStruct object.
"""
if isinstance(address, int):
self._struct = GEOSStruct()
# Integer addresses passed in, use the 'set' methods.
self.set(address)
if coordseq: self.set_coordseq(coordseq)
elif isinstance(address, GEOSPointer):
# Initializing from another GEOSPointer
self._struct = address._struct
elif isinstance(address, GEOSStruct):
# GEOSStruct passed directly in as a parameter
self._struct = address
else:
raise TypeError, 'GEOSPointer object must initialize with an integer.'
def __call__(self):
"""
Returns the address value (an integer) of the GEOS Geometry pointer.
If the pointer is NULL, a GEOSException will be raised, thus preventing
an invalid memory address from being passed to a C routine.
"""
if self.valid: return self.address
else: raise GEOSException, 'GEOS pointer no longer valid (was this geometry or the parent geometry deleted or modified?)'
def __getitem__(self, index):
"Returns a GEOSpointer object at the given child index."
n_child = len(self)
if n_child:
if index < 0 or index >= n_child:
raise IndexError, 'invalid child index'
else:
CHILD_PTR = POINTER(GEOSStruct * len(self))
return GEOSPointer(cast(self._struct.child, CHILD_PTR).contents[index])
else:
raise GEOSException, 'This GEOSPointer is not a parent'
def __iter__(self):
"""
Iterates over the children Geometry pointers, return as GEOSPointer
objects to the caller.
"""
for i in xrange(len(self)):
yield self[i]
def __len__(self):
"Returns the number of children Geometry pointers (or 0 if no children)."
return self._struct.nchild
def __nonzero__(self):
"Returns True when the GEOSPointer is valid."
return self.valid
def __str__(self):
"Returns the string representation of this GEOSPointer."
# First getting the address for the Geometry pointer.
if self.valid: geom_addr = self.address
else: geom_addr = 0
# If there's a coordinate sequence, construct accoringly.
if self.coordseq_valid:
return 'GEOSPointer(%s, %s)' % (geom_addr, self.coordseq_address)
else:
return 'GEOSPointer(%s)' % geom_addr
def __repr__(self):
return str(self)
#### GEOSPointer Properties ####
@property
def address(self):
"Returns the address of the GEOSPointer (represented as an integer)."
return self._struct.geom.contents.value
@property
def valid(self):
"Tests whether this GEOSPointer is valid."
return bool(self._struct.geom)
#### Parent & Child Properties ####
@property
def parent(self):
"Returns the GEOSPointer for the parent or None."
if self.child:
return GEOSPointer(self._struct.parent.contents)
else:
return None
@property
def child(self):
"Returns True if the GEOSPointer has a parent."
return bool(self._struct.parent)
#### Coordinate Sequence routines and properties ####
def coordseq(self):
"""
If the coordinate sequence pointer is NULL or 0, an exception will
be raised.
"""
if self.coordseq_valid: return self.coordseq_address
else: raise GEOSException, 'GEOS coordinate sequence pointer invalid (was this geometry or the parent geometry deleted or modified?)'
@property
def coordseq_address(self):
"Returns the address of the related coordinate sequence."
return self._struct.cs.contents.value
@property
def coordseq_valid(self):
"Returns True if the coordinate sequence address is valid, False otherwise."
return bool(self._struct.cs)
#### GEOSPointer Methods ####
def destroy(self):
"""
Calls GEOSGeom_destroy on the address of this pointer, and nullifies
this pointer. Use VERY carefully, as trying to destroy an address that
no longer holds a valid GEOS Geometry may crash Python.
"""
if self.valid:
# ONLY on valid geometries.
lgeos.GEOSGeom_destroy(self.address)
self.nullify()
def set(self, address):
"""
Sets the address of this pointer with the given address (represented
as an integer). Using 0 or None will set the pointer to NULL.
"""
if address in (0, None):
self._struct.geom = None
else:
self._struct.geom.contents = c_int(address)
def set_coordseq(self, address):
"""
Sets the address of the coordinate sequence associated with
this pointer.
"""
if address in (0, None):
self._struct.cs = None
else:
self._struct.cs.contents = c_int(address)
def set_children(self, ptr_list):
"""
Sets children pointers with the given pointer list (of integers).
Alternatively, a list of tuples for the geometry and coordinate
sequence pointers of the children may be used.
"""
# The number of children geometries is the number of pointers (or
# tuples) passed in via the `ptr_list`.
n_child = len(ptr_list)
# Determining whether coordinate sequences pointers were passed in.
if n_child and isinstance(ptr_list[0], (tuple, list)):
self._child_cs = True
else:
self._child_cs = False
# Dynamically creating the C types for the children array (CHILD_ARR),
# initializing with the created type, and creating a parent pointer
# for the children.
CHILD_ARR = GEOSStruct * n_child
children = CHILD_ARR()
parent = pointer(self._struct)
# Incrementing through each of the children, and setting the
# pointers in the array of GEOSStructs.
for i in xrange(n_child):
if self._child_cs:
geom_ptr, cs_ptr = ptr_list[i]
if cs_ptr is not None:
children[i].cs.contents = c_int(cs_ptr)
else:
geom_ptr = ptr_list[i]
children[i].geom.contents = c_int(geom_ptr)
children[i].parent = parent
# Casting the CHILD_ARR to the contents of the void pointer, and
# setting the number of children within the struct (used by __len__).
self._struct.child = cast(pointer(children), c_void_p)
self._struct.nchild = c_int(n_child)
def nullify(self):
"""
Nullify the geometry and coordinate sequence pointers (sets the
pointers to NULL). This does not delete any memory (destroy()
handles that), rather, it sets the GEOS pointers to a NULL address,
preventing access to deleted objects.
"""
# Nullifying both the geometry and coordinate sequence pointer.
self.set(None)
self.set_coordseq(None)
def summary(self):
"""
Returns a summary string containing information about the pointer,
including the geometry address, any associated coordinate sequence
address, and child geometry addresses.
"""
sum = '%s\n' % str(self)
for p1 in self:
sum += ' * %s\n' % p1
for p2 in p1: sum += ' - %s\n' % p2
if bool(self._struct.parent):
sum += 'Parent: %s\n' % self.parent
return sum
class NullGEOSPointer(GEOSPointer):
"The NullGEOSPointer is always NULL, and cannot be set to anything else."
def __init__(self):
self._struct = GEOSStruct()
@property
def valid(self):
return False
@property
def coordseq_valid(self):
return False
def set(self, *args):
pass
def set_coordseq(self, *args):
pass
def set_children(self, *args):
pass
NULL_GEOM = NullGEOSPointer()

View File

@ -0,0 +1,33 @@
"""
This module contains all of the GEOS ctypes function prototypes. Each
prototype handles the interaction between the GEOS library and Python
via ctypes.
"""
# Coordinate sequence routines.
from django.contrib.gis.geos.prototypes.coordseq import create_cs, get_cs, \
cs_clone, cs_getordinate, cs_setordinate, cs_getx, cs_gety, cs_getz, \
cs_setx, cs_sety, cs_setz, cs_getsize, cs_getdims
# Geometry routines.
from django.contrib.gis.geos.prototypes.geom import from_hex, from_wkb, from_wkt, \
create_point, create_linestring, create_linearring, create_polygon, create_collection, \
destroy_geom, get_extring, get_intring, get_nrings, get_geomn, geom_clone, \
geos_normalize, geos_type, geos_typeid, geos_get_srid, geos_set_srid, \
get_dims, get_num_coords, get_num_geoms, \
to_hex, to_wkb, to_wkt
# Miscellaneous routines.
from django.contrib.gis.geos.prototypes.misc import geos_area, geos_distance, geos_length
# Predicates
from django.contrib.gis.geos.prototypes.predicates import geos_hasz, geos_isempty, \
geos_isring, geos_issimple, geos_isvalid, geos_contains, geos_crosses, \
geos_disjoint, geos_equals, geos_equalsexact, geos_intersects, \
geos_intersects, geos_overlaps, geos_relatepattern, geos_touches, geos_within
# Topology routines
from django.contrib.gis.geos.prototypes.topology import \
geos_boundary, geos_buffer, geos_centroid, geos_convexhull, geos_difference, \
geos_envelope, geos_intersection, geos_pointonsurface, geos_preservesimplify, \
geos_simplify, geos_symdifference, geos_union, geos_relate

View File

@ -0,0 +1,82 @@
from ctypes import c_double, c_int, c_uint, POINTER
from django.contrib.gis.geos.libgeos import lgeos, GEOM_PTR, CS_PTR
from django.contrib.gis.geos.prototypes.errcheck import last_arg_byref, GEOSException
## Error-checking routines specific to coordinate sequences. ##
def check_cs_ptr(result, func, cargs):
"Error checking on routines that return Geometries."
if not result:
raise GEOSException('Error encountered checking Coordinate Sequence returned from GEOS C function "%s".' % func.__name__)
return result
def check_cs_op(result, func, cargs):
"Checks the status code of a coordinate sequence operation."
if result == 0:
raise GEOSException('Could not set value on coordinate sequence')
else:
return result
def check_cs_get(result, func, cargs):
"Checking the coordinate sequence retrieval."
check_cs_op(result, func, cargs)
# Object in by reference, return its value.
return last_arg_byref(cargs)
## Coordinate sequence prototype generation functions. ##
def cs_int(func):
"For coordinate sequence routines that return an integer."
func.argtypes = [CS_PTR, POINTER(c_uint)]
func.restype = c_int
func.errcheck = check_cs_get
return func
def cs_operation(func, ordinate=False, get=False):
"For coordinate sequence operations."
if get:
# Get routines get double parameter passed-in by reference.
func.errcheck = check_cs_get
dbl_param = POINTER(c_double)
else:
func.errcheck = check_cs_op
dbl_param = c_double
if ordinate:
# Get/Set ordinate routines have an extra uint parameter.
func.argtypes = [CS_PTR, c_uint, c_uint, dbl_param]
else:
func.argtypes = [CS_PTR, c_uint, dbl_param]
func.restype = c_int
return func
def cs_output(func, argtypes):
"For routines that return a coordinate sequence."
func.argtypes = argtypes
func.restype = CS_PTR
func.errcheck = check_cs_ptr
return func
## Coordinate Sequence ctypes prototypes ##
# Coordinate Sequence constructors & cloning.
cs_clone = cs_output(lgeos.GEOSCoordSeq_clone, [CS_PTR])
create_cs = cs_output(lgeos.GEOSCoordSeq_create, [c_uint, c_uint])
get_cs = cs_output(lgeos.GEOSGeom_getCoordSeq, [GEOM_PTR])
# Getting, setting ordinate
cs_getordinate = cs_operation(lgeos.GEOSCoordSeq_getOrdinate, ordinate=True, get=True)
cs_setordinate = cs_operation(lgeos.GEOSCoordSeq_setOrdinate, ordinate=True)
# For getting, x, y, z
cs_getx = cs_operation(lgeos.GEOSCoordSeq_getX, get=True)
cs_gety = cs_operation(lgeos.GEOSCoordSeq_getY, get=True)
cs_getz = cs_operation(lgeos.GEOSCoordSeq_getZ, get=True)
# For setting, x, y, z
cs_setx = cs_operation(lgeos.GEOSCoordSeq_setX)
cs_sety = cs_operation(lgeos.GEOSCoordSeq_setY)
cs_setz = cs_operation(lgeos.GEOSCoordSeq_setZ)
# These routines return size & dimensions.
cs_getsize = cs_int(lgeos.GEOSCoordSeq_getSize)
cs_getdims = cs_int(lgeos.GEOSCoordSeq_getDimensions)

View File

@ -0,0 +1,76 @@
"""
Error checking functions for GEOS ctypes prototype functions.
"""
import os
from ctypes import string_at, CDLL
from ctypes.util import find_library
from django.contrib.gis.geos.error import GEOSException
# Getting the C library, needed to free the string pointers
# returned from GEOS.
if os.name == 'nt':
libc_name = 'msvcrt'
else:
libc_name = 'libc'
libc = CDLL(find_library(libc_name))
### ctypes error checking routines ###
def last_arg_byref(args):
"Returns the last C argument's by reference value."
return args[-1]._obj.value
def check_dbl(result, func, cargs):
"Checks the status code and returns the double value passed in by reference."
# Checking the status code
if result != 1: return None
# Double passed in by reference, return its value.
return last_arg_byref(cargs)
def check_geom(result, func, cargs):
"Error checking on routines that return Geometries."
if not result:
raise GEOSException('Error encountered checking Geometry returned from GEOS C function "%s".' % func.__name__)
return result
def check_minus_one(result, func, cargs):
"Error checking on routines that should not return -1."
if result == -1:
raise GEOSException('Error encountered in GEOS C function "%s".' % func.__name__)
else:
return result
def check_predicate(result, func, cargs):
"Error checking for unary/binary predicate functions."
val = ord(result) # getting the ordinal from the character
if val == 1: return True
elif val == 0: return False
else:
raise GEOSException('Error encountered on GEOS C predicate function "%s".' % func.__name__)
def check_sized_string(result, func, cargs):
"Error checking for routines that return explicitly sized strings."
if not result:
raise GEOSException('Invalid string pointer returned by GEOS C function "%s"' % func.__name__)
# A c_size_t object is passed in by reference for the second
# argument on these routines, and its needed to determine the
# correct size.
s = string_at(result, last_arg_byref(cargs))
libc.free(result)
return s
def check_string(result, func, cargs):
"Error checking for routines that return strings."
if not result: raise GEOSException('Error encountered checking string return value in GEOS C function "%s".' % func.__name__)
# Getting the string value at the pointer address.
s = string_at(result)
# Freeing the memory allocated by the GEOS library.
libc.free(result)
return s
def check_zero(result, func, cargs):
"Error checking on routines that should not return 0."
if result == 0:
raise GEOSException('Error encountered in GEOS C function "%s".' % func.__name__)
else:
return result

View File

@ -0,0 +1,105 @@
from ctypes import c_char_p, c_int, c_size_t, c_uint, POINTER
from django.contrib.gis.geos.libgeos import lgeos, CS_PTR, GEOM_PTR
from django.contrib.gis.geos.prototypes.errcheck import \
check_geom, check_minus_one, check_sized_string, check_string, check_zero
### ctypes generation functions ###
def bin_constructor(func):
"Generates a prototype for binary construction (HEX, WKB) GEOS routines."
func.argtypes = [c_char_p, c_size_t]
func.restype = GEOM_PTR
func.errcheck = check_geom
return func
# HEX & WKB output
def bin_output(func):
"Generates a prototype for the routines that return a a sized string."
func.argtypes = [GEOM_PTR, POINTER(c_size_t)]
func.errcheck = check_sized_string
return func
def geom_output(func, argtypes):
"For GEOS routines that return a geometry."
if argtypes: func.argtypes = argtypes
func.restype = GEOM_PTR
func.errcheck = check_geom
return func
def geom_index(func):
"For GEOS routines that return geometries from an index."
return geom_output(func, [GEOM_PTR, c_int])
def int_from_geom(func, zero=False):
"Argument is a geometry, return type is an integer."
func.argtypes = [GEOM_PTR]
func.restype = c_int
if zero:
func.errcheck = check_zero
else:
func.errcheck = check_minus_one
return func
def string_from_geom(func):
"Argument is a Geometry, return type is a string."
# We do _not_ specify an argument type because we want just an
# address returned from the function.
func.argtypes = [GEOM_PTR]
func.errcheck = check_string
return func
### ctypes prototypes ###
# Creation routines from WKB, HEX, WKT
from_hex = bin_constructor(lgeos.GEOSGeomFromHEX_buf)
from_wkb = bin_constructor(lgeos.GEOSGeomFromWKB_buf)
from_wkt = geom_output(lgeos.GEOSGeomFromWKT, [c_char_p])
# Output routines
to_hex = bin_output(lgeos.GEOSGeomToHEX_buf)
to_wkb = bin_output(lgeos.GEOSGeomToWKB_buf)
to_wkt = string_from_geom(lgeos.GEOSGeomToWKT)
# The GEOS geometry type, typeid, num_coordites and number of geometries
geos_normalize = int_from_geom(lgeos.GEOSNormalize)
geos_type = string_from_geom(lgeos.GEOSGeomType)
geos_typeid = int_from_geom(lgeos.GEOSGeomTypeId)
get_dims = int_from_geom(lgeos.GEOSGeom_getDimensions, zero=True)
get_num_coords = int_from_geom(lgeos.GEOSGetNumCoordinates)
get_num_geoms = int_from_geom(lgeos.GEOSGetNumGeometries)
# Geometry creation factories
create_point = geom_output(lgeos.GEOSGeom_createPoint, [CS_PTR])
create_linestring = geom_output(lgeos.GEOSGeom_createLineString, [CS_PTR])
create_linearring = geom_output(lgeos.GEOSGeom_createLinearRing, [CS_PTR])
# Polygon and collection creation routines are special and will not
# have their argument types defined.
create_polygon = geom_output(lgeos.GEOSGeom_createPolygon, None)
create_collection = geom_output(lgeos.GEOSGeom_createCollection, None)
# Ring routines
get_extring = geom_output(lgeos.GEOSGetExteriorRing, [GEOM_PTR])
get_intring = geom_index(lgeos.GEOSGetInteriorRingN)
get_nrings = int_from_geom(lgeos.GEOSGetNumInteriorRings)
# Collection Routines
get_geomn = geom_index(lgeos.GEOSGetGeometryN)
# Cloning
geom_clone = lgeos.GEOSGeom_clone
geom_clone.argtypes = [GEOM_PTR]
geom_clone.restype = GEOM_PTR
# Destruction routine.
destroy_geom = lgeos.GEOSGeom_destroy
destroy_geom.argtypes = [GEOM_PTR]
destroy_geom.restype = None
# SRID routines
geos_get_srid = lgeos.GEOSGetSRID
geos_get_srid.argtypes = [GEOM_PTR]
geos_get_srid.restype = c_int
geos_set_srid = lgeos.GEOSSetSRID
geos_set_srid.argtypes = [GEOM_PTR, c_int]
geos_set_srid.restype = None

View File

@ -0,0 +1,27 @@
"""
This module is for the miscellaneous GEOS routines, particularly the
ones that return the area, distance, and length.
"""
from ctypes import c_int, c_double, POINTER
from django.contrib.gis.geos.libgeos import lgeos, GEOM_PTR
from django.contrib.gis.geos.prototypes.errcheck import check_dbl
### ctypes generator function ###
def dbl_from_geom(func, num_geom=1):
"""
Argument is a Geometry, return type is double that is passed
in by reference as the last argument.
"""
argtypes = [GEOM_PTR for i in xrange(num_geom)]
argtypes += [POINTER(c_double)]
func.argtypes = argtypes
func.restype = c_int # Status code returned
func.errcheck = check_dbl
return func
### ctypes prototypes ###
# Area, distance, and length prototypes.
geos_area = dbl_from_geom(lgeos.GEOSArea)
geos_distance = dbl_from_geom(lgeos.GEOSDistance, num_geom=2)
geos_length = dbl_from_geom(lgeos.GEOSLength)

View File

@ -0,0 +1,43 @@
"""
This module houses the GEOS ctypes prototype functions for the
unary and binary predicate operations on geometries.
"""
from ctypes import c_char, c_char_p, c_double
from django.contrib.gis.geos.libgeos import lgeos, GEOM_PTR
from django.contrib.gis.geos.prototypes.errcheck import check_predicate
## Binary & unary predicate functions ##
def binary_predicate(func, *args):
"For GEOS binary predicate functions."
argtypes = [GEOM_PTR, GEOM_PTR]
if args: argtypes += args
func.argtypes = argtypes
func.restype = c_char
func.errcheck = check_predicate
return func
def unary_predicate(func):
"For GEOS unary predicate functions."
func.argtypes = [GEOM_PTR]
func.restype = c_char
func.errcheck = check_predicate
return func
## Unary Predicates ##
geos_hasz = unary_predicate(lgeos.GEOSHasZ)
geos_isempty = unary_predicate(lgeos.GEOSisEmpty)
geos_isring = unary_predicate(lgeos.GEOSisRing)
geos_issimple = unary_predicate(lgeos.GEOSisSimple)
geos_isvalid = unary_predicate(lgeos.GEOSisValid)
## Binary Predicates ##
geos_contains = binary_predicate(lgeos.GEOSContains)
geos_crosses = binary_predicate(lgeos.GEOSCrosses)
geos_disjoint = binary_predicate(lgeos.GEOSDisjoint)
geos_equals = binary_predicate(lgeos.GEOSEquals)
geos_equalsexact = binary_predicate(lgeos.GEOSEqualsExact, c_double)
geos_intersects = binary_predicate(lgeos.GEOSIntersects)
geos_overlaps = binary_predicate(lgeos.GEOSOverlaps)
geos_relatepattern = binary_predicate(lgeos.GEOSRelatePattern, c_char_p)
geos_touches = binary_predicate(lgeos.GEOSTouches)
geos_within = binary_predicate(lgeos.GEOSWithin)

View File

@ -0,0 +1,35 @@
"""
This module houses the GEOS ctypes prototype functions for the
topological operations on geometries.
"""
from ctypes import c_char_p, c_double, c_int
from django.contrib.gis.geos.libgeos import lgeos, GEOM_PTR
from django.contrib.gis.geos.prototypes.errcheck import check_geom, check_string
def topology(func, *args):
"For GEOS unary topology functions."
argtypes = [GEOM_PTR]
if args: argtypes += args
func.argtypes = argtypes
func.restype = GEOM_PTR
func.errcheck = check_geom
return func
### Topology Routines ###
geos_boundary = topology(lgeos.GEOSBoundary)
geos_buffer = topology(lgeos.GEOSBuffer, c_double, c_int)
geos_centroid = topology(lgeos.GEOSGetCentroid)
geos_convexhull = topology(lgeos.GEOSConvexHull)
geos_difference = topology(lgeos.GEOSDifference, GEOM_PTR)
geos_envelope = topology(lgeos.GEOSEnvelope)
geos_intersection = topology(lgeos.GEOSIntersection, GEOM_PTR)
geos_pointonsurface = topology(lgeos.GEOSPointOnSurface)
geos_preservesimplify = topology(lgeos.GEOSTopologyPreserveSimplify, c_double)
geos_simplify = topology(lgeos.GEOSSimplify, c_double)
geos_symdifference = topology(lgeos.GEOSSymDifference, GEOM_PTR)
geos_union = topology(lgeos.GEOSUnion, GEOM_PTR)
# GEOSRelate returns a string, not a geometry.
geos_relate = lgeos.GEOSRelate
geos_relate.argtypes = [GEOM_PTR, GEOM_PTR]
geos_relate.errcheck = check_string

View File

@ -1,11 +1,12 @@
import random, unittest
import random, unittest, sys
from ctypes import ArgumentError
from django.contrib.gis.geos import \
GEOSException, GEOSGeometryIndexError, \
GEOSGeometry, Point, LineString, LinearRing, Polygon, \
MultiPoint, MultiLineString, MultiPolygon, GeometryCollection, \
fromstr, HAS_NUMPY
fromstr, geos_version, HAS_NUMPY
from django.contrib.gis.geos.base import HAS_GDAL
from geometries import *
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
@ -36,7 +37,10 @@ class GEOSTest(unittest.TestCase):
# string-based
print "\nBEGIN - expecting GEOS_ERROR; safe to ignore.\n"
for err in errors:
self.assertRaises(GEOSException, fromstr, err.wkt)
try:
g = fromstr(err.wkt)
except (GEOSException, ValueError):
pass
print "\nEND - expecting GEOS_ERROR; safe to ignore.\n"
class NotAGeometry(object):
@ -249,7 +253,7 @@ class GEOSTest(unittest.TestCase):
# Testing __getitem__ and __setitem__ on invalid indices
self.assertRaises(GEOSGeometryIndexError, poly.__getitem__, len(poly))
self.assertRaises(GEOSGeometryIndexError, poly.__setitem__, len(poly), False)
#self.assertRaises(GEOSGeometryIndexError, poly.__setitem__, len(poly), False)
self.assertRaises(GEOSGeometryIndexError, poly.__getitem__, -1)
# Testing __iter__
@ -260,24 +264,11 @@ class GEOSTest(unittest.TestCase):
# Testing polygon construction.
self.assertRaises(TypeError, Polygon.__init__, 0, [1, 2, 3])
self.assertRaises(TypeError, Polygon.__init__, 'foo')
rings = tuple(r.clone() for r in poly)
rings = tuple(r for r in poly)
self.assertEqual(poly, Polygon(rings[0], rings[1:]))
self.assertEqual(poly.wkt, Polygon(*tuple(r.clone() for r in poly)).wkt)
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)
# Setting the second point of the first ring (which should set the
# first point of the polygon).
prev = poly.clone() # Using clone() to get a copy of the current polygon
self.assertEqual(True, poly == prev) # They clone should be equal to the first
newval = (poly[0][1][0] + 5.0, poly[0][1][1] + 5.0) # really testing __getitem__ ([ring][point][tuple])
try:
poly[0][1] = ('cannot assign with', 'string values')
except TypeError:
pass
poly[0][1] = newval # setting the second point in the polygon with the newvalue (based on the old)
self.assertEqual(newval, poly[0][1]) # The point in the polygon should be the new value
self.assertEqual(False, poly == prev) # Should be different from the clone we just made
def test05b_multipolygons(self):
"Testing MultiPolygon objects."
print "\nBEGIN - expecting GEOS_NOTICE; safe to ignore.\n"
@ -321,149 +312,11 @@ class GEOSTest(unittest.TestCase):
# Deleting the polygon
del poly
# Ensuring that trying to access the deleted memory (by getting the string
# representation of the ring of a deleted polygon) raises a GEOSException
# instead of something worse..
self.assertRaises(GEOSException, str, ring1)
self.assertRaises(GEOSException, str, ring2)
# Access to these rings is OK since they are clones.
s1, s2 = str(ring1), str(ring2)
def test06b_memory_hijinks(self):
"Testing Geometry __del__() on collections."
#### Memory issues with geometries from Geometry Collections
mp = fromstr('MULTIPOINT(85 715, 235 1400, 4620 1711)')
# Getting the points
pts = [p for p in mp]
# More 'harmless' child geometry deletes
for p in pts: del p
# Cloning for comparisons
clones = [p.clone() for p in pts]
for i in xrange(len(clones)):
# Testing equivalence before & after modification
self.assertEqual(True, pts[i] == clones[i]) # before
pts[i].x = 3.14159
pts[i].y = 2.71828
self.assertEqual(False, pts[i] == clones[i]) # after
self.assertEqual(3.14159, mp[i].x) # parent x,y should be modified
self.assertEqual(2.71828, mp[i].y)
# Should raise GEOSException when trying to get geometries from the multipoint
# after it has been deleted.
parr1 = [ptr for ptr in mp._ptr]
parr2 = [p._ptr for p in mp]
del mp
for p in pts:
self.assertRaises(GEOSException, str, p) # tests p's geometry pointer
self.assertRaises(GEOSException, p.get_coords) # tests p's coordseq pointer
# Now doing this with a GeometryCollection
poly = fromstr(polygons[3].wkt) # a 'real life' polygon.
linring = fromstr(linearrings[0].wkt) # a 'real life' linear ring
# Pulling out the shell and cloning our initial geometries for later comparison.
shell = poly.shell
polyc = poly.clone()
linringc = linring.clone()
gc = GeometryCollection(poly, linring, Point(5, 23))
# Should no longer be able to access these variables
self.assertRaises(GEOSException, str, poly)
self.assertRaises(GEOSException, str, shell)
self.assertRaises(GEOSException, str, linring)
# Deep-indexing delete should be 'harmless.'
tmpr1 = gc[0][0]
del tmpr1
r1 = gc[0][0] # pulling out shell from polygon -- deep indexing
r2 = gc[1] # pulling out the ring
pnt = gc[2] # pulling the point from the geometry collection
# Now lets create a MultiPolygon from the geometry collection components
mpoly = MultiPolygon(gc[0], Polygon(gc[1]))
self.assertEqual(polyc.wkt, mpoly[0].wkt)
self.assertEqual(linringc.wkt, mpoly[1][0].wkt) # deep-indexed ring
# Should no longer be able to access the geometry collection directly
self.assertRaises(GEOSException, len, gc)
# BUT, should still be able to access the Point we obtained earlier,
# however, not the linear ring (since it is now part of the
# MultiPolygon).
self.assertEqual(5, pnt.x)
self.assertEqual(23, pnt.y)
for tmpr in [r1, r2]:
# __len__ is called on the coordinate sequence pointer --
# making sure its nullified as well
self.assertRaises(GEOSException, len, tmpr)
self.assertRaises(GEOSException, str, tmpr)
# Can't access point after deletion of parent geometry.
del gc
self.assertRaises(GEOSException, str, pnt)
# Cleaning up.
del polyc, mpoly
def test06c_memory_hijinks(self):
"Testing __init__ using other Geometries as parameters."
#### Memory issues with creating geometries from coordinate sequences within other geometries
# Creating the initial polygon from the following tuples, and then pulling out
# the individual rings.
ext_tup = ((0, 0), (0, 7), (7, 7), (7, 0), (0, 0))
itup1 = ((1, 1), (1, 2), (2, 2), (2, 1), (1, 1))
itup2 = ((4, 4), (4, 5), (5, 5), (5, 4), (4, 4))
poly1 = Polygon(LinearRing(ext_tup), LinearRing(itup1), LinearRing(itup2))
shell = poly1.shell
hole1 = poly1[1]
hole2 = poly1[2]
# Creating a Polygon from the shell and one of the holes
poly2 = Polygon(shell, hole1)
# We should no longer be able to access the original Polygon, its
# shell or its first internal ring.
self.assertRaises(GEOSException, str, poly1)
self.assertRaises(GEOSException, str, shell)
self.assertRaises(GEOSException, str, hole1)
# BUT, the second hole is still accessible.
self.assertEqual(itup2, hole2.tuple)
# Deleting the first polygon, and ensuring that
# the second hole is now gone for good.
del poly1, poly2
self.assertRaises(GEOSException, str, hole2)
#### Testing creating Geometries w/"deep-indexed" rings ####
mpoly = MultiPolygon(Polygon(LinearRing(ext_tup), LinearRing(itup1), LinearRing(itup2)),
Polygon(LinearRing(itup1)))
r0 = mpoly[0][1] # itup1, left alone
r1 = mpoly[0][0] # ext_tup
r2 = mpoly[1][0] # itup1
r3 = mpoly[0][2] # itup2
# Using the rings of the multipolygon to create a new Polygon, should
# no longer be able to access the MultiPolygon or the ring objects
# used in initialization.
p1 = Polygon(r1, r2, r3)
for g in [mpoly, r1, r2, r3]:
self.assertRaises(GEOSException, len, g)
self.assertRaises(GEOSException, str, g)
# However, the middle ring of the first Polygon was not used.
self.assertEqual(r0.tuple, itup1)
self.assertEqual(r0, p1[1])
# This deletes the leftover ring (or whenever mpoly is garbage collected)
del mpoly
# The previous hijinks tests are now moot because only clones are
# now used =)
def test08_coord_seq(self):
"Testing Coordinate Sequence objects."
@ -496,7 +349,6 @@ class GEOSTest(unittest.TestCase):
"Testing relate() and relate_pattern()."
g = fromstr('POINT (0 0)')
self.assertRaises(GEOSException, g.relate_pattern, 0, 'invalid pattern, yo')
for i in xrange(len(relate_geoms)):
g_tup = relate_geoms[i]
a = fromstr(g_tup[0].wkt)
@ -504,7 +356,7 @@ class GEOSTest(unittest.TestCase):
pat = g_tup[2]
result = g_tup[3]
self.assertEqual(result, a.relate_pattern(b, pat))
self.assertEqual(g_tup[2], a.relate(b))
self.assertEqual(pat, a.relate(b))
def test10_intersection(self):
"Testing intersects() and intersection()."
@ -569,7 +421,7 @@ class GEOSTest(unittest.TestCase):
exp_buf = fromstr(g_tup[1].wkt)
# Can't use a floating-point for the number of quadsegs.
self.assertRaises(TypeError, g.buffer, g_tup[2], float(g_tup[3]))
self.assertRaises(ArgumentError, g.buffer, g_tup[2], float(g_tup[3]))
# Constructing our buffer
buf = g.buffer(g_tup[2], g_tup[3])
@ -593,7 +445,7 @@ class GEOSTest(unittest.TestCase):
self.assertEqual(4326, pnt.srid)
pnt.srid = 3084
self.assertEqual(3084, pnt.srid)
self.assertRaises(TypeError, pnt.set_srid, '4326')
self.assertRaises(ArgumentError, pnt.set_srid, '4326')
# Testing SRID keyword on fromstr(), and on Polygon rings.
poly = fromstr(polygons[1].wkt, srid=4269)
@ -635,57 +487,56 @@ class GEOSTest(unittest.TestCase):
shell_tup = poly.shell.tuple
new_coords = []
for point in shell_tup: new_coords.append((point[0] + 500., point[1] + 500.))
shell1 = LinearRing(*tuple(new_coords))
shell2 = shell1.clone()
new_shell = LinearRing(*tuple(new_coords))
# Assigning polygon's exterior ring w/the new shell
poly.exterior_ring = shell1
self.assertRaises(GEOSException, str, shell1) # shell1 should no longer be accessible
self.assertEqual(poly.exterior_ring, shell2)
self.assertEqual(poly[0], shell2)
del poly, shell1, shell_tup # cleaning up
poly.exterior_ring = new_shell
s = str(new_shell) # new shell is still accessible
self.assertEqual(poly.exterior_ring, new_shell)
self.assertEqual(poly[0], new_shell)
### Testing the mutability of Geometry Collections
for tg in multipoints:
mp = fromstr(tg.wkt)
for i in range(len(mp)):
# Creating a random point.
pnt = mp[i].clone()
pnt = mp[i]
new = Point(random.randint(1, 100), random.randint(1, 100))
tmp = new.clone()
# Testing the assignmen
mp[i] = tmp
self.assertRaises(GEOSException, len, tmp)
# Testing the assignment
mp[i] = new
s = str(new) # what was used for the assignment is still accessible
self.assertEqual(mp[i], new)
self.assertEqual(mp[i].wkt, new.wkt)
self.assertNotEqual(pnt, mp[i])
del mp
# MultiPolygons involve much more memory management because each
# Polygon w/in the collection has its own rings.
# Polygon w/in the collection has its own rings.
for tg in multipolygons:
mpoly = fromstr(tg.wkt)
for i in xrange(len(mpoly)):
poly = mpoly[i].clone()
poly = mpoly[i]
old_poly = mpoly[i]
# Offsetting the each ring in the polygon by 500.
tmp = poly.clone()
for r in tmp:
for j in xrange(len(r)): r[j] = (r[j][0] + 500., r[j][1] + 500.)
self.assertNotEqual(poly, tmp)
new = tmp.clone() # a 'reference' copy of the geometry used in assignment
for j in xrange(len(poly)):
r = poly[j]
for k in xrange(len(r)): r[k] = (r[k][0] + 500., r[k][1] + 500.)
poly[j] = r
self.assertNotEqual(mpoly[i], poly)
# Testing the assignment
mpoly[i] = tmp
self.assertRaises(GEOSException, str, tmp)
self.assertEqual(mpoly[i], new)
self.assertNotEqual(poly, mpoly[i])
mpoly[i] = poly
s = str(poly) # Still accessible
self.assertEqual(mpoly[i], poly)
self.assertNotEqual(mpoly[i], old_poly)
# Extreme (!!) __setitem__
mpoly[0][0][0] = (3.14, 2.71)
self.assertEqual((3.14, 2.71), mpoly[0][0][0])
# Doing it more slowly..
self.assertEqual((3.14, 2.71), mpoly[0].shell[0])
del mpoly
# Extreme (!!) __setitem__ -- no longer works, have to detect
# in the first object that __setitem__ is called in the subsequent
# objects -- maybe mpoly[0, 0, 0] = (3.14, 2.71)?
#mpoly[0][0][0] = (3.14, 2.71)
#self.assertEqual((3.14, 2.71), mpoly[0][0][0])
# Doing it more slowly..
#self.assertEqual((3.14, 2.71), mpoly[0].shell[0])
#del mpoly
def test17_threed(self):
"Testing three-dimensional geometries."