Несогласованное декартово произведение с использованием NumPy

Question

У меня есть следующие определения трех векторов осей:

arena = np.array([[-0.1, 0.1], [-0.1, 0.1], [0.06, 0.14]], dtype=np.float32)
resolution = np.array([0.01, 0.01, 0.01], dtype=np.float32)

axes_size = np.array((arena[:, 1] - arena[:, 0]) / resolution).astype(np.int32)
axes = [np.linspace(*arena[i], axes_size[i], endpoint=False, dtype=np.float32)
        for i in np.arange(axes_size.shape[0])]

Обратите внимание, что axes_size = array([20, 20, 8], dtype=int32).

Я пытаюсь векторизовать и обобщить следующий фрагмент кода, который создать декартово произведение axes:

points = np.zeros([*axes_size, 3])
points[:, :, :, 0] = axes[0][:, None, None]
points[:, :, :, 1] = axes[1][None, :, None]
points[:, :, :, 2] = axes[2][None, None, :]

Обратите внимание, что points.shape = (20, 20, 8, 3).

Я пробовал следующее:

points = np.stack(np.meshgrid(*axes), axis=-1)
points = np.array(np.meshgrid(*axes)).T.reshape(-1, 3).T.reshape(*axes_size, 3)

Обе строки, кажется, не работают, когда сравнивая результаты, используя np.testing.assert_allclose(), и я не могу понять почему.

Answer 1

Проблема в том, что np.meshgrid имеет настройку по умолчанию, которая не соответствует ожиданиям большинства людей.

Добавьте indexing='ij' для создания разумного поведения:

points = np.stack(np.meshgrid(*axes, indexing='ij'), axis=-1)

Это пройдет тест np.allclose.

По умолчанию np.meshgrid равно indexing='xy', что меняет порядок первых двух осей. Предположительно это удобно для данных изображений, где изображение 800x400 (ширина 800, высота 400) должно быть массивом shape (400, 800).

>>> np.meshgrid(np.arange(20), np.arange(21), np.arange(3))[0].shape
(21, 20, 3)

>>> np.meshgrid(np.arange(20), np.arange(21), np.arange(3), indexing='ij')[0].shape
(20, 21, 3)

Если я имею дело с многомерными массивами, я стараюсь иметь разные размеры по разным осям, например (20, 21, 8, 3), а не (20, 20, 8, 3), чтобы я мог уловить этот тип проблемы.