numpy: непрерывный объем 3d смежных 2d-фрагментов?

#python #numpy

#python #numpy

Вопрос:

Я хочу загрузить несколько смежных массивов изображений 2D C в объем 3D и получить к ним доступ с помощью интуитивного индексирования, например slice3 = volume[:,:,2] , что требует некоторого изменения формы исходного объединенного 1D-представления.

Теперь, поскольку я загружаю много изображений и вычисляю из него еще больше новых объемов, которые должны работать на вашем обычном ПК, меня беспокоит использование памяти, а также производительность вычислений.

В: Как я могу добиться того, чтобы объем 3D, а также фрагменты 2D оставались непрерывными, чтобы я мог эффективно выполнять некоторые операции с объемом, а некоторые — с отдельными фрагментами.

Вот несколько примеров, с которыми можно поиграть:

 import numpy as np

# dimensions:
rows   = 2
cols   = 2
slices = 3

# create array
a = np.arange(rows*cols*slices)
print(a)
# [ 0  1  2  3  4  5  6  7  8  9 10 11]
# this is the original concatenated input of the slices [[0,1],[2,3]], [[4,5],[6,7]], [[8,9],[10,11]]

# a contiguous?
print(a.flags['C_CONTIGUOUS'])
# True

a = a.reshape(rows,cols,slices)
print(a)
# a still contiguous?
print(a.flags['C_CONTIGUOUS'])
# True

# what about a slice?
print(a[:,:,0])
# [[0 3]
#  [6 9]]
# ouch! that's not the slice I wanted! I wanted to keep [[0,1],[2,3]] together
# this slice is of course also not contiguous:
print(a[:,:,0].flags['C_CONTIGUOUS'])
# False

# ok, let's start over
a = a.ravel()
print(a)
# [ 0  1  2  3  4  5  6  7  8  9 10 11]
a = a.reshape(slices,rows,cols)
a = a.swapaxes(0,1)
a = a.swapaxes(1,2)
# what about a slice?
print(a[:,:,0])
# [[0 1]
#  [2 3]]
# now that's the kind of slice I wanted!

# a still contiguous?
print(a.flags['C_CONTIGUOUS'])
# False

# slice contiguous?
print(a[:,:,0].flags['C_CONTIGUOUS'])
# True
# only halfway there.. :(
  

Опять же: есть ли способ добиться желаемой индексации фрагментов при сохранении объема, а также отдельных фрагментов C смежными?

Ответ №1:

Вариант 1

Не зацикливайтесь на том, какое измерение соответствует slices :

 a = a.reshape(slices, rows, cols)
  

Оставьте все как есть. Приятным побочным эффектом является то, что для получения, скажем, фрагмента 1 вы просто делаете

 a[1]
  

Нет необходимости в a[1, :, :] or even a[1, ...] .

Вариант 2

Если вы зацикливаетесь на том, какое измерение slices соответствует, вам необходимо скопировать данные для поддержания непрерывности. Однако вы можете получить его только с одной копией, а не использовать swapaxes дважды.

Один из способов — использовать transpose :

 a = a.reshape(slices, rows, cols).transpose(1, 2, 0)
  

Или вы можете использовать np.moveaxis :

 a = np.moveaxis(a.reshape(slices, rows, cols), 0, -1)
  

Комментарии:

1. Хороший ответ. Способ, которым NumPy хранит свои элементы, находится в «порядке C», где измерение в «правом» изменяется быстрее, чем измерение в «левом». Вот почему более естественно, чтобы измерение, представляющее срез, было первым. Это противоположно, например, MATLAB, где использование последнего измерения в качестве среза более естественно из-за того, что MATLAB использует «порядок Fortran».

2. Спасибо, что показали, как этого можно добиться более элегантно, хотя я пытаюсь избежать копирования. Я думаю, я заставлю своих коллег индексировать [срез строки col]. Поскольку строка (измерение y) и столбец (измерение x) уже обратны декартовым координатам, наличие среза (z-измерение) спереди, по крайней мере, согласовано (с приятным побочным эффектом адресации всего среза только с одним индексом, как упоминалось выше).

Ответ №2:

Нет, невозможно, чтобы 3D-массив a был C-непрерывным, а также a[:,:,0] чтобы он был C-непрерывным.

Почему? По определению, непрерывный массив C имеет единичный шаг в последнем измерении, например

 >>> a = np.arange(rows*cols*slices)
>>> a = a.reshape(slices,rows,cols)
>>> print([stride // a.itemsize for stride in a.strides])
[4, 2, 1]
  

Для любого C-смежного массива нарезка последнего измерения создаст не C-непрерывный массив, поскольку он не может иметь единичного шага в последнем измерении:

 >>> b = a[:, :, 0]
>>> print([stride // b.itemsize for stride in b.strides])
[4, 2]
  

Если вы хотите, чтобы ваши фрагменты были C-непрерывными, вам нужно будет нарезать в первом измерении и только в первом измерении:

 >>> c = a[0, :, :]
>>> print([stride // c.itemsize for stride in c.strides])
[2, 1]
  

Комментарии:

1. Спасибо, что внесли изменения в мое понимание.