GDAL和OGR Python绑定中的Python Gotchas

本页列出了GDAL和OGR的Python绑定的一些方面,这些方面可能会让Python程序员大吃一惊。如果你发现了新的东西,可以随意添加到列表中,但是可以考虑在 gdal-dev mailing list 首先,要确保您完全理解这个问题,并且其他人都同意它是意外的,“非Pythonic”,或者会让许多Python程序员大吃一惊的东西。确保引用电子邮件线程、GitHub票证和其他附加信息源。

这个列表不是报告错误的地方。如果你相信有东西是虫子,请 open a ticket 然后将问题报告给gdal-dev,如果它是与Python相关的,那么可以考虑在这里列出它。如果它通常与GDAL或OGR相关,而不是特定的Python绑定,请不要在这里列出它。

并不是这里列出的所有项目都是bug。其中一些就是GDAL和OGR的工作方式,如果不破坏现有代码,就很难修复。如果你不喜欢某件事情的工作方式,认为它应该改变,可以在gdal dev上自由讨论,看看能做些什么。

是故意的。。。或按历史记录

这些是不可预料的行为,GDAL和OGR团队不认为它们是bug,也不太可能由于所需的工作而更改,或者其修复可能会影响向后兼容性等。

Python绑定不会引发异常,除非显式调用 UseExceptions()

默认情况下,当发生错误时,GDAL和OGR Python绑定不会引发异常。相反,它们返回一个错误值,如 None 并将错误消息写入 sys.stdout . 例如,当尝试使用GDAL打开不存在的数据集时:

>>> from osgeo import gdal
>>> gdal.Open('C:\\foo.img')
ERROR 4: 'C:\foo.img does not exist in the file system,
and is not recognized as a supported dataset name.

>>>

在Python中,传统的方法是通过引发异常来报告错误。通过调用 UseExceptions() 功能:

>>> from osgeo import gdal
>>> gdal.UseExceptions()    # Enable exceptions
>>> gdal.open('C:\\foo.img')
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
RuntimeError: 'C:\foo.img' does not exist in the file system,
and is not recognized as a supported dataset name.

>>>

GDAL团队承认Python程序员希望在默认情况下启用异常,但表示在默认情况下禁用异常是为了 preserve backward compatibility .

如果在删除与之有关系的对象后使用该对象,则Python将崩溃

举个例子:

>>> from osgeo import gdal
>>> dataset = gdal.Open('C:\\RandomData.img')
>>> band = dataset.GetRasterBand(1)
>>> print(band.Checksum())
31212

在这个例子中, band 与…有关系 dataset 这需要 dataset 为了 band 去工作。如果我们删除 dataset 然后尝试使用 band ,Python将崩溃:

>>> from osgeo import gdal
>>> dataset = gdal.Open('C:\\RandomData.img')
>>> band = dataset.GetRasterBand(1)
>>> del dataset           # This will cause the Python garbage collector to deallocate dataset
>>> band.GetChecksum()    # This will now crash Python because the band's dataset is gone
< Python crashes >

这个问题可以用微妙的方式表现出来。例如,如果尝试在一行代码中实例化临时数据集实例,则可能会发生这种情况:

>>> from osgeo import gdal
>>> print(gdal.Open('C:\\RandomData.img').GetRasterBand(1).Checksum())
< Python crashes >

在本例中,在调用 GetRasterBand() 所以Python取消了它 之前 打电话 Checksum() .

>>> from osgeo import gdal
>>>
>>> def load_band(data_filename):
>>>     dataset = gdal.Open(data_filename)
>>>     return dataset.GetRasterBand(1)
>>>
>>> band = load_band('c:\\RandomData.img')
>>> print(band.Checksum())
< Python crashes >

this example is the same case with above but it looks different. the dataset object is only available in the load_band function and it'll be deleted right after leaving the function.

这个问题是因为GDAL和OGR对象是用C++实现的,它们之间的关系用指针来保持在C++中。当您在Python中删除DataSet实例时,它会导致C++对象在其上被释放。但是,在带实例后面的C++对象不知道这是发生的,因此它包含指向不再存在的C++数据集对象的指针。当带区尝试访问不存在的对象时,进程将崩溃。

GDAL团队知道这种设计不是Python程序员所期望的。不幸的是,这个设计很难纠正,所以它可能会保留一段时间。请咨询GDAL团队以获取更多信息。

该问题不限于GDAL带和数据集对象。它发生在对象之间有关系的其他区域。不幸的是没有完整的清单,所以你必须自己看。另一个已知的地方涉及OGR GetGeometryRef() 功能:

>>> feat = lyr.GetNextFeature()
>>> geom = feat.GetGeometryRef()     # geom contains a reference into the C++ geometry object maintained by the C++ feature object
>>> del feat                         # This deallocates the C++ feature object, and its C++ geometry
>>> print(geom.ExportToWkt())        # Crash here. The C++ geometry no longer exists
< Python crashes >

如果仔细阅读GDAL和OGR API文档,您将看到以“Ref”结尾的函数获得对内部对象的引用,而不是生成新的副本。这是问题可能发生的线索。使用“Ref”函数时要小心。还要注意以“直接”结尾的函数,例如 SetGeometryDirectly() ,它转移内部对象的所有权:

>>> point = ogr.Geometry(ogr.wkbPoint)
>>> feature = ogr.Feature(layer_defn)
>>> feature.SetGeometryDirectly(point)    # Transfers ownership of the C++ geometry from point to feature
>>> del feature                           # point becomes implicitly invalid, because feature owns the C++ geometry
>>> print(point.ExportToWkt())            # Crash here
< Python crashes >

“Ref”和“Directly”函数的优点是它们提供更快的性能,因为不需要创建重复的对象。缺点是你得小心这个问题。

当从这个层定义派生的特性仍然处于活动状态时,如果向OGR层添加一个新字段,Python就会崩溃

例如:

>>> feature = lyr.GetNextFeature()
>>> field_defn = ogr.FieldDefn("foo", ogr.OFTString)
>>> lyr.CreateField(field_defn)                       # now, existing features deriving from this layer are invalid
>>> feature.DumpReadable()                            # segfault
< Python crashes >

有关详细信息,请参见 #3552 .

带属性过滤器的图层 (SetAttributeFilter() )仅在使用时返回筛选的功能 GetNextFeature()

如果您阅读了 SetAttributeFilter() 小心你会看到警告的 OGR_L_GetNextFeature() . 这意味着如果你使用 GetFeature() ,而不是 GetNextFeature() ,则仍然可以访问和使用过滤器未覆盖的图层中的要素。 GetFeatureCount() 将尊重筛选器并显示筛选的正确功能数。但是,与 GetFeatureCount() 在循环中可能会导致一些微妙的混乱。在层对象上迭代或使用 GetNextFeature() 应为访问功能的默认方法:

>>> lyr = inDataSource.GetLayer()
>>> lyr.SetAttributeFilter("PIN = '0000200001'")      # this is a unique filter for only one record
>>> for i in range( 0, lyr.GetFeatureCount() ):
...    feat = lyr.GetFeature( i )
...    print(feat)                                    # this will print one feat, but it's the first feat in the Layer and not the filtered feat
...

某些对象包含 Destroy() 方法,但你不应该使用它

您可能会遇到称为 Destroy() method. This tutorial 甚至在第12页给出了具体的建议 Destroy. But according to email from Even RouaultDestroy() 不需要叫:

>I have some Python code that uses OGR geometry objects internally, creating
> them like this:
>
> point = ogr.Geometry(ogr.wkbPoint)
>
> Does this code need to explicitly destroy these geometries, like the
> following, to avoid leaks, or can it simply allow them to go out of scope
> and have Python's reference counting and garbage collector clean them up?
>
> point.Destroy()

There's no reason to call Destroy(), at all. Native object gets destroyed when
Python object goes out of scope, or when they are assigned to None. So replace
foo.Destroy() by foo = None if you really want to control when the underlying
C++ object is destroyed.

> I'm sorry for my ignorance here. I found a nice GDAL tutorial that seems to
> say they *should* be explicitly destroyed in certain circumstances (see
> http://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_slides2.pdf, page
> 12). But I have not really seen any other examples of this.
>

Destroy() was perhaps necessary with old-gen bindings, but I'm not even sure
of that... Perhaps this shouldn't have been exposed at all... But, as
mentioned in the slides, it is true that there are situations where you
shouldn't call Destroy() at all.

保存和关闭数据集/数据源

要保存和关闭GDAL栅格数据集或OGR矢量数据源,需要取消对该对象的引用,例如将其设置为 None ,其他值,或删除对象。如果数据集或数据源对象有多个副本,则需要取消对每个副本的引用。

例如,创建和保存栅格数据集:

>>> from osgeo import gdal
>>> driver = gdal.GetDriverByName('GTiff')
>>> dst_ds = driver.Create('new.tif', 10, 15)
>>> band = dst_ds.GetRasterBand(1)
>>> arr = band.ReadAsArray()  # raster values are all zero
>>> arr[2, 4:] = 50  # modify some data
>>> band.WriteArray(arr)  # raster file still unmodified
>>> band = None  # dereference band to avoid gotcha described previously
>>> dst_ds = None  # save, close

最后一次取消对栅格数据集的引用将写入数据修改并关闭栅格文件。 WriteArray(arr) 不会将阵列写入磁盘,除非GDAL块缓存已满(通常为40 MB)。

使用某些驱动程序,可以在不关闭的情况下间歇保存栅格数据集 FlushCache() . 类似地,向量数据集可以使用 SyncToDisk() . 但是,这两种方法都不能保证数据写入磁盘,因此首选的方法是如上所示解除分配。

在自定义错误处理程序中引发的异常不会被捕获

python绑定允许您指定一个可调用的python作为错误处理程序 (#4993 _). 但是,这些错误处理程序似乎是在单独的线程中调用的,并且引发的任何异常都不会传播回主线程 (#5186 ②)

所以如果你想 catch warnings as well as errors ,像这样的方法行不通:

from osgeo import gdal

def error_handler(err_level, err_no, err_msg):
    if err_level >= gdal.CE_Warning:
        raise RuntimeError(err_level, err_no, err_msg)  # this exception does not propagate back to main thread!

if __name__ == '__main__':
    # Test custom error handler
    gdal.PushErrorHandler(error_handler)
    gdal.Error(gdal.CE_Warning, 2, 'test warning message')
    gdal.PopErrorHandler()

但你可以这样做:

from osgeo import gdal

class GdalErrorHandler(object):
    def __init__(self):
        self.err_level = gdal.CE_None
        self.err_no = 0
        self.err_msg = ''

    def handler(self, err_level, err_no, err_msg):
        self.err_level = err_level
        self.err_no = err_no
        self.err_msg = err_msg

if __name__ == '__main__':
    err = GdalErrorHandler()
    gdal.PushErrorHandler(err.handler)
    gdal.UseExceptions()  # Exceptions will get raised on anything >= gdal.CE_Failure

    assert err.err_level == gdal.CE_None, 'the error level starts at 0'

    try:
        # Demonstrate handling of a warning message
        try:
            gdal.Error(gdal.CE_Warning, 8675309, 'Test warning message')
        except Exception:
            raise AssertionError('Operation raised an exception, this should not happen')
        else:
            assert err.err_level == gdal.CE_Warning, (
                'The handler error level should now be at warning')
            print('Handled error: level={}, no={}, msg={}'.format(
                err.err_level, err.err_no, err.err_msg))

        # Demonstrate handling of an error message
        try:
            gdal.Error(gdal.CE_Failure, 42, 'Test error message')
        except Exception as e:
            assert err.err_level == gdal.CE_Failure, (
                'The handler error level should now be at failure')
            assert err.err_msg == e.args[0], 'raised exception should contain the message'
            print('Handled warning: level={}, no={}, msg={}'.format(
                err.err_level, err.err_no, err.err_msg))
        else:
            raise AssertionError('Error message was not raised, this should not happen')

    finally:
        gdal.PopErrorHandler()

由其他软件的错误或行为引起的问题

升级或降级numpy时,GDAL函数中的Python会崩溃

大部分的GDAL的Python绑定都是用C++实现的。NUMPI的大部分核心是在C中实现的。GDAL的Python绑定的C++部分与NoMPY的C部分通过NUMPY的ABI(应用二进制接口)交互。这要求GDAL的Python绑定使用定义numpy C数据结构的numpy头文件进行编译。这些数据结构有时会在numpy版本之间发生变化。当这种情况发生时,新版本的numpy在二进制级别上与旧版本不兼容,必须重新编译GDAL Python绑定才能与新版本的numpy一起工作。当它们被重新编译时,它们可能无法与旧版本一起工作。

如果您获得了GDAL的Python绑定的预编译版本,例如 http://gisinternals.com/sdk.php 确保您查找了用于编译它们的numpy版本,并在您的计算机上安装了该版本的numpy。

无法从ArcGIS进程内地理处理工具(ArcGIS 9.3及更高版本)成功使用Python绑定

ArcGIS允许创建基于Python的自定义地理处理工具。直到ArcGIS 10,没有一种简单的方法可以将栅格数据读入内存。GDAL提供了这样一种机制。

从ArcGIS 9.3开始,地理处理工具可以在ArcGIS进程本身(ArcCatalog.exe或ArcMap.exe)中运行,也可以在单独的python.exe工作进程中运行。不幸的是,ArcGIS在如何运行进程内工具方面存在缺陷。因此,如果您使用进程内工具中的GDAL,它将在第一次运行时运行良好,但之后可能会失败 TypeError 在重新启动ArcGIS进程之前出现异常。例如,band.ReadAsArray()失败,原因是:

TypeError: in method 'BandRasterIONumpy', argument 1 of type 'GDALRasterBandShadow *

这是ArcGIS中的一个错误。请看 #3672 有关解决方案的完整详细信息和建议。