GeoDjango教程[1]

**引子**

最初接触GeoDjango应该是很早的事情了，一直希望将其完整的学下来
期间也曾向公司申请了一个小项目作为学习平台，但左右还是一些阻碍与懒惰，终究没有
完成
http://www.geodjango.org

说来惭愧，进公司管了三年多的IT管理，竟连一段码都没有写过（虽说也不需要写，但纸
上谈兵的管理终究有点心虚，最近工作之中受累技术短板的感觉也越发强烈，也希望通过
GeoDjango能为自己找些补强的动力）

大约是去年在某师弟的博客里看到部分GeoDjango教程翻译时，就有一股冲动想通过翻译教
程的办法督促自己完成这一拖了很久的工作
http://www.geoinformatics.cn/2009/07/geodjango中文教程/

最近又一次偶然看到，既然不积跬步无以至千里，不妨集腋成裘便利于你我

ps：先从已有的译稿开始，比较v1.3不同之处略作更新



**简介**

GeoDjango是一款支持地理数据存储和空间查询的Django衍生项目。开发者力图使
GeoDjango做到尽可能简单的创建地理网络应用，比如基于位置的服务。GeoDjango的一些
特征包括：支持OGC地物类型的数据模型，面向空间数据查询与操作的Django ORM增强，针
对地物分析与数据格式的松耦合、高效python接口，直接通过Admin界面的地物数据编辑功
能等

这份教程假定你熟知Django，如果你是刚刚接触Django，那么请您先阅读“Django教程”
以获得Django的一些基础知识。

同时，GeoDjango运行除了Django之外，仍需要一些其他条件——请参考“GeoDjango安装
文档”以获详情。

这份教程将利用“世界地图”这个GIS应用来指导你的学习。这份教程中部分代码受到了或
者直接来自“GeoDjango 基础实例应用”这个项目的启发。

后续部分建议按序阅读。



**配置**

*建立空间数据库*

注：MySQL和Oracle数据库用户可以略过这一部分，因为这两种数据库在创建时已有空间特
征。

首先，项目需要创建一个空间数据库。如果利用PostgreSQL和PostGIS，则可利用如下的命
令从空间数据库模板(spatial database template) 来建立数据库:

$ createdb -T template_postgis geodjango

注：请确保用户有权限创建数据库，如下命令可创建授权用户

$ sudo su - postgres
$ createuser --createdb geo
$ exit

选项 geo 指定 createdb 利用 geo 用户作为数据库的拥有者。当然您也可替换使用其他 用户名。 此外，使用SQLite和SpatiaLite的用户请使用相关文档创建上述过程。 *创建GeoDjango项目* 利用 django-admin.py 脚本来创建一个 geodjango 项目:

$ django-admin.py startproject geodjango

当 geodjango 项目初始化后，在里面建立 world Django 实例:

$ cd geodjango
$ python manage.py startapp world

*配置 settings.py* geodjango 的项目设置存储在文件 settings.py 中，编辑数据库连接:

DATABASE_ENGINE ='django.contrib.gis.db.backends.postgis'
DATABASE_NAME ='geodjango'
DATABASE_USER ='geo'

注：上述配置仅适用于1.2版本以上的Django。 另外，编辑 INSTALLED_APPS 设置，使之包含 django.contrib.admin 、 django. contrib.gis 、以及 geodjango.world (刚刚新建的项目):

INSTALLED_APPS =(
    'django.contrib.auth',
    'django.contrib.contenttypes',
    'django.contrib.sessions',
    'django.contrib.sites',
    'django.contrib.admin',
    'django.contrib.gis',
    'geodjango.world'
)

**地理数据** *世界地图* 世界地图数据保存于 zip 压缩文件 . 在 world 实例中创建一个数据目录，下载并解压世 界地图数据:

$ 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 ../..

