2
This module contains the 'base' GEOSGeometry object -- all GEOS Geometries
3
inherit from this object.
5
# Python, ctypes and types dependencies.
7
from ctypes import addressof, byref, c_double, c_size_t
9
# super-class for mutable list behavior
10
from django.contrib.gis.geos.mutable_list import ListMixin
12
# GEOS-related dependencies.
13
from django.contrib.gis.geos.base import GEOSBase, gdal
14
from django.contrib.gis.geos.coordseq import GEOSCoordSeq
15
from django.contrib.gis.geos.error import GEOSException, GEOSIndexError
16
from django.contrib.gis.geos.libgeos import GEOM_PTR, GEOS_PREPARE
17
from django.contrib.gis.geos.mutable_list import ListMixin
19
# All other functions in this module come from the ctypes
20
# prototypes module -- which handles all interaction with
21
# the underlying GEOS library.
22
from django.contrib.gis.geos import prototypes as capi
24
# Regular expression for recognizing HEXEWKB and WKT. A prophylactic measure
25
# to prevent potentially malicious input from reaching the underlying C
26
# library. Not a substitute for good web security programming practices.
27
hex_regex = re.compile(r'^[0-9A-F]+$', re.I)
28
wkt_regex = re.compile(r'^(SRID=(?P<srid>\d+);)?(?P<wkt>(POINT|LINESTRING|LINEARRING|POLYGON|MULTIPOINT|MULTILINESTRING|MULTIPOLYGON|GEOMETRYCOLLECTION)[ACEGIMLONPSRUTY\d,\.\-\(\) ]+)$', re.I)
30
class GEOSGeometry(GEOSBase, ListMixin):
31
"A class that, generally, encapsulates a GEOS geometry."
33
# Raise GEOSIndexError instead of plain IndexError
34
# (see ticket #4740 and GEOSIndexError docstring)
35
_IndexError = GEOSIndexError
39
#### Python 'magic' routines ####
40
def __init__(self, geo_input, srid=None):
42
The base constructor for GEOS geometry objects, and may take the
47
- HEXEWKB (a PostGIS-specific canonical form)
48
- GeoJSON (requires GDAL)
52
The `srid` keyword is used to specify the Source Reference Identifier
53
(SRID) number for this Geometry. If not set, the SRID will be None.
55
if isinstance(geo_input, basestring):
56
if isinstance(geo_input, unicode):
57
# Encoding to ASCII, WKT or HEXEWKB doesn't need any more.
58
geo_input = geo_input.encode('ascii')
60
wkt_m = wkt_regex.match(geo_input)
63
if wkt_m.group('srid'): srid = int(wkt_m.group('srid'))
64
g = wkt_r.read(wkt_m.group('wkt'))
65
elif hex_regex.match(geo_input):
66
# Handling HEXEWKB input.
67
g = wkb_r.read(geo_input)
68
elif gdal.GEOJSON and gdal.geometries.json_regex.match(geo_input):
69
# Handling GeoJSON input.
70
g = wkb_r.read(gdal.OGRGeometry(geo_input).wkb)
72
raise ValueError('String or unicode input unrecognized as WKT EWKT, and HEXEWKB.')
73
elif isinstance(geo_input, GEOM_PTR):
74
# When the input is a pointer to a geomtry (GEOM_PTR).
76
elif isinstance(geo_input, buffer):
77
# When the input is a buffer (WKB).
78
g = wkb_r.read(geo_input)
79
elif isinstance(geo_input, GEOSGeometry):
80
g = capi.geom_clone(geo_input.ptr)
82
# Invalid geometry type.
83
raise TypeError('Improper geometry input type: %s' % str(type(geo_input)))
86
# Setting the pointer object with a valid pointer.
89
raise GEOSException('Could not initialize GEOS Geometry with given input.')
91
# Post-initialization setup.
94
def _post_init(self, srid):
95
"Helper routine for performing post-initialization setup."
96
# Setting the SRID, if given.
97
if srid and isinstance(srid, int): self.srid = srid
99
# Setting the class type (e.g., Point, Polygon, etc.)
100
self.__class__ = GEOS_CLASSES[self.geom_typeid]
102
# Setting the coordinate sequence for the geometry (will be None on
103
# geometries that do not have coordinate sequences)
108
Destroys this Geometry; in other words, frees the memory used by the
111
if self._ptr: capi.destroy_geom(self._ptr)
115
Returns a clone because the copy of a GEOSGeometry may contain an
116
invalid pointer location if the original is garbage collected.
120
def __deepcopy__(self, memodict):
122
The `deepcopy` routine is used by the `Node` class of django.utils.tree;
123
thus, the protocol routine needs to be implemented to return correct
124
copies (clones) of these GEOS objects, which use C pointers.
129
"WKT is used for the string representation."
133
"Short-hand representation because WKT may be very large."
134
return '<%s object at %s>' % (self.geom_type, hex(addressof(self.ptr)))
137
def __getstate__(self):
138
# The pickled state is simply a tuple of the WKB (in string form)
140
return str(self.wkb), self.srid
142
def __setstate__(self, state):
143
# Instantiating from the tuple state that was pickled.
145
ptr = capi.from_wkb(wkb, len(wkb))
146
if not ptr: raise GEOSException('Invalid Geometry loaded from pickled state.')
148
self._post_init(srid)
150
# Comparison operators
151
def __eq__(self, other):
153
Equivalence testing, a Geometry may be compared with another Geometry
154
or a WKT representation.
156
if isinstance(other, basestring):
157
return self.wkt == other
158
elif isinstance(other, GEOSGeometry):
159
return self.equals_exact(other)
163
def __ne__(self, other):
164
"The not equals operator."
165
return not (self == other)
167
### Geometry set-like operations ###
168
# Thanks to Sean Gillies for inspiration:
169
# http://lists.gispython.org/pipermail/community/2007-July/001034.html
171
def __or__(self, other):
172
"Returns the union of this Geometry and the other."
173
return self.union(other)
176
def __and__(self, other):
177
"Returns the intersection of this Geometry and the other."
178
return self.intersection(other)
181
def __sub__(self, other):
182
"Return the difference this Geometry and the other."
183
return self.difference(other)
186
def __xor__(self, other):
187
"Return the symmetric difference of this Geometry and the other."
188
return self.sym_difference(other)
190
#### Coordinate Sequence Routines ####
193
"Returns True if this Geometry has a coordinate sequence, False if not."
194
# Only these geometries are allowed to have coordinate sequences.
195
if isinstance(self, (Point, LineString, LinearRing)):
201
"Sets the coordinate sequence for this Geometry."
203
self._cs = GEOSCoordSeq(capi.get_cs(self.ptr), self.hasz)
209
"Returns a clone of the coordinate sequence for this Geometry."
211
return self._cs.clone()
213
#### Geometry Info ####
216
"Returns a string representing the Geometry type, e.g. 'Polygon'"
217
return capi.geos_type(self.ptr)
220
def geom_typeid(self):
221
"Returns an integer representing the Geometry type."
222
return capi.geos_typeid(self.ptr)
226
"Returns the number of geometries in the Geometry."
227
return capi.get_num_geoms(self.ptr)
230
def num_coords(self):
231
"Returns the number of coordinates in the Geometry."
232
return capi.get_num_coords(self.ptr)
235
def num_points(self):
236
"Returns the number points, or coordinates, in the Geometry."
237
return self.num_coords
241
"Returns the dimension of this Geometry (0=point, 1=line, 2=surface)."
242
return capi.get_dims(self.ptr)
245
"Converts this Geometry to normal form (or canonical form)."
246
return capi.geos_normalize(self.ptr)
248
#### Unary predicates ####
252
Returns a boolean indicating whether the set of points in this Geometry
255
return capi.geos_isempty(self.ptr)
259
"Returns whether the geometry has a 3D dimension."
260
return capi.geos_hasz(self.ptr)
264
"Returns whether or not the geometry is a ring."
265
return capi.geos_isring(self.ptr)
269
"Returns false if the Geometry not simple."
270
return capi.geos_issimple(self.ptr)
274
"This property tests the validity of this Geometry."
275
return capi.geos_isvalid(self.ptr)
277
#### Binary predicates. ####
278
def contains(self, other):
279
"Returns true if other.within(this) returns true."
280
return capi.geos_contains(self.ptr, other.ptr)
282
def crosses(self, other):
284
Returns true if the DE-9IM intersection matrix for the two Geometries
285
is T*T****** (for a point and a curve,a point and an area or a line and
286
an area) 0******** (for two curves).
288
return capi.geos_crosses(self.ptr, other.ptr)
290
def disjoint(self, other):
292
Returns true if the DE-9IM intersection matrix for the two Geometries
295
return capi.geos_disjoint(self.ptr, other.ptr)
297
def equals(self, other):
299
Returns true if the DE-9IM intersection matrix for the two Geometries
302
return capi.geos_equals(self.ptr, other.ptr)
304
def equals_exact(self, other, tolerance=0):
306
Returns true if the two Geometries are exactly equal, up to a
309
return capi.geos_equalsexact(self.ptr, other.ptr, float(tolerance))
311
def intersects(self, other):
312
"Returns true if disjoint returns false."
313
return capi.geos_intersects(self.ptr, other.ptr)
315
def overlaps(self, other):
317
Returns true if the DE-9IM intersection matrix for the two Geometries
318
is T*T***T** (for two points or two surfaces) 1*T***T** (for two curves).
320
return capi.geos_overlaps(self.ptr, other.ptr)
322
def relate_pattern(self, other, pattern):
324
Returns true if the elements in the DE-9IM intersection matrix for the
325
two Geometries match the elements in pattern.
327
if not isinstance(pattern, basestring) or len(pattern) > 9:
328
raise GEOSException('invalid intersection matrix pattern')
329
return capi.geos_relatepattern(self.ptr, other.ptr, pattern)
331
def touches(self, other):
333
Returns true if the DE-9IM intersection matrix for the two Geometries
334
is FT*******, F**T***** or F***T****.
336
return capi.geos_touches(self.ptr, other.ptr)
338
def within(self, other):
340
Returns true if the DE-9IM intersection matrix for the two Geometries
343
return capi.geos_within(self.ptr, other.ptr)
345
#### SRID Routines ####
347
"Gets the SRID for the geometry, returns None if no SRID is set."
348
s = capi.geos_get_srid(self.ptr)
349
if s == 0: return None
352
def set_srid(self, srid):
353
"Sets the SRID for the geometry."
354
capi.geos_set_srid(self.ptr, srid)
355
srid = property(get_srid, set_srid)
357
#### Output Routines ####
360
"Returns the EWKT (WKT + SRID) of the Geometry."
361
if self.get_srid(): return 'SRID=%s;%s' % (self.srid, self.wkt)
362
else: return self.wkt
366
"Returns the WKT (Well-Known Text) of the Geometry."
367
return wkt_w.write(self)
372
Returns the HEX of the Geometry -- please note that the SRID is not
373
included in this representation, because the GEOS C library uses
374
-1 by default, even if the SRID is set.
376
# A possible faster, all-python, implementation:
377
# str(self.wkb).encode('hex')
378
return wkb_w.write_hex(self)
383
Returns GeoJSON representation of this Geometry if GDAL 1.5+
389
raise GEOSException('GeoJSON output only supported on GDAL 1.5+.')
394
"Returns the WKB of the Geometry as a buffer."
395
return wkb_w.write(self)
399
"Returns the KML representation of this Geometry."
400
gtype = self.geom_type
401
return '<%s>%s</%s>' % (gtype, self.coord_seq.kml, gtype)
406
Returns a PreparedGeometry corresponding to this geometry -- it is
407
optimized for the contains, intersects, and covers operations.
410
return PreparedGeometry(self)
412
raise GEOSException('GEOS 3.1+ required for prepared geometry support.')
414
#### GDAL-specific output routines ####
417
"Returns the OGR Geometry for this Geometry."
420
return gdal.OGRGeometry(self.wkb, self.srid)
422
return gdal.OGRGeometry(self.wkb)
424
raise GEOSException('GDAL required to convert to an OGRGeometry.')
428
"Returns the OSR SpatialReference for SRID of this Geometry."
431
return gdal.SpatialReference(self.srid)
435
raise GEOSException('GDAL required to return a SpatialReference object.')
439
"Alias for `srs` property."
442
def transform(self, ct, clone=False):
444
Requires GDAL. Transforms the geometry according to the given
445
transformation object, which may be an integer SRID, and WKT or
446
PROJ.4 string. By default, the geometry is transformed in-place and
447
nothing is returned. However if the `clone` keyword is set, then this
448
geometry will not be modified and a transformed clone will be returned
452
if gdal.HAS_GDAL and srid:
453
# Creating an OGR Geometry, which is then transformed.
454
g = gdal.OGRGeometry(self.wkb, srid)
456
# Getting a new GEOS pointer
457
ptr = wkb_r.read(g.wkb)
459
# User wants a cloned transformed geometry returned.
460
return GEOSGeometry(ptr, srid=g.srid)
462
# Reassigning pointer, and performing post-initialization setup
463
# again due to the reassignment.
464
capi.destroy_geom(self.ptr)
466
self._post_init(g.srid)
468
raise GEOSException('Transformed WKB was invalid.')
470
#### Topology Routines ####
471
def _topology(self, gptr):
472
"Helper routine to return Geometry from the given pointer."
473
return GEOSGeometry(gptr, srid=self.srid)
477
"Returns the boundary as a newly allocated Geometry object."
478
return self._topology(capi.geos_boundary(self.ptr))
480
def buffer(self, width, quadsegs=8):
482
Returns a geometry that represents all points whose distance from this
483
Geometry is less than or equal to distance. Calculations are in the
484
Spatial Reference System of this Geometry. The optional third parameter sets
485
the number of segment used to approximate a quarter circle (defaults to 8).
486
(Text from PostGIS documentation at ch. 6.1.3)
488
return self._topology(capi.geos_buffer(self.ptr, width, quadsegs))
493
The centroid is equal to the centroid of the set of component Geometries
494
of highest dimension (since the lower-dimension geometries contribute zero
495
"weight" to the centroid).
497
return self._topology(capi.geos_centroid(self.ptr))
500
def convex_hull(self):
502
Returns the smallest convex Polygon that contains all the points
505
return self._topology(capi.geos_convexhull(self.ptr))
507
def difference(self, other):
509
Returns a Geometry representing the points making up this Geometry
510
that do not make up other.
512
return self._topology(capi.geos_difference(self.ptr, other.ptr))
516
"Return the envelope for this geometry (a polygon)."
517
return self._topology(capi.geos_envelope(self.ptr))
519
def intersection(self, other):
520
"Returns a Geometry representing the points shared by this Geometry and other."
521
return self._topology(capi.geos_intersection(self.ptr, other.ptr))
524
def point_on_surface(self):
525
"Computes an interior point of this Geometry."
526
return self._topology(capi.geos_pointonsurface(self.ptr))
528
def relate(self, other):
529
"Returns the DE-9IM intersection matrix for this Geometry and the other."
530
return capi.geos_relate(self.ptr, other.ptr)
532
def simplify(self, tolerance=0.0, preserve_topology=False):
534
Returns the Geometry, simplified using the Douglas-Peucker algorithm
535
to the specified tolerance (higher tolerance => less points). If no
536
tolerance provided, defaults to 0.
538
By default, this function does not preserve topology - e.g. polygons can
539
be split, collapse to lines or disappear holes can be created or
540
disappear, and lines can cross. By specifying preserve_topology=True,
541
the result will have the same dimension and number of components as the
542
input. This is significantly slower.
544
if preserve_topology:
545
return self._topology(capi.geos_preservesimplify(self.ptr, tolerance))
547
return self._topology(capi.geos_simplify(self.ptr, tolerance))
549
def sym_difference(self, other):
551
Returns a set combining the points in this Geometry not in other,
552
and the points in other not in this Geometry.
554
return self._topology(capi.geos_symdifference(self.ptr, other.ptr))
556
def union(self, other):
557
"Returns a Geometry representing all the points in this Geometry and other."
558
return self._topology(capi.geos_union(self.ptr, other.ptr))
560
#### Other Routines ####
563
"Returns the area of the Geometry."
564
return capi.geos_area(self.ptr, byref(c_double()))
566
def distance(self, other):
568
Returns the distance between the closest points on this Geometry
569
and the other. Units will be in those of the coordinate system of
572
if not isinstance(other, GEOSGeometry):
573
raise TypeError('distance() works only on other GEOS Geometries.')
574
return capi.geos_distance(self.ptr, other.ptr, byref(c_double()))
579
Returns the extent of this geometry as a 4-tuple, consisting of
580
(xmin, ymin, xmax, ymax).
583
if isinstance(env, Point):
584
xmin, ymin = env.tuple
585
xmax, ymax = xmin, ymin
587
xmin, ymin = env[0][0]
588
xmax, ymax = env[0][2]
589
return (xmin, ymin, xmax, ymax)
594
Returns the length of this Geometry (e.g., 0 for point, or the
595
circumfrence of a Polygon).
597
return capi.geos_length(self.ptr, byref(c_double()))
600
"Clones this Geometry."
601
return GEOSGeometry(capi.geom_clone(self.ptr), srid=self.srid)
603
# Class mapping dictionary. Has to be at the end to avoid import
604
# conflicts with GEOSGeometry.
605
from django.contrib.gis.geos.linestring import LineString, LinearRing
606
from django.contrib.gis.geos.point import Point
607
from django.contrib.gis.geos.polygon import Polygon
608
from django.contrib.gis.geos.collections import GeometryCollection, MultiPoint, MultiLineString, MultiPolygon
609
GEOS_CLASSES = {0 : Point,
616
7 : GeometryCollection,
619
# Similarly, import the GEOS I/O instances here to avoid conflicts.
620
from django.contrib.gis.geos.io import wkt_r, wkt_w, wkb_r, wkb_w
622
# If supported, import the PreparedGeometry class.
624
from django.contrib.gis.geos.prepared import PreparedGeometry