Нарязване на 2D масиви с помощта на индекси от масиви в python

Работя с парчета от 2D numpy масив. За да избера срезовете, имам индексите, съхранени в масиви. Например имам:

mat = np.zeros([xdim,ydim], float)
xmin = np.array([...]) # Array of minimum indices in x
xmax = np.array([...]) # Array of maximum indices in x
ymin = np.array([...]) # Array of minimum indices in y
ymax = np.array([...]) # Array of maximum indices in y
value = np.array([...]) # Values

Където ... просто обозначава някои цели числа, изчислени преди това. Всички масиви са добре дефинирани и имат дължини ~265000. Това, което искам да направя, е нещо като:

mat[xmin:xmax, ymin:ymax] += value

По такъв начин, че за първите елементи ще имам:

mat[xmin[0]:xmax[0], ymin[0]:ymax[0]] += value[0]
mat[xmin[1]:xmax[1], ymin[1]:ymax[1]] += value[1]

и така нататък, за ~265000 елемента на масива. За съжаление това, което току-що написах, не работи и извежда грешката: IndexError: invalid slice.

Опитвам се да използвам np.meshgrid, както е предложено тук: NumPy: използвайте 2D индексен масив от argmin в 3D срез, но все още не ми работи. Освен това търся pythonic начин да направя това, като избягвам циклите for.

Всяка помощ ще бъде много оценена!

Благодаря!


python numpy arrays slice
person fbecerra    schedule 03.03.2014    source източник
comment
np.array() е невалиден   -  person zhangxaochen    schedule 03.03.2014
comment
Съжалявам за объркването, добавих малко текст, за да стане по-ясно.   -  person fbecerra    schedule 03.03.2014
comment
Припокриват ли се резените ви? Колко голям е mat?   -  person Jaime    schedule 03.03.2014
comment
Да, може да се припокриват. Размерът на постелката се определя от потребителя, така че може да бъде или доста малък, или доста голям.   -  person fbecerra    schedule 03.03.2014
comment
Може да искате да опитате нещо като Cython или Numba, вместо да търсите трикове за нарязване.   -  person user2357112 supports Monica    schedule 03.03.2014
comment
Благодаря за съвета! Имам предвид да напиша Cython/Numba версия на кода, но за момента просто трябва да използвам numpy.   -  person fbecerra    schedule 03.03.2014

Отговори (1)


arrow_upward
2
arrow_downward

Не мисля, че има задоволителен начин да векторизирате проблема си, без да прибягвате до Cython или други подобни. Позволете ми да очертая как би могло да изглежда чисто numpy решение, което трябва да изясни защо това вероятно не е много добър подход.

Първо, нека да разгледаме 1D случай. Не можете да направите много с куп срезове в numpy, така че първата задача е да ги разширите в отделни индекси. Кажете, че вашите масиви са били:

mat = np.zeros((10,))
x_min = np.array([2, 5, 3, 1])
x_max = np.array([5, 9, 8, 7])
value = np.array([0.2, 0.6, 0.1, 0.9])

Тогава следният код разширява ограниченията на срезовете в списъци с (евентуално повтарящи се) индекси и стойности, обединява ги заедно с bincount и ги добавя към оригиналния mat:

x_len = x_max - x_min
x_cum_len = np.cumsum(x_len)
x_idx = np.arange(x_cum_len[-1])
x_idx[x_len[0]:] -= np.repeat(x_cum_len[:-1], x_len[1:])
x_idx += np.repeat(x_min, x_len)
x_val = np.repeat(value, x_len)
x_cumval = np.bincount(x_idx, weights=x_val)
mat[:len(x_cumval)] += x_cumval

>>> mat
array([ 0. ,  0.9,  1.1,  1.2,  1.2,  1.6,  1.6,  0.7,  0.6,  0. ])

Възможно е да разширите това до вашия 2D случай, въпреки че е всичко друго, но не и тривиално и нещата започват да стават трудни за следване:

mat = np.zeros((10, 10))
x_min = np.array([2, 5, 3, 1])
x_max = np.array([5, 9, 8, 7])
y_min = np.array([1, 7, 2, 6])
y_max = np.array([6, 8, 6, 9])
value = np.array([0.2, 0.6, 0.1, 0.9])

x_len = x_max - x_min
y_len = y_max - y_min
total_len = x_len * y_len
x_cum_len = np.cumsum(x_len)
x_idx = np.arange(x_cum_len[-1])
x_idx[x_len[0]:] -= np.repeat(x_cum_len[:-1], x_len[1:])
x_idx += np.repeat(x_min, x_len)
x_val = np.repeat(value, x_len)
y_min_ = np.repeat(y_min, x_len)
y_len_ = np.repeat(y_len, x_len)
y_cum_len = np.cumsum(y_len_)
y_idx = np.arange(y_cum_len[-1])
y_idx[y_len_[0]:] -= np.repeat(y_cum_len[:-1], y_len_[1:])
y_idx += np.repeat(y_min_, y_len_)
x_idx_ = np.repeat(x_idx, y_len_)
xy_val = np.repeat(x_val, y_len_)
xy_idx = np.ravel_multi_index((x_idx_, y_idx), dims=mat.shape)
xy_cumval = np.bincount(xy_idx, weights=xy_val)
mat.ravel()[:len(xy_cumval)] += xy_cumval

Което произвежда:

>>> mat
array([[ 0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ],
       [ 0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0.9,  0.9,  0.9,  0. ],
       [ 0. ,  0.2,  0.2,  0.2,  0.2,  0.2,  0.9,  0.9,  0.9,  0. ],
       [ 0. ,  0.2,  0.3,  0.3,  0.3,  0.3,  0.9,  0.9,  0.9,  0. ],
       [ 0. ,  0.2,  0.3,  0.3,  0.3,  0.3,  0.9,  0.9,  0.9,  0. ],
       [ 0. ,  0. ,  0.1,  0.1,  0.1,  0.1,  0.9,  1.5,  0.9,  0. ],
       [ 0. ,  0. ,  0.1,  0.1,  0.1,  0.1,  0.9,  1.5,  0.9,  0. ],
       [ 0. ,  0. ,  0.1,  0.1,  0.1,  0.1,  0. ,  0.6,  0. ,  0. ],
       [ 0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0.6,  0. ,  0. ],
       [ 0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. ]])

Но ако имате 265 000 двуизмерни среза с произволен размер, тогава индексиращите масиви ще влязат в многото милиони елементи много бързо. Необходимостта да се справите с четенето и записването на толкова много данни може да отмени подобренията на скоростта, които идват с използването на numpy. Честно казано, съмнявам се, че това изобщо е добър вариант, ако не друго, заради това колко загадъчен ще стане кодът ви.

person Jaime    schedule 03.03.2014