project logo

GDAL/OGR 快速入门

  本教程只涉及命令行操作。对于数据的可视化,请使用 OSGeo-Live 所提供的各款桌面 GIS 系统,例如 Quantum GIS (QGIS)

  本教程分为 GDAL(栅格数据处理)和 OGR(矢量数据处理)两个部分。第一部分为 GDAL 。

  涉及的内容有:

GDAL
  • 使用 gdalinfo 浏览影象
  • 使用 gdal_translate 转换数据格式
  • 使用 gdalwarp 重投影
  • 使用 gdal_warp 或 gdal_merge.py 拼接影象
  • 使用 gdaltindex 生成作为栅格切片索引的 shp 文件
OGR
  • 使用 ogrinfo 浏览元数据
  • 使用 ogr2ogr 转换数据格式

GDAL

  OSGeo Live 的演示文件位于 /usr/local/share/data 。这里使用 Natural Earth data 来进行演示。请先将其拷贝至 home 目录。即:

cd /home/user
cp -R /usr/local/share/data/natural_earth/ ./gdal_natural_earth

以下位置将出现一个 NaturalEarth 栅格文件和一个 tfw 文件:

cd /home/user/gdal_natural_earth/HYP_50M_SR_W

Tip

建议使用 QGIS 等程序浏览数据。

使用 gdalinfo 获取元数据

