Закрепите определенную область в кубе IRIS с помощью шейп-файла

Question

Я обрабатываю кубики ириса, содержащие метеорологические данные (lon, lat, осадки, температуру, ...), и мне интересно рассчитать статистику в определенных областях (например, в стране).

В этом post объясняется, как обрезать куб с помощью поля (мин., Мин. Лат, макс. Лон, макс. Лат), но я хотел бы сделать еще один шаг и выбрать точную область с помощью шейп-файла.

Это post объясняет, что можно обрезать изображение с помощью файла формы, связанного с маской, но я не знаю, как я могу заставить его работать для моих кубов ириса.

Если кто-нибудь может дать мне пример или объяснить мне, как это сделать, это было бы очень полезно.

PS: Я довольно Noobie с питоном

Answer 1

Прочитав шейп, используя, например, Fiona что-то, как это должно работать:

from shapely.geometry import MultiPoint 

# Create a mask for the data 
mask = np.ones(cube.shape, dtype=bool) 

# Create a set of x,y points from the cube 
x, y = np.meshgrid(cube.coord(axis='X').points, cube.coord(axis='Y').points) 
lat_lon_points = np.vstack([x.flat, y.flat]) 
points = MultiPoint(lat_lon_points.T) 

# Find all points within the region of interest (a Shapely geometry) 
indices = [i for i, p in enumerate(points) if region.contains(p)] 

mask[np.unravel_index(indices)] = False 

# Then apply the mask 
if isinstance(cube.data, np.ma.MaskedArray): 
    cube.data.mask &= mask 
else: 
    cube.data = np.ma.masked_array(cube.data, mask) 

это работает только для 2D-кубов, а просто нуждается в тонкой настройки для более высоких размерностей, так что маска только по Lat/Lon размеров.

Я фактически реализовал это поведение в CIS недавно, так что вы можете сделать cube.subset(shape=region), что может быть проще для вас.