Интерполяция трехмерного массива в Python.Как избежать петель? - PullRequest
3 голосов
/ 13 октября 2011

У меня есть массив, который я хочу интерполировать по 1-м осям.Сейчас я делаю это, как в следующем примере:

import numpy as np
from scipy.interpolate import interp1d

array = np.random.randint(0, 9, size=(100, 100, 100))
new_array = np.zeros((1000, 100, 100))
x = np.arange(0, 100, 1)
x_new = np.arange(0, 100, 0.1)

for i in x:
    for j in x:
        f = interp1d(x, array[:, i, j])
        new_array[:, i, j] = f(xnew)

Используемые мной данные представляют собой 10-летние усредненные значения за 5 дней для каждой широты и долготы в домене.Я хочу создать массив ежедневных значений.

Я также пытался использовать сплайны.Я действительно не знаю, как они работают, но это было не намного быстрее.

Есть ли способ сделать это без использования циклов for?Если необходимо использовать циклы for, есть ли другие способы ускорить это?

Заранее благодарим вас за любые предложения.

Ответы [ 2 ]

7 голосов
/ 16 октября 2011

Вы можете указать аргумент оси для interp1d: import numpy as np from scipy.interpolate import interp1d array = np.random.randint(0, 9, size=(100, 100, 100)) x = np.linspace(0, 100, 100) x_new = np.linspace(0, 100, 1000) new_array = interp1d(x, array, axis=0)(x_new) new_array.shape # -> (1000, 100, 100)

6 голосов
/ 13 октября 2011

Поскольку вы интерполируете данные с регулярной сеткой, обратите внимание на использование scipy.ndimage.map_coordinates.

В качестве быстрого примера:

import numpy as np
import scipy.ndimage as ndimage

interp_factor = 10
nx, ny, nz = 100, 100, 100
array = np.random.randint(0, 9, size=(nx, ny, nz))

# If you're not familiar with mgrid: 
# http://docs.scipy.org/doc/numpy/reference/generated/numpy.mgrid.html
new_indicies = np.mgrid[0:nx:interp_factor*nx*1j, 0:ny, 0:nz]

# order=1 indicates bilinear interpolation. Default is 3 (cubic interpolation)
# We're also indicating the output array's dtype should be the same as the 
# original array's. Otherwise, a new float array would be created.
interp_array = ndimage.map_coordinates(array, new_indicies, 
                                       order=1, output=array.dtype)
interp_array = interp_array.reshape((interp_factor * nx, ny, nz))
...