cómo enmascarar los datos de la matriz específica en función del shapefile

Aquí está mi pregunta:

  • los datos de la matriz numpy 2-d representan alguna propiedad de cada espacio de cuadrícula
  • el shapefile como la división administrativa del área de estudio (como una ciudad).

Por ejemplo:

http://i4.tietuku.com/84ea2afa5841517a.png

Toda el área tiene una red de 40×40 grids, y quiero extraer los datos dentro del área púrpura. En otras palabras, quiero enmascarar los datos fuera del límite administrativo en np.nan.

Mi primer bash

Etiqueto el número de la cuadrícula y selecciono los datos de la matriz específica en np.nan.

http://i4.tietuku.com/523df4783bea00e2.png

value[0,:] = np.nan value[1,:] = np.nan . . . . 

¿Puede alguien mostrarme un método más fácil para lograr el objective?

Añadir

Encontré una respuesta aquí que puede trazar los datos ráster en shapefile, pero los datos en sí no cambian.

Actualización -2016-01-16

Ya he resuelto este problema inspirado en algunas respuestas.
Alguien que esté interesado en este objective, revise estos dos mensajes que le he preguntado:
1. Punto de prueba con entrada / salida de un vector shapefile
2. Cómo utilizar la ruta de acceso recortada para el polígono de mapa base

El paso clave fue probar el punto dentro / fuera del shapefile que ya he transformado en shapely.polygon.

Paso 1. Rasterizar shapefile

Cree una función que pueda determinar si un punto en las coordenadas (x, y) está o no en el área. Consulte aquí para obtener más detalles sobre cómo rasterizar su shapefile en una matriz de las mismas dimensiones que su máscara de destino

 def point_is_in_mask(mask, point): # this is just pseudocode return mask.contains(point) 

Paso 2. Crea tu máscara.

 mask = np.zeros((height, width)) value = np.zeros((height, width)) for y in range(height): for x in range(width): if not point_is_in_mask(mask, (x, y)): value[y][x] = np.nan 

Lo mejor es usar matplotlib :

 def outline_to_mask(line, x, y): """Create mask from outline contour Parameters ---------- line: array-like (N, 2) x, y: 1-D grid coordinates (input for meshgrid) Returns ------- mask : 2-D boolean array (True inside) """ import matplotlib.path as mplp mpath = mplp.Path(line) X, Y = np.meshgrid(x, y) points = np.array((X.flatten(), Y.flatten())).T mask = mpath.contains_points(points).reshape(X.shape) return mask 

alternativamente, puede usar el método de shapely como se sugiere en la respuesta anterior. Puede acelerar los cálculos subdividiendo el espacio de forma recursiva, como se indica en este resumen (pero la solución matplotlib fue 1.5 veces más rápida en mis pruebas):

https://gist.github.com/perrette/a78f99b76aed54b6babf3597e0b331f8