这份世界地图zip压缩文件包含了一整套文件，诸如 ESRI Shapefile ，它是众所周知的地 理数据格式。解压过后，您会发现还包括: .shp: 存储世界地图边界的矢量几何数据。 .shx: 是存储在 .shp 中的几何数据的空间索引文件。 .dbf: 用来存储和几何数据无关的属性数据的数据库文件 (例如整形和字符型的字段)。 .prj: 包含了shapefile文件中空间参考系统的信息。 *利用 ogrinfo 检验空间数据* GDAL的 ogrinfo 是用来检验shapefile或者其他矢量数据的元数据 (metadata)的工具:

   $ 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)

``ogrinfo`` 提示说shapefile只有一个图层，并且该图层含有多边形数据。为了获得更多 信息，制定图层名称，并利用选项 ``-so`` 来获取概述信息::

   $ 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)

这份详细的概述信息表明图层中地理要素的数量为246，以及其地理几何边界，空间参考系 统(“SRS WKT”) 同每一个属性字段的详细信息。例如， FIPS: String (2.0) 表明有一 个名为 FIPS 的字符串字段，最大字符长度为2，类似的 LON: Real (8.3) 是一个浮点型 字段，最大可容纳8个字节以及最多3个浮点位。尽管这些信息可以很容易的在 世界地图 的网站上找到，但是该命令可以在没有提供元数据时获取这些信息。 **地理模型** *定义地理模型* 现在利用 ogrinfo 来检验世界地图数据，创建一个GeoDjango模型来表达数据:

from django.contrib.gis.db import models

classWorldBorders(models.Model):
    # 基于世界地图shapefile文件中的属性字段而构建的Django字段
    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制定的一个地理几何字段 (MultiPolygonField),
    # 并用GeoManager实例覆盖默认的管理器。
    mpoly = models.MultiPolygonField()
    objects = models.GeoManager()

    classMeta:
        # 制定复数名.
        verbose_name_plural ="World Borders"

    # 返回表达模型的字符串。
    def __unicode__(self):
        return self.name

两个需要注意的问题： models 模块必须由 django.contrib.gis.db 导入。 模型必须用 GeoManager ; 这样做是 必须的 ，因为 GeoManager 是用来进行空间查询的 。 用户在模型中声明的地理几何字段，其默认的空间参考系统是WGS84 (其 SRID 是 4326) – 换句话说，经度和纬度是用”度”做单位。如果你想用其他的坐标系统，你可以制定地 理几何字段的 srid 。 *运行 syncdb* 定义好模型之后，需要同步空间数据库。首先，让我们看看 WorldBorders 模型产生的SQ L语句:

$ python manage.py sqlall world

这个命令将输出如下的SQL语句:

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;

如果这段SQL语句符合要求，你可以通过执行 syncdb 命令来建立数据表:

$ python manage.py syncdb
Creating table world_worldborders
Installing custom SQL for world.WorldBorders model

利用 syncdb 命令可以方便的创建一个管理员用户；如果感兴趣，可以参考 createsuperuser 命令。 *导入空间数据* 这一节将向您展示如何利用 LayerMapping 工具将世界地图的shapefile文件导入 GeoDajngo模型。导入数据的方式有很多种，除去使用刚才提到的方法，您还可以利用: ogr2ogr: 该命令行工具属于GDAL，支持将多种矢量格式的数据导入 PostGIS、 MySQL、 以及Oracle 空间数据库。 shp2pgsql: 这个工具包含在PostGIS中，但是只支持ESRI Shapefile文件。 GDAL 接口 前文中谈到利用 ogrinfo 来获取shapefile文件的内容。GeoDjango内置了GDAL强大的OGR 库 – 换句话说，您可以利用支持OGR的Python API来处理矢量数据。 首先，开启Django shell:

$ python manage.py shell

如果 世界地图 数据像前文所述的那样下载后，便可以用python的内置模块 os 来确定其 路径:

>>>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'))

现在，世界地图shapefile文件会被GeoDjango的 DataSource 接口打开:

>>>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) 可以包含多个空间数据图层，然后，shapefile文件只允许有 一个图层:

>>>print len(ds)
1
>>> lyr = ds[0]
>>>print lyr
TM_WORLD_BORDERS-0.3

