==================

GeoDjango Tutorial

==================

Introduction

============

GeoDjango is an add-on for Django that turns it into a world-class geographic

Web framework.  GeoDjango strives to make at as simple as possible to create

geographic Web applications, like location-based services.  Some features include:

* Django model fields for `OGC`_ geometries.

* Extensions to Django's ORM for the querying and manipulation of spatial data.

* Loosely-coupled, high-level Python interfaces for GIS geometry operations and

  data formats.

* Editing of geometry fields inside the admin.

This tutorial assumes a familiarity with Django; thus, if you're brand new to

Django please read through the :doc:`regular tutorial </intro/tutorial01>` to introduce

yourself with basic Django concepts.

.. note::

    GeoDjango has special prerequisites overwhat is required by Django --

    please consult the :ref:`installation documentation <ref-gis-install>`

    for more details.

This tutorial will guide you through the creation of a geographic Web

application for viewing the `world borders`_. [#]_ Some of the code

used in this tutorial is taken from and/or inspired by the `GeoDjango

basic apps`_ project. [#]_

.. note::

    Proceed through the tutorial sections sequentially for step-by-step

    instructions.

.. _OGC: http://www.opengeospatial.org/

.. _world borders: http://thematicmapping.org/downloads/world_borders.php

.. _GeoDjango basic apps: http://code.google.com/p/geodjango-basic-apps/

Setting Up

==========

Create a Spatial Database

-------------------------

.. note::

    MySQL and Oracle users can skip this section because spatial types

    are already built into the database.

First, a spatial database needs to be created for our project.  If using

PostgreSQL and PostGIS, then the following commands will

create the database from a :ref:`spatial database template <spatialdb_template>`::

    $ createdb -T template_postgis geodjango

.. note::

    This command must be issued by a database user that has permissions to

    create a database.  Here is an example set of commands to create such

    a user::

        $ sudo su - postgres

        $ createuser --createdb geo

        $ exit

    Replace ``geo`` to correspond to the system login user name will be

    connecting to the database.  For example, ``johndoe`` if that is the

    system user that will be running GeoDjango.

Users of SQLite and SpatiaLite should consult the instructions on how

to create a :ref:`SpatiaLite database <create_spatialite_db>`.

Create GeoDjango Project

------------------------

Use the ``django-admin.py`` script like normal to create a ``geodjango`` project::

    $ django-admin.py startproject geodjango

With the project initialized, now create a ``world`` Django application within

the ``geodjango`` project::

    $ cd geodjango

    $ python manage.py startapp world

Configure ``settings.py``

-------------------------

The ``geodjango`` project settings are stored in the ``settings.py`` file.  Edit

the database connection settings appropriately::

    DATABASES = {

        'default': {

             'ENGINE': 'django.contrib.gis.db.backends.postgis',

             'NAME': 'geodjango',

             'USER': 'geo',

         }

    }

.. note::

    These database settings are for Django 1.2 and above.

In addition, modify the :setting:`INSTALLED_APPS` setting to include

:mod:`django.contrib.admin`, :mod:`django.contrib.gis`,

and ``world`` (our newly created application)::

    INSTALLED_APPS = (

        'django.contrib.auth',

        'django.contrib.contenttypes',

        'django.contrib.sessions',

        'django.contrib.sites',

        'django.contrib.admin',

        'django.contrib.gis',

        'world'

    )

Geographic Data

===============

.. _worldborders:

World Borders

-------------

The world borders data is available in this `zip file`__.  Create a data directory

in the ``world`` application, download the world borders data, and unzip.

On GNU/Linux platforms the following commands should do it::

    $ mkdir world/data

    $ cd world/data

    $ wget http://thematicmapping.org/downloads/TM_WORLD_BORDERS-0.3.zip

    $ unzip TM_WORLD_BORDERS-0.3.zip

    $ cd ../..

The world borders ZIP file contains a set of data files collectively known as

an `ESRI Shapefile`__, one of the most popular geospatial data formats.  When

unzipped the world borders data set includes files with the following extensions:

* ``.shp``: Holds the vector data for the world borders geometries.

* ``.shx``: Spatial index file for geometries stored in the ``.shp``.

* ``.dbf``: Database file for holding non-geometric attribute data

  (e.g., integer and character fields).

* ``.prj``: Contains the spatial reference information for the geographic

  data stored in the shapefile.

__ http://thematicmapping.org/downloads/TM_WORLD_BORDERS-0.3.zip

__ http://en.wikipedia.org/wiki/Shapefile

Use ``ogrinfo`` to examine spatial data

---------------------------------------

The GDAL ``ogrinfo`` utility is excellent for examining metadata about

shapefiles (or other vector data sources)::

    $ ogrinfo world/data/TM_WORLD_BORDERS-0.3.shp

    INFO: Open of `world/data/TM_WORLD_BORDERS-0.3.shp'

          using driver `ESRI Shapefile' successful.

    1: TM_WORLD_BORDERS-0.3 (Polygon)

Here ``ogrinfo`` is telling us that the shapefile has one layer, and that

layer contains polygon data.  To find out more we'll specify the layer name

and use the ``-so`` option to get only important summary information::

    $ ogrinfo -so world/data/TM_WORLD_BORDERS-0.3.shp TM_WORLD_BORDERS-0.3

    INFO: Open of `world/data/TM_WORLD_BORDERS-0.3.shp'

          using driver `ESRI Shapefile' successful.

    Layer name: TM_WORLD_BORDERS-0.3

    Geometry: Polygon

    Feature Count: 246

    Extent: (-180.000000, -90.000000) - (180.000000, 83.623596)

    Layer SRS WKT:

    GEOGCS["GCS_WGS_1984",

        DATUM["WGS_1984",

            SPHEROID["WGS_1984",6378137.0,298.257223563]],

        PRIMEM["Greenwich",0.0],

        UNIT["Degree",0.0174532925199433]]

    FIPS: String (2.0)

    ISO2: String (2.0)

    ISO3: String (3.0)

    UN: Integer (3.0)

    NAME: String (50.0)

    AREA: Integer (7.0)

    POP2005: Integer (10.0)

    REGION: Integer (3.0)

    SUBREGION: Integer (3.0)

    LON: Real (8.3)

    LAT: Real (7.3)

This detailed summary information tells us the number of features in the layer

(246), the geographical extent, the spatial reference system ("SRS WKT"),

as well as detailed information for each attribute field.  For example,

``FIPS: String (2.0)`` indicates that there's a ``FIPS`` character field

with a maximum length of 2; similarly, ``LON: Real (8.3)`` is a floating-point

field that holds a maximum of 8 digits up to three decimal places.  Although

this information may be found right on the `world borders`_ Web site, this shows

you how to determine this information yourself when such metadata is not

provided.

Geographic Models

=================

Defining a Geographic Model

---------------------------

Now that we've examined our world borders data set using ``ogrinfo``, we can

create a GeoDjango model to represent this data::

    from django.contrib.gis.db import models

    class WorldBorders(models.Model):

        # Regular Django fields corresponding to the attributes in the

	# world borders shapefile.

        name = models.CharField(max_length=50)

        area = models.IntegerField()

        pop2005 = models.IntegerField('Population 2005')

        fips = models.CharField('FIPS Code', max_length=2)

        iso2 = models.CharField('2 Digit ISO', max_length=2)

        iso3 = models.CharField('3 Digit ISO', max_length=3)

        un = models.IntegerField('United Nations Code')

        region = models.IntegerField('Region Code')

        subregion = models.IntegerField('Sub-Region Code')

    	lon = models.FloatField()

    	lat = models.FloatField()

	# GeoDjango-specific: a geometry field (MultiPolygonField), and

        # overriding the default manager with a GeoManager instance.

	mpoly = models.MultiPolygonField()

	objects = models.GeoManager()

        # So the model is pluralized correctly in the admin.

        class Meta:

            verbose_name_plural = "World Borders"

        # Returns the string representation of the model.

        def __unicode__(self):

            return self.name

Two important things to note:

1. The ``models`` module is imported from :mod:`django.contrib.gis.db`.

2. The model overrides its default manager with

   :class:`~django.contrib.gis.db.models.GeoManager`; this is *required*

   to perform spatial queries.

When declaring a geometry field on your model the default spatial reference system

is WGS84 (meaning the `SRID`__ is 4326) -- in other words, the field coordinates are in

longitude/latitude pairs in units of degrees.  If you want the coordinate system to be

different, then SRID of the geometry field may be customized by setting the ``srid``

with an integer corresponding to the coordinate system of your choice.

__ http://en.wikipedia.org/wiki/SRID

Run ``syncdb``

--------------

After you've defined your model, it needs to be synced with the spatial database.

First, let's look at the SQL that will generate the table for the ``WorldBorders``

model::

    $ python manage.py sqlall world

This management command should produce the following output::

    BEGIN;

    CREATE TABLE "world_worldborders" (

        "id" serial NOT NULL PRIMARY KEY,

        "name" varchar(50) NOT NULL,

        "area" integer NOT NULL,

        "pop2005" integer NOT NULL,

        "fips" varchar(2) NOT NULL,

        "iso2" varchar(2) NOT NULL,

        "iso3" varchar(3) NOT NULL,

        "un" integer NOT NULL,

        "region" integer NOT NULL,

        "subregion" integer NOT NULL,

        "lon" double precision NOT NULL,

        "lat" double precision NOT NULL

    )

    ;

    SELECT AddGeometryColumn('world_worldborders', 'mpoly', 4326, 'MULTIPOLYGON', 2);

    ALTER TABLE "world_worldborders" ALTER "mpoly" SET NOT NULL;

    CREATE INDEX "world_worldborders_mpoly_id" ON "world_worldborders" USING GIST ( "mpoly" GIST_GEOMETRY_OPS );

    COMMIT;

If satisfied, you may then create this table in the database by running the

``syncdb`` management command::

    $ python manage.py syncdb

    Creating table world_worldborders

    Installing custom SQL for world.WorldBorders model

The ``syncdb`` command may also prompt you to create an admin user; go ahead and

do so (not required now, may be done at any point in the future using the

``createsuperuser`` management command).

Importing Spatial Data

======================

This section will show you how to take the data from the world borders

shapefile and import it into GeoDjango models using the :ref:`ref-layermapping`.

There are many different different ways to import data in to a

spatial database -- besides the tools included within GeoDjango, you

may also use the following to populate your spatial database:

* `ogr2ogr`_: Command-line utility, included with GDAL, that

  supports loading a multitude of vector data formats into

  the PostGIS, MySQL, and Oracle spatial databases.

* `shp2pgsql`_: This utility is included with PostGIS and only supports

  ESRI shapefiles.

.. _ogr2ogr: http://www.gdal.org/ogr2ogr.html

.. _shp2pgsql: http://postgis.refractions.net/documentation/manual-1.5/ch04.html#shp2pgsql_usage

.. _gdalinterface:

GDAL Interface

--------------

Earlier we used the the ``ogrinfo`` to explore the contents of the world borders

shapefile.  Included within GeoDjango is an interface to GDAL's powerful OGR

library -- in other words, you'll be able explore all the vector data sources

that OGR supports via a Pythonic API.

First, invoke the Django shell::

    $ python manage.py shell

If the :ref:`worldborders` data was downloaded like earlier in the

tutorial, then we can determine the path using Python's built-in

``os`` module::

    >>> import os

    >>> from geodjango import world

    >>> world_shp = os.path.abspath(os.path.join(os.path.dirname(world.__file__),

    ...                             'data/TM_WORLD_BORDERS-0.3.shp'))

Now, the world borders shapefile may be opened using GeoDjango's

:class:`~django.contrib.gis.gdal.DataSource` interface::

    >>> from django.contrib.gis.gdal import *

    >>> ds = DataSource(world_shp)

    >>> print ds

    / ... /geodjango/world/data/TM_WORLD_BORDERS-0.3.shp (ESRI Shapefile)

Data source objects can have different layers of geospatial features; however,

shapefiles are only allowed to have one layer::

    >>> print len(ds)

    1

    >>> lyr = ds[0]

    >>> print lyr

    TM_WORLD_BORDERS-0.3

You can see what the geometry type of the layer is and how many features it

contains::

    >>> print lyr.geom_type

    Polygon

    >>> print len(lyr)

    246

.. note::

    Unfortunately the shapefile data format does not allow for greater

    specificity with regards to geometry types.  This shapefile, like

    many others, actually includes ``MultiPolygon`` geometries in its

    features.  You need to watch out for this when creating your models

    as a GeoDjango ``PolygonField`` will not accept a ``MultiPolygon``

    type geometry -- thus a ``MultiPolygonField`` is used in our model's

    definition instead.

The :class:`~django.contrib.gis.gdal.Layer` may also have a spatial reference

system associated with it -- if it does, the ``srs`` attribute will return a

:class:`~django.contrib.gis.gdal.SpatialReference` object::

    >>> srs = lyr.srs

    >>> print srs

    GEOGCS["GCS_WGS_1984",

        DATUM["WGS_1984",

            SPHEROID["WGS_1984",6378137.0,298.257223563]],

        PRIMEM["Greenwich",0.0],

        UNIT["Degree",0.0174532925199433]]

    >>> srs.proj4 # PROJ.4 representation

    '+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs '

Here we've noticed that the shapefile is in the popular WGS84 spatial reference

system -- in other words, the data uses units of degrees longitude and latitude.

In addition, shapefiles also support attribute fields that may contain

additional data.  Here are the fields on the World Borders layer:

    >>> print lyr.fields

    ['FIPS', 'ISO2', 'ISO3', 'UN', 'NAME', 'AREA', 'POP2005', 'REGION', 'SUBREGION', 'LON', 'LAT']

Here we are examining the OGR types (e.g., whether a field is an integer or

a string) associated with each of the fields:

    >>> [fld.__name__ for fld in lyr.field_types]

    ['OFTString', 'OFTString', 'OFTString', 'OFTInteger', 'OFTString', 'OFTInteger', 'OFTInteger', 'OFTInteger', 'OFTInteger', 'OFTReal', 'OFTReal']

You can iterate over each feature in the layer and extract information from both

the feature's geometry (accessed via the ``geom`` attribute) as well as the

feature's attribute fields (whose **values** are accessed via ``get()``

method)::

    >>> for feat in lyr:

    ...    print feat.get('NAME'), feat.geom.num_points

    ...

    Guernsey 18

    Jersey 26

    South Georgia South Sandwich Islands 338

    Taiwan 363

:class:`~django.contrib.gis.gdal.Layer` objects may be sliced::

    >>> lyr[0:2]

    [<django.contrib.gis.gdal.feature.Feature object at 0x2f47690>, <django.contrib.gis.gdal.feature.Feature object at 0x2f47650>]

And individual features may be retrieved by their feature ID::

    >>> feat = lyr[234]

    >>> print feat.get('NAME')

    San Marino

Here the boundary geometry for San Marino is extracted and looking

exported to WKT and GeoJSON::

    >>> geom = feat.geom

    >>> print geom.wkt

    POLYGON ((12.415798 43.957954,12.450554 ...

    >>> print geom.json

    { "type": "Polygon", "coordinates": [ [ [ 12.415798, 43.957954 ], [ 12.450554, 43.979721 ], ...

``LayerMapping``

----------------

We're going to dive right in -- create a file called ``load.py`` inside the

``world`` application, and insert the following::

    import os

    from django.contrib.gis.utils import LayerMapping

    from models import WorldBorders

    world_mapping = {

        'fips' : 'FIPS',

        'iso2' : 'ISO2',

        'iso3' : 'ISO3',

        'un' : 'UN',

        'name' : 'NAME',

        'area' : 'AREA',

        'pop2005' : 'POP2005',

        'region' : 'REGION',

        'subregion' : 'SUBREGION',

        'lon' : 'LON',

        'lat' : 'LAT',

        'mpoly' : 'MULTIPOLYGON',

    }

    world_shp = os.path.abspath(os.path.join(os.path.dirname(__file__), 'data/TM_WORLD_BORDERS-0.3.shp'))

    def run(verbose=True):

        lm = LayerMapping(WorldBorders, world_shp, world_mapping,

                          transform=False, encoding='iso-8859-1')

        lm.save(strict=True, verbose=verbose)

A few notes about what's going on:

* Each key in the ``world_mapping`` dictionary corresponds to a field in the

  ``WorldBorders`` model, and the value is the name of the shapefile field

  that data will be loaded from.

* The key ``mpoly`` for the geometry field is ``MULTIPOLYGON``, the

  geometry type we wish to import as.  Even if simple polygons are encountered

  in the shapefile they will automatically be converted into collections prior

  to insertion into the database.

* The path to the shapefile is not absolute -- in other words, if you move the

  ``world`` application (with ``data`` subdirectory) to a different location,

  then the script will still work.

* The ``transform`` keyword is set to ``False`` because the data in the

  shapefile does not need to be converted -- it's already in WGS84 (SRID=4326).

* The ``encoding`` keyword is set to the character encoding of string values in

  the shapefile. This ensures that string values are read and saved correctly

  from their original encoding system.

Afterwards, invoke the Django shell from the ``geodjango`` project directory::

   $ python manage.py shell

Next, import the ``load`` module, call the ``run`` routine, and watch ``LayerMapping``

do the work::

   >>> from world import load

   >>> load.run()

.. _ogrinspect-intro:

Try ``ogrinspect``

------------------

Now that you've seen how to define geographic models and import data with the

:ref:`ref-layermapping`, it's possible to further automate this process with

use of the :djadmin:`ogrinspect` management command.  The :djadmin:`ogrinspect`

command  introspects a GDAL-supported vector data source (e.g., a shapefile) and

generates a model definition and ``LayerMapping`` dictionary automatically.

The general usage of the command goes as follows::

    $ python manage.py ogrinspect [options] <data_source> <model_name> [options]

Where ``data_source`` is the path to the GDAL-supported data source and

``model_name`` is the name to use for the model.  Command-line options may

be used to further define how the model is generated.

For example, the following command nearly reproduces the ``WorldBorders`` model

and mapping dictionary created above, automatically::

    $ python manage.py ogrinspect world/data/TM_WORLD_BORDERS-0.3.shp WorldBorders --srid=4326 --mapping --multi

A few notes about the command-line options given above:

* The ``--srid=4326`` option sets the SRID for the geographic field.

* The ``--mapping`` option tells ``ogrinspect`` to also generate a

  mapping dictionary for use with :class:`~django.contrib.gis.utils.LayerMapping`.