修改

下面是做到这一点的适当方式,而文档

import random
from osgeo import gdal, ogr    

RASTERIZE_COLOR_FIELD = "__color__"

def rasterize(pixel_size=25)
    # Open the data source
    orig_data_source = ogr.Open("test.shp")
    # Make a copy of the layer's data source because we'll need to 
    # modify its attributes table
    source_ds = ogr.GetDriverByName("Memory").CopyDataSource(
            orig_data_source, "")
    source_layer = source_ds.GetLayer(0)
    source_srs = source_layer.GetSpatialRef()
    x_min, x_max, y_min, y_max = source_layer.GetExtent()
    # Create a field in the source layer to hold the features colors
    field_def = ogr.FieldDefn(RASTERIZE_COLOR_FIELD, ogr.OFTReal)
    source_layer.CreateField(field_def)
    source_layer_def = source_layer.GetLayerDefn()
    field_index = source_layer_def.GetFieldIndex(RASTERIZE_COLOR_FIELD)
    # Generate random values for the color field (it's here that the value
    # of the attribute should be used, but you get the idea)
    for feature in source_layer:
        feature.SetField(field_index, random.randint(0, 255))
        source_layer.SetFeature(feature)
    # Create the destination data source
    x_res = int((x_max - x_min) / pixel_size)
    y_res = int((y_max - y_min) / pixel_size)
    target_ds = gdal.GetDriverByName('GTiff').Create('test.tif', x_res,
            y_res, 3, gdal.GDT_Byte)
    target_ds.SetGeoTransform((
            x_min, pixel_size, 0,
            y_max, 0, -pixel_size,
        ))
    if source_srs:
        # Make the target raster have the same projection as the source
        target_ds.SetProjection(source_srs.ExportToWkt())
    else:
        # Source has no projection (needs GDAL >= 1.7.0 to work)
        target_ds.SetProjection('LOCAL_CS["arbitrary"]')
    # Rasterize
    err = gdal.RasterizeLayer(target_ds, (3, 2, 1), source_layer,
            burn_values=(0, 0, 0),
            options=["ATTRIBUTE=%s" % RASTERIZE_COLOR_FIELD])
    if err != 0:
        raise Exception("error rasterizing layer: %s" % err)

<强>原始问题

我在寻找关于如何使用osgeo.gdal.RasterizeLayer()(文档字符串是非常简洁,并在C或C ++的API文档我无法找到它的信息,我只发现了的 Java绑定的)。

我适于一个单元测试并试图在由多边形的的.shp:

import os
import sys
from osgeo import gdal, gdalconst, ogr, osr

def rasterize():
    # Create a raster to rasterize into.
    target_ds = gdal.GetDriverByName('GTiff').Create('test.tif', 1280, 1024, 3,
            gdal.GDT_Byte)
    # Create a layer to rasterize from.
    cutline_ds = ogr.Open("data.shp")
    # Run the algorithm.
    err = gdal.RasterizeLayer(target_ds, [3,2,1], cutline_ds.GetLayer(0),
            burn_values=[200,220,240])
    if err != 0:
        print("error:", err)

if __name__ == '__main__':
    rasterize()

它运行良好,但所有我得到为黑色的.tif。

什么是burn_values参数?可以RasterizeLayer()被用于光栅化基于属性的值与着色不同特征的层?

如果不能,我应该用什么?是 AGG 适于使地理数据(我想要的没有反锯齿和一个非常坚固的渲染器,能够正确地绘制非常大和非常小的特征,有可能从“脏数据”(退化多边形等),有时在大坐标指定)?

在这里,多边形是由一个属性的值区分(颜色并不重要,我只想有一个不同的属性的每个值)。

有帮助吗?

解决方案

编辑:我想我会使用QGIS python绑定: http://www.qgis.org/维基/ Python_Bindings

这是我能想到的最简单的方法。我记得前手滚动的东西,但它的丑陋。 QGIS会更容易,即使你不得不做出一个单独的Windows安装(获得蟒与它的工作),然后建立一个XML-RPC服务器在一个单独的Python程序运行它。

我就可以得到GDAL光栅化妥善这是伟大的。

我没有用GDAL了一段时间,但这里是我的猜测:

burn_values是假的颜色,如果你不使用Z值。您的多边形内部的一切都是[255,0,0](红色)如果你使用burn=[1,2,3],burn_values=[255,0,0]。我不知道会发生什么点 - 他们可能不积

。 如果您想使用Z值

使用gdal.RasterizeLayer(ds,bands,layer,burn_values, options = ["BURN_VALUE_FROM=Z"])

我只是从你在看测试这个拉http://svn.osgeo.org/gdal/trunk/autotest/alg/rasterize.py

另一种方法 - 拉多边形对象出来,并使用匀称,这可能不是有吸引力的绘制它们。或可考虑GeoDjango内置(我认为它使用的OpenLayers绘制成使用JavaScript的浏览器)。

此外,你需要栅格化?一个PDF的出口可能会更好,如果你真的想要的精度。

其实,我觉得我找到了使用Matplotlib(提取和突出的特征后)比光栅化更容易,我能得到更多的控制。

编辑:

一个较低的水平的方法是在这里:

http://svn.osgeo.org/ GDAL /中继/ GDAL /痛饮/蟒/样品/ gdal2grd.py \

最后,可以遍历多边形(它们转化成一个局部投影之后),并直接绘制它们。但是你最好不要有复杂的多边形,否则你将有一点悲伤。如果您有复杂的多边形......你可能是最好关闭用匀称和R树从 HTTP://trac.gispython如果你想推出自己的绘图仪.ORG /实验室

GeoDjango内置可能是一个好地方,问..他们会认识很多比我多。他们有一个邮件列表?还有很多蟒蛇测绘专家身边,但他们都不担心这个问题。我想他们只是绘制它在QGIS或草或东西。

说真的,我希望有人谁知道他们在做什么回答。

许可以下: CC-BY-SA归因
不隶属于 StackOverflow
scroll top