======== GEOS API ======== Background ========== What is GEOS? ------------- `GEOS`__ stands for **G**\ eometry **E**\ ngine - **O**\ pen **S**\ ource, and is a C++ port of the `Java Topology Suite`__, implementing the OpenGIS `Simple Features for SQL`__ spatial predicate functions and spatial operators. GEOS, now an OSGeo project, was initially developed and maintained by `Refractions Research`__ of Victoria, Canada. __ http://trac.osgeo.org/geos/ __ http://sourceforge.net/projects/jts-topo-suite/ __ http://www.opengeospatial.org/standards/sfs __ http://www.refractions.net/ Why the new API? ---------------- 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. Cross-platform compatibility. Thus, the Python ``ctypes`` [#]_ package was used to wrap the `GEOS C API`__ to bring the rich capabilities of GEOS to Python and GeoDjango. Features: * A BSD-licensed interface to the GEOS geometry routines, implemented purely in Python using ``ctypes``. * Loosely-coupled to GeoDjango. For example, GEOS geometry objects may be used outside a django project/application (no need to have ``DJANGO_SETTINGS_MODULE`` set, etc.) * Mutability. GEOS Geometry objects may be modified. * Cross-platform and tested. GeoDjango GEOS geometries are well-tested and compatible with Windows, Linux, Solaris, and Mac OS X platforms. __ http://www.swig.org/ __ http://trac.gispython.org/projects/PCL/browser/PCL/trunk/PCL-Core/cartography/geometry __ http://trac.osgeo.org/geos/browser/trunk/capi/geos_c.h.in Related Work ------------ The ``GEOSGeometry`` interface was first committed to the Django SVN by Justin Bronn (the lead developer of GeoDjango) in April, 2007. [#]_ After learning of GeoDjango's interface [#]_, Sean Gillies created his own ``ctypes`` interface, "Shapely." [#]_ Justin Bronn declined to merge the projects because Shapely was initially licensed under the LGPL, a license incompatible with GeoDjango. [#]_ While Shapely was later relicensed under the BSD license (the same as GeoDjango) [#]_, differences in design and implementation prevent the projects from merging at the moment, but are not irreconcilable. Geometry Objects ================ .. _point: .. class:: Point The ``Point`` object may be initialized with either a tuple, or individual parameters. For example:: >>> from django.contrib.gis.geos import Point >>> p = Point((5, 23)) # 2D point, passed in as a tuple >>> p = Point(5, 23) # Same, passed in with individual parameters The ``srid`` keyword may be used to specify the SRID for the point:: >>> pnt = Point(5, 23, srid=4326) Additionally, 3D geometries may be created by specifing a Z value:: >>> pnt_3d = Point(5, 23, 17) # Also: Point( (5, 23, 17) ) >>> pnt_3d.hasz True >>> pnt_3d.z 17 Properties ---------- .. attribute:: x This property sets or retrieves the X coordinate for the point. .. attribute:: y This property sets or retrieves the Y coordinate for the point. .. attribute:: z This property sets or retrieves the Z (3D) coordinate for the point. If the point has no Z coordinate, ``None`` is returned. .. note: The ``Point`` must have been created as 3D before it is possible to set the ``z`` property. .. _linestring: .. class:: LineString ``LineString`` objects initialize on a given sequence. For example, the constructor may take lists, tuples, `NumPy`__ arrays of X,Y[,Z] pairs, or ``Point`` objects. If ``Point`` objects are used, ownership of the points is *not* transferred to the ``LineString`` object. Examples:: >>> from django.contrib.gis.geos import LineString, Point >>> ls = LineString((1, 1), (2, 2)) >>> ls = LineString([(1, 1), (2, 2)]) >>> ls = LineString(Point(1, 1), Point(2, 2)) >>> from numpy import array >>> ls = LineString(array([(1, 1), (2, 2)])) .. class:: LinearRing ``LinearRing`` objects are subclasses of ``LineString``; however, an error will be raised if the points used during initialization are not closed (the first point is equal to the last point):: >>> from django.contrib.gis.geos import LinearRing >>> lr = LinearRing((0, 0), (0, 1), (1, 1), (1, 0), (0, 0)) >>> lr = LinearRing((0, 0), (0, 1)) GEOS_ERROR: IllegalArgumentException: points must form a closed linestring __ http://numpy.scipy.org/ Methods ------- .. method:: LineString(coords) .. attribute:: x This property returns the X values in a list, or a NumPy array (if installed). .. attribute:: y This property returns the Y values in a list, or a NumPy array (if installed). .. attribute:: z This property returns the Z values in a list, or a NumPy array (if installed). .. _polygon: .. class:: Polygon Polygons are composed of an exterior ring (the shell), and may also have interior rings that denote areas excluded from the exterior. Methods ^^^^^^^ .. method:: Polygon(rings) The ``Polygon`` initializes on arguments that are either ``LinearRing`` instances or may be accepted by the ``LinearRing`` constructor. .. method:: from_bbox(bbox) .. versionadded:: 1.1 Returns a new ``Polygon`` object for the given bounding box. The bounding box should be a four-tuple comprising the X and Y minimum values followed by the X and Y maximum values. For example:: >>> from django.contrib.gis.geos import Polygon >>> bbox = (0, 0, 5, 5) >>> poly = Polygon.from_bbox(bbox) >>> print poly POLYGON ((0.0000000000000000 0.0000000000000000, 0.0000000000000000 5.0000000000000000, 5.0000000000000000 5.0000000000000000, 5.0000000000000000 0.0000000000000000, 0.0000000000000000 0.0000000000000000)) .. method:: shell() .. method:: exterior_ring() Returns a ``LinearRing`` corresponding to the exterior ring, or shell, of the ``Polygon``. .. method:: num_interior_rings Returns the number of interior rings contained in the ``Polygon``. Geometry Collections ==================== .. _multipoint: .. class:: MultiPoint >>> from django.contrib.gis.geos import Point >>> from django.contrib.gis.geos import MultiPoint >>> p1 = Point((0,0)) >>> p2 = Point((1,2)) >>> mp = MultiPoint(p1, p2) >>> mp >>> mp[0].wkt 'POINT (0.0000000000000000 0.0000000000000000)' >>> len(mp) 2 >>> [ p.wkt for p in mp ] ['POINT (0.0000000000000000 0.0000000000000000)', 'POINT (1.0000000000000000 2.0000000000000000)'] >>> mp.ring False .. _multilinestring: .. class:: MultiLineString .. _multipolygon: .. class:: MultiPolygon Methods ------- .. method:: cascaded_union .. versionadded:: 1.1 .. note:: Use of this method requires at least GEOS 3.1. .. _geometrycollection: .. class:: GeometryCollection .. _prepared: Prepared Geometries =================== .. versionadded: 1.1 In order to obtain a prepared geometry, just access the ``prepared`` property on any ``GEOSGeometry`` object. Once you have a ``PreparedGeometry`` instance its spatial predicate methods, listed below, may be used with other ``GEOSGeometry`` objects. An operation with a prepared geometry can be orders of magnitude faster -- the more complex the geometry that is prepared, the larger the speedup in the operation. .. note:: GEOS 3.1 is *required* in order to use prepared geometries. For example:: >>> from django.contrib.gis.geos import Point, Polygon >>> poly = Polygon.from_bbox((0, 0, 5, 5)) >>> prep_poly = poly.prepared >>> prep_poly.contains(Point(2.5, 2.5)) True .. class:: PreparedGeometry All methods on ``PreparedGeometry`` take an ``other`` argument, which must be a ``GEOSGeometry`` object. .. method:: contains(other) .. method:: contains_properly(other) .. method:: covers(other) .. method:: intersects(other) I/O Objects =========== .. versionadded: 1.1 Inside GEOS, I/O classes are used to output geometries to WKT, WKB, and EWKB. GeoDjango allows access to these I/O classes if finer-grained control of serialization is required. Typically, Reader Objects -------------- The reader I/O classes simply return a ``GEOSGeometry`` instance from the WKB and/or WKT input given to their ``read(geom)`` method. .. class:: WKBReader Example:: >>> from django.contrib.gis.geos import WKBReader >>> wkb_r = WKBReader() >>> wkb_r.read('0101000000000000000000F03F000000000000F03F') .. class:: WKTReader Example:: >>> from django.contrib.gis.geos import WKTReader >>> wkt_r = WKTReader() >>> wkt_r.read('POINT(1 1)') Writer Objects -------------- All writer objects have a ``write(geom)`` method that returns either the WKB or WKT of the given geometry. In addition, ``WKBWriter`` objects also have properties that may be used to change the byte order, and or include the SRID and 3D values (in other words, EWKB). .. class:: WKBWriter The ``WKBWriter`` provides the most control over its output. By default it returns OGC-compliant WKB when it's ``write`` method is called. However, it has properties that .. function:: write(geom) Returns the WKB of the given geometry as a Python ``buffer`` object. Example:: >>> from django.contrib.gis.geos import Point, WKBWriter >>> pnt = Point(1, 1) >>> wkb_w = WKBWriter() >>> wkb_w.write(pnt) .. function:: write_hex(geom) Returns WKB of the geometry in hexadecimal. Example:: >>> from django.contrib.gis.geos import Point, WKBWriter >>> pnt = Point(1, 1) >>> wkb_w = WKBWriter() >>> wkb_w.write_hex(pnt) '0101000000000000000000F03F000000000000F03F' .. function:: byteorder This property may be be set to change the byte-order of the geometry representation. =============== ================================================= Byteorder Value Description =============== ================================================= 0 Big Endian (e.g., compatible with RISC systems) 1 Little Endian (e.g., compatible with x86 systems) =============== ================================================= Example:: >>> from django.contrib.gis.geos import Point, WKBWriter >>> wkb_w = WKBWriter() >>> pnt = Point(1, 1) >>> wkb_w.write_hex(pnt) '0101000000000000000000F03F000000000000F03F' >>> wkb_w.byteorder = 0 '00000000013FF00000000000003FF0000000000000' .. attribute:: outdim This property may be set to change the output dimension of the geometry representation. In other words, if you have a 3D geometry then set to 3 so that the Z value is included in the WKB. ============ =========================== Outdim Value Description ============ =========================== 2 The default, output 2D WKB. 3 Output 3D EWKB. ============ =========================== Example:: >>> from django.contrib.gis.geos import Point, WKBWriter >>> wkb_w = WKBWriter() >>> wkb_w.outdim 2 >>> pnt = Point(1, 1, 1) >>> wkb_w.write_hex(pnt) # By default, no Z value included: '0101000000000000000000F03F000000000000F03F' >>> wkb_w.outdim = 3 # Tell writer to include Z values >>> wkb_w.write_hex(pnt) '0101000080000000000000F03F000000000000F03F000000000000F03F' .. attribute:: srid Set this property with a boolean to indicate whether the SRID of the geometry should be included with the WKB representation. Example:: >>> from django.contrib.gis.geos import Point, WKBWriter >>> wkb_w = WKBWriter() >>> pnt = Point(1, 1, srid=4326) >>> wkb_w.write_hex(pnt) # By default, no SRID included: '0101000000000000000000F03F000000000000F03F' >>> wkb_w.srid = True # Tell writer to include SRID >>> wkb_w.write_hex(pnt) '0101000020E6100000000000000000F03F000000000000F03F' .. class:: WKTWriter .. method:: write(geom) Returns the WKT of the given geometry. Example:: >>> from django.contrib.gis.geos import Point, WKTWriter >>> pnt = Point(1, 1) >>> wkt_w = WKTWriter() >>> wkt_w.write(pnt) 'POINT (1.0000000000000000 1.0000000000000000)' API === Creation -------- .. class:: GEOSGeometry(geo_input, srid=None) Geometries may be created using the ``GEOSGeometry`` constructor; ``geo_input`` may be any of the following types of data: ============ ====================== ========================================================================================= Input Type Python Type Example ============ ====================== ========================================================================================= WKT ``str`` / ``unicode`` ``'POINT(5 23)'`` EWKT ``str`` / ``unicode`` ``'SRID=4326;POINT(5 23)'`` HEX ``str`` / ``unicode`` ``'010100000000000000000014400000000000003740'`` HEXEWKB ``str`` / ``unicode`` ``''0101000020E610000000000000000014400000000000003740'`` GeoJSON ``str`` / ``unicode`` ``'{ "type": "Point", "coordinates": [ 5.000000, 23.000000 ] }'`` WKB ``buffer`` ``buffer('\x01\x01\x00\x00\x00\x00\x00\x00\x00\x00\x00\x14@\x00\x00\x00\x00\x00\x007@')`` ============ ====================== ========================================================================================= While this may seem like an acronym soup, all are standardized geospatial data formats or extensions thereof. The ``srid`` keyword may be used to set the spatial reference system identifier number for the geometry. This will be used to conduct any needed transformations for spatial lookups and geographic model creation. .. method:: fromfile(file_h) This factory creates a GEOS geometry from the given file name or an open file handle (1.1 only):: >>> from django.contrib.gis.geos import fromfile >>> g = fromfile('/home/bob/geom.wkt') .. method:: fromstr(string) GEOS geometry objects may be created from strings using the ``fromstr()`` factory, or using the constructor for each geometry object (as described above):: >>> from django.contrib.gis.geos import fromstr >>> pnt = fromstr('POINT(-90.5 29.5)', srid=4326) It should be noted that ``fromstr`` is a shortcut to the constructor for the base ``GEOSGeometry`` object. Geometry Properties ------------------- .. attribute:: empty Returns whether or not the set of points in the geometry is empty. .. attribute:: geom_type Returns a string corresponding to the type of geometry. For example:: >>> pnt = GEOSGeometry('POINT(5 23)') >>> pnt.geom_type 'Point' .. attribute:: geom_typeid Returns the GEOS geometry type identification number. The following table shows the value for each geometry type: ====================== ======== Geometry ID ====================== ======== ``Point`` 0 ``LineString`` 1 ``LinearRing`` 2 ``Polygon`` 3 ``MultiPoint`` 4 ``MultiLineString`` 5 ``MultiPolygon`` 6 ``GeometryCollection`` 7 ====================== ======== .. attribute:: hasz Returns a boolean indicating whether the geometry is three-dimensional. .. attribute:: num_coords .. attribute:: num_points Returns the total number of coordinates in this geometry. For polygons and collections, this is the cumulative number of all coordinates from the component geometries. .. attribute:: ring Returns a boolean indicating whether the geometry is a ``LinearRing``. .. attribute:: simple A Geometry is simple if and only if the only self-intersections are at boundary points. For example, a ``LineString`` object is not simple if it intersects itself. Thus, ``LinearRing`` and ``Polygon`` objects are always simple because they do not intersect themselves. .. attribute:: valid Returns a boolean indicating whether the geometry is valid. Output Properties ----------------- .. _ewkt: .. attribute:: ewkt Returns the "extended" Well-Known Text of the geometry. This representation is specific to PostGIS and is a super set of the OGC WKT standard. [#]_ Essentially the SRID is prepended to the WKT representation, for example ``SRID=4326;POINT(5 23)``. Please note that this does not include the 3dm, 3dz, and 4d information that PostGIS supports in its EWKT representations. .. _hex: .. attribute:: hex Returns the WKB of this Geometry in hexadecimal form. Please note that the SRID and Z values are not included in this representation because it is not a part of the OGC specification (use the :ref:`hexewkb` property instead). .. _hexewkb: .. versionadded:: 1.2 .. attribute:: hexewkb Returns the EWKB of this Geometry in hexadecimal form. This is an extension of the WKB specification that includes SRID and Z values that are a part of this geometry. .. note:: GEOS 3.1 is *required* if you want valid 3D HEXEWKB. .. attribute:: json .. attribute:: geojson Returns the GeoJSON representation of the geometry. Requires GDAL. .. attribute:: kml Returns a `KML`__ (Keyhole Markup Language) representation of the geometry. This should only be used for geometries with an SRID of 4326 (WGS84), but this restriction is not enforced. .. attribute:: ogr Returns an OGR ``OGRGeometry`` object correspondg to the GEOS geometry. Consult the `OGRGeometry`_ documentation for more information. .. _OGRGeometry: gdal.html#ogrgeometry .. note:: Requires GDAL. .. _wkb: .. attribute:: wkb Returns the WKB (Well-Known Binary) representation of this Geometry as a Python buffer. SRID and Z values are not included, use the :ref:`ewkb` property instead. .. _ewkb: .. attribute:: ewkb .. versionadded:: 1.2 Return the EWKB representation of this Geometry as a Python buffer. This is an extension of the WKB specification that includes any SRID and Z values that are a part of this geometry. .. note:: GEOS 3.1 is *required* if you want valid 3D EWKB. .. attribute:: wkt Returns the Well-Known Text of the geometry (an OGC standard). __ http://code.google.com/apis/kml/documentation/ Spatial Predicate Methods ------------------------- All of the following spatial predicate methods take another GEOS Geometry instance (``other``) as an argument. .. method:: contains(other) Returns True if ``within(other)`` is False. .. method:: crosses(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). .. method:: disjoint(other) Returns true if the DE-9IM intersection matrix for the two Geometries is ``FF*FF****``. .. method:: equals(other) Returns true if the DE-9IM intersection matrix for the two Geometries is ``T*F**FFF*``. .. method:: equals_exact(other, tolerance=0) Returns true if the two Geometries are exactly equal, up to a specified tolerance. The ``tolerance`` value should be a floating point number representing the error tolerance in the comparison, e.g., ``poly1.equals_exact(poly2, 0.001)`` will compare equality to within one thousandth of a unit. .. method:: intersects(other) Returns True if ``disjoint(other)`` is False. .. method:: overlaps(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). .. method:: relate_pattern(other, pattern) Returns true if the elements in the DE-9IM intersection matrix for this geometry and the other matches the given ``pattern`` -- a string of nine characters from the alphabet: {``T``, ``F``, ``*``, ``0``}. .. method:: touches(other) Returns true if the DE-9IM intersection matrix for the two Geometries is ``FT*******``, ``F**T*****`` or ``F***T****``. .. method:: within(other) Returns true if the DE-9IM intersection matrix for the two Geometries is ``T*F**F***``. Topological Methods ------------------- .. method:: buffer(width, quadsegs=8) Returns a geometry that represents all points whose distance from this geometry is less than or equal to the given ``width``. The optional ``quadsegs`` keyword sets the number of segments used to approximate a quarter circle (defaults is 8). .. method:: difference(other) Returns a geometry representing the points making up this geometry that do not make up other. .. method:: intersection(other) Returns a geometry representing the points shared by this geometry and other. .. method:: relate(other) Returns the DE-9IM intersection matrix for this geometry and the other. .. method:: simplify(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. 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. .. method:: sym_difference(other) Returns a set combining the points in this geometry not in other, and the points in other not in this geometry. .. method:: union(other) Returns a Geometry representing all the points in this Geometry and other. Topological Properties ---------------------- .. attribute:: boundary Returns the boundary as a newly allocated Geometry object. .. attribute:: centroid 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). .. attribute:: convex_hull Returns the smallest convex Polygon that contains all the points in the Geometry. .. attribute:: envelope Returns a ``Polygon`` that represents the bounding envelope of this geometry. .. attribute:: point_on_surface Computes and returns a ``Point`` guaranteed to be on the interior of this geometry. Other Properties & Methods -------------------------- .. attribute:: extent This property returns the extent of this geometry as a 4-tuple, consisting of (xmin, ymin, xmax, ymax). .. attribute:: area This property returns the area of the Geometry. .. method:: distance(geom) Returns the distance between the closest points on this Geometry and the given ``geom`` (another ``GEOSGeometry`` object). .. note:: GEOS distance calculations are linear -- in other words, GEOS will not perform a spherical calculation even if the SRID specifies a geographic coordinate system. .. attribute:: length Returns the length of this Geometry (e.g., 0 for point or the circumference of a Polygo .. attribute:: prepared .. versionadded:: 1.1 .. note:: Support for prepared geometries requires GEOS 3.1. Returns a GEOS ``PreparedGeometry`` for the contents of this geometry. ``PreparedGeometry`` objects are optimized for the contains, intersects, and covers operations. Refer to the :ref:`prepared` documentation for more information. .. attribute:: srs Returns an OGR ``SpatialReference`` object corresponding to the SRID of the geometry or ``None``. Consult the `SpatialReference documentation`_ for more information about these objects. .. _SpatialReference documentation: gdal.html#spatialreference .. note:: Requires GDAL. .. method:: transform(ct, clone=False) Transforms the geometry according to the given transformation object, which may be an integer SRID, spatial reference WKT string, a PROJ.4 string, or a ``SpatialReference`` object. By default, the geometry is transformed in-place and nothing is returned. However if the `clone` keyword is set, then the geometry is not modified and a transformed clone is returned instead. .. note:: GDAL and PROJ.4 are required to perform coordinate system transformations. .. rubric:: Footnotes .. [#] *See* Sean Gillies, `Geometries for Python `_ (blog post explaining rationale for abandoning GEOS support); *see also* Sean's message on the `geos-devel mailing list `_, Mar. 5, 2007 .. [#] *See generally* `Python's ctypes documentation `_, at Ch. 14.14. .. [#] Specifically, ``GEOSGeometry`` was introduced in `revision 5008 `_ on April 15th, 2007. .. [#] *See* Sean Gillies, `Geometries for Python Update `_, April 16th 2007 ("The ctypes-based geometry module in r5008 looks kickass. I'm checking it out now."). .. [#] Sean Gillies, `Proposal to launch the Shapely Project `_, May 1, 2007. .. [#] Justin Bronn, `RE: Proposal to launch the Shapely project `_, May 1, 2007. .. [#] Sean Gillies, `Proposal to change Shapely license from LGPL to BSD `_, Nov. 20, 2007. .. [#] *See* `PostGIS EWKB, EWKT and Canonical Forms `_, PostGIS documentation at Ch. 4.1.2.