您同样可看到该图层的地理类型以及所包含的地理要数个数:

>>>print lyr.geom_type
Polygon
>>>print len(lyr)
246

注: shapefile并不允许对地理几何类型的定义。这份shapefile文件和大多数 shapefile文件一样实际上包含的是 MultiPolygon 作为其地理几何要素类型。您需要注意 的是 GeoDajngo的 PolygonField 字段并不会接受 MultiPolygon 类型的数据 – 这也是 在模型中使用 MultiPolygonField 的原因。 Layer 也可以指定空间参考系对象。 – 如果有的话, 则 srs 属性将返回 SpatialReference 对象:

>>> 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 '

在此，我们注意到shapfile文件使用的是流行的WGS84空间参考系统 – 也就是说利用度作 为单位表达经纬度。 另外，shapefile文件中属性数据也可能含有数据。在世界地图图层中，将包含如下属性数 据

>>>print lyr.fields
['FIPS','ISO2','ISO3','UN','NAME','AREA','POP2005','REGION','SUBREGION
','LON','LAT']

检测每个字段的OGR字段类型 (例如，一个字段是整形还是字符串):

>>>[fld.__name__ for fld in lyr.field_types]
['OFTString','OFTString','OFTString','OFTInteger','OFTString','OFTInteger
','OFTInteger','OFTInteger','OFTInteger','OFTReal','OFTReal']

您可以读取图层中的每一个要素并提取其地理几何字段 (通过 geom 字段获得) 以及其属 性字段 (该字段 值 可通过 get() 方法) 读取:

>>>for feat in lyr:
...    print feat.get('NAME'), feat.geom.num_points
...
Guernsey18
Jersey26
SouthGeorgiaSouthSandwichIslands338
Taiwan363

Layer objects may be sliced:

>>> lyr[0:2]
[<django.contrib.gis.gdal.feature.Featureobject at 0x2f47690>,<django.
contrib.gis.gdal.feature.Featureobject at 0x2f47650>]

每个不同的要素可以通过其要素ID获取:

>>> feat = lyr[234]
>>>print feat.get('NAME')
San Marino

例如，San Marino的地理几何信息被提取并解压为WKT格式或者GeoJSON格式:

>>> geom = feat.geom
>>>print geom.wkt
POLYGON ((12.41579843.957954,12.450554...
>>>print geom.json
{"type":"Polygon","coordinates":[[[12.415798,43.957954],[12.450554,
 43.979721],...

*LayerMapping* 我们将进一步在 world 实例中创建一个python脚本文件 load.py ，加入如下语句:

import os
from django.contrib.gis.utils importLayerMapping
from models importWorldBorders

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_W
ORLD_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)

解释: 每一个键(key)是 world_mapping 字典类型中所对应的 WorldBorders 模型的字段， 而值 (value)是模型读取的shapefile文件中属性字段的名称。 键 mpoly 代表的地理几何类型是 MULTIPOLYGON ，这也是我们想输入的数据类型。即使在 shapefile文件中遇到了简单多边形也将自动转换为集合来插入数据库。 shapefile文件的路径不是绝对路径 – 也就是说，如果你将 world 应用实例 (包含着 data 子目录) 移动到另一个地址, 这份脚本依然会有效工作。 transform 关键字被设置为 否(False) ，因为shapefile文件的数据不需要转换 – 它已 基于WGS84坐标系统了 (SRID=4326)。 encoding 关键字将设置为shapefile文件的字符编码类型。这将保证字符串可读并正确的 保存。 紧接着，从 geodjango 项目目录运行Django shell:

$ python manage.py shell

然后，导入 load 模块，运行 run 让 LayerMapping 做剩下的工作:

>>>from world import load
>>> load.run()

*ogrinspect 检测* 至此，您已经看到如何定义地理模型和如何利用 LayerMapping 工具导入shapefile数据, 其实您还可以进一步利用 ogrinspect 管理命令来自动化您的工作. ogrinspect 命令将 支持自省某个GDAL支持的矢量数据源(例如shapefile文件)并自动产生一个模型定义和 LayerMapping 字典数据。 基本使用方法如下:

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

data_source 指示了GDAL支持的矢量数据的路径，``model_name`` 是模型的名称。 options将用于进一步定义模型。 举例而言，如下的命令将自动产生 WorldBorders 模型和先前的映射字典数据:

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

解释: --srid=4326 选项设置地理几何字段的 SRID。 --mapping 选项指明 ogrinspect 要同时产生给 LayerMapping 使用的映射字典数据. --multi 选项指明地理数据类型是 MultiPolygonField 而不是 PolygonField 。 该命令将输出如下的语句，该语句将被复制到GeoDjango应用实例的 models.py 脚本文件 中:

# 这是一份由 ogrinspect 自动产生的 Django 模型文件。
from django.contrib.gis.db import models

classWorldBorders(models.Model):
    fips = models.CharField(max_length=2)
    iso2 = models.CharField(max_length=2)
    iso3 = models.CharField(max_length=3)
    un = models.IntegerField()
    name = models.CharField(max_length=50)
    area = models.IntegerField()
    pop2005 = models.IntegerField()
    region = models.IntegerField()
    subregion = models.IntegerField()
    lon = models.FloatField()
    lat = models.FloatField()
    geom = models.MultiPolygonField(srid=4326)
    objects = models.GeoManager()

# 为 WorldBorders 模型自动产生的 `LayerMapping` 字典数据
worldborders_mapping ={
    'fips':'FIPS',
    'iso2':'ISO2',
    'iso3':'ISO3',
    'un':'UN',
    'name':'NAME',
    'area':'AREA',
    'pop2005':'POP2005',
    'region':'REGION',
    'subregion':'SUBREGION',
    'lon':'LON',
    'lat':'LAT',
    'geom':'MULTIPOLYGON',
}

**空间查询** *空间查找* GeoDjango 扩展了Django 的ORM模型，并允许使用空间查找。 让我们做一个实验，看看 WorldBorder 模型是否包含某个兴趣点。首先，激活管理shell:

$ python manage.py shell

定义一个兴趣点:

>>> pnt_wkt ='POINT(-95.3385 29.7245)'

pnt_wkt 字符串代表一个点在经度 -95.3385 纬度 29.7245. 地理数据类型定义为 Well Known Text (WKT), 它是 Open Geospatial Consortium (OGC)所推荐的地理数据类型。 导入 WorldBorders 模型, 并利用 pnt_wkt 作为参数执行 contains 查找:

>>>from world.models importWorldBorders
>>> qs =WorldBorders.objects.filter(mpoly__contains=pnt_wkt)
>>> qs
[<WorldBorders:UnitedStates>]

这里我们将获得一个 GeoQuerySet 对象，该对象只包含一个要素，这也正是我们所需要寻 找的。 同样的，一个 GEOS 地理对象 也可以能被使用 – 比如用 intersects 空间查找 来获取 WorldBorders 对象的实例 San Marino ，该空间查找是和 get() 方法一同使用:

>>>from django.contrib.gis.geos importPoint
>>> pnt =Point(12.4604,43.9420)
>>> sm =WorldBorders.objects.get(mpoly__intersects=pnt)
>>> sm
<WorldBorders:SanMarino>

空间查找 contains 和 intersects 属于 GeoDjango 数据库接口 ，查看文档获得更多详 细资料。 *空间投影转换* 当查询空间数据时，如果地理数据拥有不同的坐标投影，GeoDjango将自动转换。在接下来 的的例子中， 坐标系统将表达为 EPSG SRID 32140 ， 地图投影坐标系统 仅 在南德州使 用，其单位为米而非度。

>>>from django.contrib.gis.geos import*
>>> pnt =Point(954158.1,4215137.1, srid=32140)

注意 pnt 对象也可以使用EWKT进行构建，”extended” 代表WKT可以包含SRID:

>>> pnt =GEOSGeometry('SRID=32140;POINT(954158.1 4215137.1)')

当使用 GeoDjango 的 ORM 模型时，将自动把地理要素值转化为SQL，这样可以使开发者在 更高层次上进行工作:

>>> qs =WorldBorders.objects.filter(mpoly__intersects=pnt)
>>> qs.query.as_sql()# Generating the SQL
('SELECT "world_worldborders"."id", "world_worldborders"."name", "world_
worldborders"."area",
"world_worldborders"."pop2005", "world_worldborders"."fips", "world_
worldborders"."iso2",
"world_worldborders"."iso3", "world_worldborders"."un", "world_worldborders"."
region",
"world_worldborders"."subregion", "world_worldborders"."lon", "world_
worldborders"."lat",
"world_worldborders"."mpoly" FROM "world_worldborders"
WHERE ST_Intersects("world_worldborders"."mpoly", ST_Transform(%s, 4326))',
(<django.contrib.gis.db.backend.postgis.adaptor.PostGISAdaptorobject at 0x256
41b0>,))
>>> qs # printing evaluates the queryset
[<WorldBorders:UnitedStates>]

*怠惰地理数据* GeoDjango可以将对地理数据进行标准化的表达。除了访问地理字段，GeoDjango可以创建 GEOS 地理对象 ， 并给地理对象强大的功能，例如序列化地理数据:

>>> sm =WorldBorders.objects.get(name='San Marino')
>>> sm.mpoly
<MultiPolygonobject at 0x24c6798>
>>> sm.mpoly.wkt # WKT
MULTIPOLYGON (((12.415798000000000643.9579540000000009,12.45055400000000034
3.9797209999999978,...
>>> sm.mpoly.wkb # WKB (as Python binary buffer)
<read-only buffer for0x1fe2c70, size -1, offset 0 at 0x2564c40>
>>> sm.mpoly.geojson # GeoJSON (requires GDAL)
'{ "type": "MultiPolygon", "coordinates": [ [ [ [ 12.415798, 43.957954 ], [ 12
.450554, 43.979721 ], ...

包括所有GEOS库提供的高级地理几何操作:

>>> pnt =Point(12.4604,43.9420)
>>> sm.mpoly.contains(pnt)
True
>>> pnt.contains(sm.mpoly)
False

*GeoQuerySet 方法* **在地图上表现数据** *谷歌地图(Google Map)* *地理数据管理(Geographic Admin)* 基础知识 GeoDjango 提供 Django Admin 以使用户通过 OpenLayers 构建并修改地理要素。 让我们再进一步，在 world 应用实例中创建 admin.py 文件，加入如下代码:

from django.contrib.gis import admin
from models importWorldBorders

admin.site.register(WorldBorders, admin.GeoModelAdmin)

然后，编辑 geodjango 项目目录下的 urls.py 文件:

from django.conf.urls.defaults import*
from django.contrib.gis import admin

admin.autodiscover()

urlpatterns = patterns('',
                       (r'^admin/(.*)', admin.site.root),
                       )

启动 Django 开发服务器:

$ python manage.py runserver

最后，您能通过 http://localhost:8000/admin/ 进行浏览，通过用 sync 所生成的管理 用户名进行登录。您可以浏览任何一个 WorldBorders 实例 – 您可以点选任何国家的边 界甚至拖动节点。 *OSMGeoAdmin* 通过 OSMGeoAdmin ，GeoDjango 将在管理模块中使用 Open Street Map 图层。这将提供 更多的内容 (包括道路和街道详情) 而不是 OpenLayers 默认提供的 Vector Map Level 0 WMS服务 (由 Metacarta 支持)。 请注意: OSMGeoAdmin 需要在 PostgreSQL/PostGIS 数据库中加入 球面莫卡托投影 到 spatial_ ref_sys 数据表中。 必须安装 PROJ.4 转换数据 (详见 PROJ.4 安装说明 )。 万事俱备则在 admin.py 文件中用如下代码替换 OSMGeoAdmin:

admin.site.register(WorldBorders, admin.OSMGeoAdmin)