Python：将点X，Y和Z网格化以提取统计属性

Question

对不起，这个问题很简单，但是我是Python新手，我需要同样的帮助。

我的数据采用点格式：X，Y，Z。其中X和Y为坐标，z为值。

我的问题是：用0.5 m x 0.5 m（或1 x 1 m）创建一个栅格（TIF或ASCII），其中每个像素的值是Z的平均值。如果我在像素i中没有指向点该值必须为NAN。

我很高兴为我可以学习和实现的一些代码提供帮助，

在此先感谢您的帮助，我真的很需要。

我试图学习并编写代码：

from osgeo import gdal, osr, ogr
import math, numpy
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.mlab as ml
import matplotlib.delaunay

tiff = 'test.tif'

gridSize = 0.5
# my area boundary
xmax, xmin = 640000.06, 636999.83
ymax, ymin = 6070000.3, 6066999.86

# number of row and columns
nx = int(math.ceil(abs(xmax - xmin)/gridSize))
ny = int(math.ceil(abs(ymax - ymin)/gridSize))

# Plot the points
plt.scatter(x,y,c=z)
plt.axis([xmin, xmax, ymin, ymax])
plt.colorbar()
plt.show()

# Generate a regular grid.
xi = np.linspace(xmin, xmax, nx)
yi = np.linspace(ymin, ymax, ny)
xi, yi = np.meshgrid(xi, yi)

从这一点上，我很难理解如何索引x，y，z点以知道它们在哪里下降。我的第一个想法是给索引一个网格网格并标记点。之后，我可以对像素内的点取平均值。空像素（没有当前点的像素）为NAN。

但是我不知道这是处理数据的正确方法。

之后，我编写了以下代码以通过GDAL以TIFF格式保存

target_ds = gdal.GetDriverByName('GTiff').Create(tiff, nx,ny, 1, gdal.GDT_Byte) #gdal.GDT_Int32
target_ds.SetGeoTransform((xmin, gridSize, 0,ymax, 0, -gridSize,))

if EPSG == None:
    proj = osr.SpatialReference()
    proj.ImportFromEPSG(EPSG)
    # Make the target raster have the same projection as the source
    target_ds.SetProjection(proj.ExportToWkt())
else:
    # Source has no projection (needs GDAL >= 1.7.0 to work)
    target_ds.SetProjection('LOCAL_CS["arbitrary"]')

target_ds.GetRasterBand(1).WriteArray(numpy.zeros((ny,nx)))
target_ds = None

真的非常感谢您的帮助

Answer 1

要走的路：

• 定义网格spacing（浮点数），即相同尺寸中两个像素/体素中点之间的距离
• 弄清楚你所需要的电网，即在网格点数量的大小x和y尺寸，N_x以及N_y。
• numpy使用例如创建两个大小均为零的数组np.zeros([N_x, N_y])
• 遍历您的（x，y，v）点集，然后
  • 将每个（x，y）对投影到其相应的像素中，该像素通过两个（整数）索引标识： x_i, y_i = tuple([int(c//spacing) for c in (x, y)])
  • 将一个1添加到一个数组中(x_i, y_i)（保留“ count”）
  • 将值添加v到其他数组(x_i, y_i)（保留值的总和）
• 填充完两个数组后，将value-of-values-array除以count-array。0 / 0.0将自动分配给NaN，而c / 0.0将自动分配给Inf。