使用GDAL将数据写入GeoTiff而不创建数据数组？

Question

有没有办法用gdal的WriteArray一个一个地写数据，而不是先创建一个完整的数组再写入？

我在创建一个大小为(50539,98357)的numpy数组时遇到了MemoryError，这意味着我的电脑内存不够用。我想我可以通过使用PyTables来解决这个问题，但我不想让事情变得更复杂。

>>> import numpy
>>> cols = 50539
>>> rows = 98357
>>> a1 = np.zeros((cols,rows))
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
NameError: name 'np' is not defined
>>> import numpy as np
>>> a1 = np.zeros((cols,rows))
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
MemoryError

更新： 最后我用了Abudis的解决方案，结合了一个稀疏矩阵来保存数据点，然后把每一行提取出来作为一个“密集”或标准矩阵。欢迎大家评论。

dataset = driver.Create(filename, cols, rows, number_of_bands, band_type)
band = dataset.GetRasterBand(1)
offset = 1 # i.e. the number of rows to write with each iteration

# values, pixel_x, and pixel_y defined earlier in my script
data = scipy.sparse.csr_matrix((values,(pixel_y,pixel_x)),dtype=float)

# iterate data writing for each row in data sparse array
for i in range(data.shape[0]):
    data_row = data[i,:].todense() # output row of sparse array as standard array
    band.WriteArray(data_row,0,offset*i)
band.SetNoDataValue(NULL_VALUE)
band.FlushCache()

Answer 1

我不是gdal方面的专家，但看起来这个方法有效。所以，诀窍在于使用WriteArray方法中的yoff参数。这样就可以分块写数据到文件里。我们基本上就是在写下一个数据块时设置偏移量。

import numpy as np
import gdal

cols = 50539
rows = 10000
offset = 1000

dst_filename = 'test.tif'
format = 'GTiff'
driver = gdal.GetDriverByName(format)

dst_ds = driver.Create(dst_filename, cols, rows, 1, gdal.GDT_Byte)

for i in range(10):
    # generate random integers from 1 to 10
    a = np.random.random_integers(1, 10, size=(offset, cols))
    # write data to band 1
    dst_ds.GetRasterBand(1).WriteArray(a, 0, offset * i)

dst_ds = None

这样可以生成一个50539乘10000像素的tiff文件。如果你真的需要98357行数据，我想你可以把offset的值设置为1，把rows设置为98357。从理论上讲，这应该是可行的，因为在处理较小的数组（10乘20像素）时是有效的。

编辑：别忘了把range(10)改成类似range(98357)或xrange(98357)的东西。