3

我有 3 个 NumPy 数组,其中包含 UTM-X(256) 和 UTM-Y(256) 坐标,以及 UTM 中天气雷达 256x256 (km) 的累积降雨量(65536)。

我在网格边界内还有一个多边形,它是 UTM 中的集水边界。

我需要确定流域多边形(RADAR 数据的裁剪子集)的平均降雨量、最大值和最大值的位置。我已经确定了整个雷达网格的平均值。

所以问题是:如何对由 Polygon 确定的 NumPy 数组的子集执行分析?我原以为这将是一个非常常见的操作,但还没有找到任何 Python 脚本来执行此操作。

下面是数据集的示意图:

UTM 中的雷达网格和白色虚线集水多边形的图像注意红色星星是雨量计位置

4

2 回答 2

1

这是一种可能的方法的概述。

首先找到界定流域边界的多边形。假设你知道你的全套点的 UTM 坐标中的哪一个形成了流域边界,假设是这样的,

catchment = (UTM_X, UTM_Y) 点元组的 np.array

您可以使用scipy.spatial.ConvexHull找到该点集的边界

边界= scipy.spatial.ConvexHull(集水区)

接下来,对于您的降雨数据数组,您必须测试坐标是落在凸包边界之内还是之外。

This previous SO question有一些很好的答案,解释了进行坐标测试的方法。

最后,您将收集那些通过边界内测试的降雨数据点,并使用适当的 NumPy/SciPy 统计函数执行您想要执行的任何统计计算。

于 2015-07-15T01:14:17.897 回答
1

假设边界是作为多边形顶点的列表给出的,您可以让 matplotlib 在数据坐标上为您生成一个掩码,然后使用该掩码仅对轮廓内的值求和。

换句话说,当您有一系列坐标来定义标记感兴趣区域的多边形边界时,让 matplotlib 生成一个布尔掩码,指示该多边形内的所有坐标。然后可以使用此掩码仅提取轮廓内有限的降雨数据集。

以下简单示例向您展示了这是如何完成的:

import numpy as np
from matplotlib.patches import PathPatch
from matplotlib.path import Path
import matplotlib.pyplot as plt

# generate some fake data
xmin, xmax, ymin, ymax = -10, 30, -4, 20 
y,x = np.mgrid[ymin:ymax+1,xmin:xmax+1]
z = (x-(xmin+xmax)/2)**2 + (y-(ymin + ymax)/2)**2
extent = [xmin-.5, xmax+.5, ymin-.5, ymax+.5]
xr, yr = [np.random.random_integers(lo, hi, 3) for lo, hi
         in ((xmin, xmax), (ymin, ymax))] # create a contour

coordlist = np.vstack((xr, yr)).T  # create an Nx2 array of coordinates
coord_map = np.vstack((x.flatten(), y.flatten())).T # create an Mx2 array listing all the coordinates in field

polypath = Path(coordlist)
mask = polypath.contains_points(coord_map).reshape(x.shape) # have mpl figure out which coords are within the contour

f, ax = plt.subplots(1,1)
ax.imshow(z, extent=extent, interpolation='none', origin='lower', cmap='hot')
ax.imshow(mask, interpolation='none', extent=extent, origin='lower', alpha=.5, cmap='gray')
patch = PathPatch(polypath, facecolor='g', alpha=.5)
ax.add_patch(patch)
plt.show(block=False)
print(z[mask].sum())  # prints out the total accumulated

用 matplotlib 掩盖轮廓 在此示例中,xy表示您的UTM-XUTM-Y数据范围。z表示天气降雨量数据,但在这种情况下是一个矩阵,与平均降雨量的单列视图不同(它很容易重新映射到网格上)。

在最后一行,我总结了z轮廓内的所有值。如果你想要平均值,只需替换summean.

于 2015-07-15T02:22:01.227 回答