gdalinfo HYP_50M_SR_W.tif
  Driver: GTiff/GeoTIFF
  Files: HYP_50M_SR_W.tif
         HYP_50M_SR_W.tfw
  Size is 10800, 5400
  Coordinate System is `'
  Origin = (-179.999999999999972,90.000000000000000)
  Pixel Size = (0.033333333333330,-0.033333333333330)
  Metadata:
    TIFFTAG_SOFTWARE=Adobe Photoshop CS3 Macintosh
    TIFFTAG_DATETIME=2009:09:19 10:13:17
    TIFFTAG_XRESOLUTION=342.85699
    TIFFTAG_YRESOLUTION=342.85699
    TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  Image Structure Metadata:
    INTERLEAVE=PIXEL
  Corner Coordinates:
  Upper Left  (-180.0000000,  90.0000000)
  Lower Left  (-180.0000000, -90.0000000)
  Upper Right ( 180.0000000,  90.0000000)
  Lower Right ( 180.0000000, -90.0000000)
  Center      (  -0.0000000,   0.0000000)
  Band 1 Block=10800x1 Type=Byte, ColorInterp=Red
  Band 2 Block=10800x1 Type=Byte, ColorInterp=Green
  Band 3 Block=10800x1 Type=Byte, ColorInterp=Blue
注:
  • 数据驱动(Driver)为“GTiff/GeoTIFF”
  • 尺寸为 10800x5400
  • 含有 3 个单字节通道
  • 显示了数据坐标信息
  • 没有地理坐标信息

简单数据格式转换

  首先,应当了解所需的数据驱动。使用 –formats 开关可以显示 gdal_translate 所支持的格式。

  各个格式的读写能力记为:
  • 只读(ro)
  • 读写(rw)
  • 读写及更新(rw+)
gdal_translate --formats

  –format 开关等够显示有关数据驱动的许多信息,包括文件创建操作和允许的格式等。

gdalinfo --format jpeg
gdal_translate --format png

转换操作

  使用 gdal_translate 完成转换操作。默认输出是 GeoTIFF:

gdal_translate HYP_50M_SR_W.tif HYP_50M_SR_W.png

  使用 -of 开关控制输出格式,使用 -co 控制文件创建参数(此处为 jpg 文件质量):

gdal_translate -of JPEG -co QUALITY=40 HYP_50M_SR_W.tif HYP_50M_SR_W.jpg

  另外,-ot 用于控制输出数据类别(此处为转换至 16 位的 tif)

gdal_translate -ot Int16 HYP_50M_SR_W.tif HYP_50M_SR_W_Int16.tif

  使用 gdalinfo 验证输出数据的格式信息。

缩放

  使用 -outsize 可以控制输出文件的尺寸。

gdal_translate -outsize 50% 50% HYP_50M_SR_W.tif  HYP_50M_SR_W_small.tif

  使用 gdalinfo 验证尺寸信息。

  使用 -scale 可以控制输入/输出范围和比例。使用 -mm 开关可以显示像元值极限。

  使用 -srcwin 可以通过定位参数(xoff yoff xsize ysize)将影象切分成两部分。使用 -projwin 可定义四角地理坐标(ulx uly lrx lry)。

gdalinfo -mm HYP_50M_SR_W.tif
gdal_translate -srcwin 0 0 5400 5400 HYP_50M_SR_W.tif  west.tif
gdal_translate -srcwin 5400 0 5400 5400 HYP_50M_SR_W.tif  east.tif

使用 gdaltindex 创建栅格切片索引

  可以建立 shp 文件显示栅格切片外框,作为数据索引。对于每幅栅格影象,将生成一个多边形显示其边界,并包含其路径。

gdaltindex index_natural_earth.shp *st.tif

  可以使用 QGIS 和 ogrinfo(相见后述)查看输出的 shp 文件。

../../_images/gdal_gdaltindex13.png
ogrinfo ../HYP_50M_SR_W/ index
INFO: Open of `../HYP_50M_SR_W/'
    using driver `ESRI Shapefile' successful.

Layer name: index
Geometry: Polygon
Feature Count: 2
Extent: (-180.000000, -90.000000) - (180.000000, 90.000000)
Layer SRS WKT: (unknown)
location: String (255.0)
OGRFeature(index):0
  location (String) = east.tif
  POLYGON ((-0.00000000001796 90.0,179.999999999964047 90.0,179.999999999964047 -89.999999999982009,-0.00000000001796 -89.999999999982009,-0.00000000001796 90.0))

OGRFeature(index):1
  location (String) = west.tif
  POLYGON ((-179.999999999999972 90.0,-0.00000000001796 90.0,-0.00000000001796 -89.999999999982009,-179.999999999999972 -89.999999999982009,-179.999999999999972 90.0))

重投影

  这里假设 HYP_50M_SR_W.tif 已被恰当地的创建,具备适当的边界。有前述可见,原始影象没有地理坐标信息,故这里假设使用的是 WGS84 地理坐标。

gdal_translate -a_srs WGS84 HYP_50M_SR_W.tif HYP_50M_SR_W_4326.tif

  使用 gdalwarp 进行投影变换。这里将影象重采样到莫卡脱投影:

gdalwarp -t_srs '+proj=merc +datum=WGS84' HYP_50M_SR_W_4326.tif mercator.tif

  使用 gdalinfo 验证或显示查看。

../../_images/gdal_mercator13.png

  这里再将影象重采样到正交投影:

gdalwarp -t_srs '+proj=ortho +datum=WGS84' HYP_50M_SR_W_4326.tif ortho.tif
../../_images/gdal_ortho13.png

  请注意,地球两极被部分裁切了,这是因为非常靠近两极的边缘是难以重投影的,gdalwarp 放弃了这部分的数据。通过强制 gdalwarp 读取两极数据可以改善这一问题。相见栅格处理教程中的有关内容:http://trac.osgeo.org/gdal/wiki/UserDocs/RasterProcTutorial

影象拼接

  gdal_merge.py 是用于拼接影象的脚本。这里将 east.tif 和 west.tif 拼合成一副影象:

gdal_merge.py  east.tif west.tif -o merged.tif

  拼接也可通过 gdalwarp 实现。它功能更强,但处理速度较慢:

gdalwarp east.tif west.tif warpmerged.tif

OGR

cd /home/usr/gdal_natural_earth/

Tip

可以先用 QGIS 等浏览矢量数据。

使用 ogrinfo 查看元数据

ogrinfo ./natural_earth
INFO: Open of `../natural_earth/'
    using driver `ESRI Shapefile' successful.
1: 10m_lakes (Polygon)
2: 10m_land (Polygon)
3: 10m_rivers_lake_centerlines (Line String)
4: 10m-admin-0-countries (Polygon)
5: 10m_ocean (Polygon)
6: 10m-urban-area (Polygon)
7: 10m_populated_places_simple (Point)

  使用 -so 获取具体文件的信息:

ogrinfo -so ../natural_earth/ 10m-admin-0-countries
INFO: Open of `../natural_earth/'
      using driver `ESRI Shapefile' successful.

Layer name: 10m-admin-0-countries
Geometry: Polygon
Feature Count: 251
Extent: (-179.999783, -89.999828) - (180.000258, 83.633811)
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]]
OBJECTID: Integer (9.0)
COUNTRY: String (100.0)
FEATURECLA: String (32.0)
SOV: String (100.0)
SHAPE_LENG: Real (19.11)
SHAPE_AREA: Real (19.11)

  若不使用以上参数,输出为综述及各个数据集的信息。

ogrinfo ../natural_earth/ 10m-admin-0-countries

  将结果递交给 grep 可过滤结果,例如根据 COUNTRY 字段:

ogrinfo ../natural_earth/ 10m-admin-0-countries | grep COUNTRY

COUNTRY: String (100.0)
COUNTRY (String) = Afghanistan
COUNTRY (String) = Akrotiri Sovereign Base Area
COUNTRY (String) = Aland
COUNTRY (String) = Albania
COUNTRY (String) = Algeria
COUNTRY (String) = American Samoa
COUNTRY (String) = Andorra
etc.

  数据可转换至其它格式。格式信息开关为 –formats

使用 ogr2ogr 转换数据格式

  ogr2ogr 用于转换矢量数据的文件。–formats 开关可用于显示各个格式的读写支持。

  将 countries 文件转换至 GML:

ogr2ogr --formats
ogr2ogr -f GML countries.xml 10m-admin-0-countries.shp

其它示例

  你可以尝试完成以下操作:

  1. 使用 gdalwarp 或 gdal_merge.py 拼接数据。
  2. 使用 gdaladdo 建立缩略图。
  3. QGIS 使用 GDAL/OGR 来支持多种数据格式。同时,它使用 GdalTools 插件进行栅格数据处理。
  4. 使用 ogr2ogr 将栅格数据输入/输出到 PostGIS 。该模块提供了很多选项。
  5. 使用 QGIS 的 OGR-Layer-Konverter 插件。

更多信息

  完成了以上最初的尝试后,以下资源将帮助你进一步学习 GDAL/OGR:

 GDAL 官方网站

 OGR 官方网站:

 GDAL 教程: