python - vdal的gdal ReadAsarray非常慢

我试图简单地减去两个光栅并将结果保存到另一个光栅中一个输入图像是tif文件，另一个是vrt文件。（输出为TIF）
这些文件非常大，所以我打开它们，把它们分成几块，然后逐一运行，然后减去问题是速度太慢了！

import gdal
import numpy as np

rA = gdal.Open(raFileName)
rB = gdal.Open(rbFileName)

nodataA = rA.GetRasterBand(1).GetNoDataValue()
nodataB = rB.GetRasterBand(1).GetNoDataValue()


raster = gdal.GetDriverByName('GTiff').Create(outputFileName, ncols, nrows, 1 ,gdal.GDT_Float32,['COMPRESS=LZW','BigTIFF=YES'])

# tile size
trows = 5000
tcols = 5000
# number of tiles in output file   (ceil)
ntrows = (nrows-1)/trows+1
ntcols = (ncols-1)/tcols+1

# tiling because of memory problems
for r in range(ntrows):
    for c in range(ntcols):
        # number of rows/cols for tile r (in case of edge)
        rtrows = min(trows,nrows-r*trows)
        ctcols = min(tcols,ncols-c*tcols)

        # the data from the input files
        rA_arr = rA.GetRasterBand(1).ReadAsArray(c*tcols,r*trows,ctcols,rtrows)
        rB_arr = rB.GetRasterBand(1).ReadAsArray(c*tcols,r*trows,ctcols,rtrows)

        mask = np.logical_or(rA_arr==nodataA,rB_arr==nodataB)

        subtraction = rA_arr-rB_arr
        subtraction[mask] = nodata

        # writing this tile to the output raster
        rasterBand.WriteArray(subtraction,c*tcols,r*trows)
        raster.FlushCache()

我正在尝试减去的raster有16*16个tiles（5000*5000 pxls），大约两个小时后，它只经过了3行！
有没有办法提高性能？?

最佳答案

你应该做的第一件事是调查文件的布局。根据定义，VRT由128x128像素的块组成对于geotiff来说，它可以是任何东西。
您可以使用gdalinfo实用程序进行检查，或者：
rA.GetRasterBand(1).GetBlockSize()
还要检查VRT底层的文件。
读取块时，最好使用本机块大小或其倍数所以试着找出你使用的所有块大小的“交集”。因此，如果VRT是128x128，而Geotiff是512x1，那么以512x128（或一个倍数）的块读取将是最有效的。
如果这是不可能的，它可以帮助将gdal的缓存设置为尽可能高的：
gdal.SetCacheMax(2**30)
该值以字节为单位，因此2**30将是GiB这可以防止对磁盘进行不必要的I/O（速度较慢）。
它是什么类型的vrt，一个“简单”的马赛克/堆栈，还是包含各种计算？
您可以使用geotiff的两倍作为输入运行一次，以测试its是否是导致延迟的vrt（或者相反）。
如果你有很多nodata，你可以优化一下你的计算。但这似乎很简单，我不认为这是任何接近成为瓶颈的情况。
编辑：
基准
我做了一个快速的测试，我故意用一个低效的块大小来阅读。使用块大小为86400x1的86400x43200光栅。我用你的代码读了一个光栅（没有写）。maxcache设置为1MIB以避免大量缓存，这将降低效率。
使用86400*1的块需要：
1 loop, best of 3: 9.49 s per loop
使用5000*5000块需要：
1 loop, best of 3: 32.6 s per loop