Editar: Actualización cosas que reflejan sus aclaraciones anteriormente. Tu pregunta es mucho más clara ahora, ¡gracias!
Básicamente, solo quiere interpolar una matriz 2D en un punto arbitrario.
scipy.ndimage.map_coordinates es lo que quiere ....
Como entiendo su pregunta, tiene una matriz 2D de los valores "z" que va desde algunos xmin a xmax, ymin y ymax que en cada dirección.
Cualquier cosa fuera de esas coordenadas de límite desea devolver valores desde los bordes de la matriz.
map_coordinates tiene varias opciones para manejar puntos fuera de los límites de la grilla, pero ninguno de ellos hace exactamente lo que usted desea. En su lugar, podemos establecer cualquier cosa fuera de los límites para ubicarnos en el borde, y usar map_coordinates como de costumbre.
lo tanto, para utilizar map_coordinates, que necesita para convertir su coodinates reales:
| <1 2 3 4 5+
-------|----------------------------
<10000 | 3.6 6.5 9.1 11.5 13.8
20000 | 3.9 7.3 10.0 13.1 15.9
20000+ | 4.5 9.2 12.2 14.8 18.2
en el índice Coordenadas:
| 0 1 2 3 4
-------|----------------------------
0 | 3.6 6.5 9.1 11.5 13.8
1 | 3.9 7.3 10.0 13.1 15.9
2 | 4.5 9.2 12.2 14.8 18.2
Nota que sus límites se comportan de manera diferente en cada dirección ... En el x-direction, las cosas se comportan sin problemas, pero en la dirección y, estás pidiendo un break "duro", donde y = 20000 -> 3.9, pero y = 20000.000001 -> 4.5.
A modo de ejemplo:
import numpy as np
from scipy.ndimage import map_coordinates
#-- Setup ---------------------------
z = np.array([ [3.6, 6.5, 9.1, 11.5, 13.8],
[3.9, 7.3, 10.0, 13.1, 15.9],
[4.5, 9.2, 12.2, 14.8, 18.2] ])
ny, nx = z.shape
xmin, xmax = 1, 5
ymin, ymax = 10000, 20000
# Points we want to interpolate at
x1, y1 = 1.3, 25000
x2, y2 = 0.2, 50000
x3, y3 = 2.5, 15000
# To make our lives easier down the road, let's
# turn these into arrays of x & y coords
xi = np.array([x1, x2, x3], dtype=np.float)
yi = np.array([y1, y2, y3], dtype=np.float)
# Now, we'll set points outside the boundaries to lie along an edge
xi[xi > xmax] = xmax
xi[xi < xmin] = xmin
# To deal with the "hard" break, we'll have to treat y differently,
# so we're ust setting the min here...
yi[yi < ymin] = ymin
# We need to convert these to (float) indicies
# (xi should range from 0 to (nx - 1), etc)
xi = (nx - 1) * (xi - xmin)/(xmax - xmin)
# Deal with the "hard" break in the y-direction
yi = (ny - 2) * (yi - ymin)/(ymax - ymin)
yi[yi > 1] = 2.0
# Now we actually interpolate
# map_coordinates does cubic interpolation by default,
# use "order=1" to preform bilinear interpolation instead...
z1, z2, z3 = map_coordinates(z, [yi, xi])
# Display the results
for X, Y, Z in zip((x1, x2, x3), (y1, y2, y3), (z1, z2, z3)):
print X, ',', Y, '-->', Z
Esto produce:
1.3 , 25000 --> 5.1807375
0.2 , 50000 --> 4.5
2.5 , 15000 --> 8.12252371652
Esperemos que ayuda ...
Gracias por las aclaraciones! He actualizado mi respuesta a continuación. Creo que hace exactamente lo que quieres, ahora. –