Интерполяция Scipy на массиве numpy
У меня есть справочная таблица, которая определяется следующий путь:
| <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
TR_ua1 = 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] ])
- Элементы строки заголовка (гд) <1,2,3,4,5 +
- Столбец заголовка (inc) элементы <10000, 20000, 20001 +
Пользователь введет пример значения (1.3, 25,000), (0.2, 50,000), и так далее. scipy.interpolate() должен интерполировать для определения правильного значения.
В настоящее время единственный способ, которым я могу сделать это, с набором if/elifs как иллюстрируется ниже. Я вполне уверен существует лучший, более эффективный способ сделать это
Вот то, что я имею до сих пор:
import numpy as np
from scipy import interpolate
if (ua == 1):
if (inc <= low_inc): # low_inc = 10,000
if (hh <= 1):
return TR_ua1[0][0]
elif (hh >= 1 & hh < 2):
return interpolate( (1, 2), (TR_ua1[0][1], TR_ua1[0][2]) )
1 ответ
Изменить: обновлено, чтобы отразить ваши пояснения выше. Теперь ваш вопрос стал более ясным, спасибо!
По сути, вы просто хотите интерполировать 2D-массив в произвольной точке.
scipy.ndimage.map_coordinates - это то, что вам нужно ....
Насколько я понимаю ваш вопрос, у вас есть 2D-массив значений «z», который варьируется от некоторого xmin до xmax и от ymin до ymax в каждом направлении.
Все, что находится за пределами этих ограничивающих координат, вы хотите вернуть с краев массива.
map_coordinates имеет несколько параметров для обработки точек за пределами сетки, но ни один из них не делает именно то, что вы хотите. Вместо этого мы можем просто установить что-либо за пределами границ так, чтобы оно лежало на краю, и использовать map_coordinates как обычно.
Итак, чтобы использовать map_coordinates, вам нужно превратить ваши фактические координаты:
| <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
в координаты индекса:
| 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
Обратите внимание, что ваши границы ведут себя по-разному в каждом направлении ... В направлении x все работает плавно, но в направлении y вы просите «жесткого» разрыва, где y = 20000 -> 3,9, но y = 20000.000001 -> 4,5.
В качестве примера:
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
Это дает:
1.3 , 25000 --> 5.1807375
0.2 , 50000 --> 4.5
2.5 , 15000 --> 8.12252371652
Надеюсь, это поможет ...