Где два двумерных массива начинают перекрывать друг друга? - PullRequest
6 голосов
/ 11 ноября 2010

В данный момент я работаю с выводом модели и не могу придумать хороший способ объединения двух массивов данных.Массивы A и B хранят разные данные, и записи в каждой соответствуют некоторой пространственной (x, y) точке - A содержит некоторый параметр, а B содержит выходные данные модели.Проблема в том, что B является пространственным подразделом A, то есть, если бы модель была для всего мира, A сохранял бы параметр в каждой точке Земли, а B сохранял бы выходные данные модели только для этих точек в Африке..

Так что мне нужно выяснить, на сколько B смещено относительно A - другими словами, мне нужно найти индексы, с которых они начинают перекрываться.Итак, если A.shape = (1000,1500), является ли B (750: 850, 200: 300) частью этого или подразделом (783: 835, 427: 440)?У меня есть массивы, связанные как с A, так и с B, которые хранят позиции (x, y) точек сетки для каждого.

Казалось бы, простая проблема - найти, где эти два массива перекрываются.И я могу решить эту проблему с KDTree от scipy.spatial достаточно просто, но это очень медленно.У кого-нибудь есть идеи получше?

Ответы [ 3 ]

1 голос
/ 16 ноября 2010

У меня есть массивы, связанные как с A, так и с B, в которых хранятся (x, y) позиции точек сетки для каждого.

В этом случае ответ должен быть довольно простым...

Находятся ли две сетки строго по одной и той же схеме сетки?Предполагая, что это так, вы можете просто сделать что-то вроде:

np.argwhere((Ax == Bx.min()) & (Ay == By.min())) 

Предполагая, что мировые координаты двух сеток увеличиваются в том же направлении, что и индикаторы сеток, это дает нижний левый угол заданной сетки,(И если они не увеличиваются в одном и том же направлении (т.е. отрицательно dx или dy), это просто дает один из других углов)

В приведенном ниже примере, очевидно, мы могли бы просто вычислитьправильные указания из ix = (Bxmin - Axmin) / dx и т. д., но при условии, что у вас более сложная система построения сетки, это все равно будет работать.Однако это предполагает, что две сетки находятся на одной и той же схеме сетки !Это немного сложнее, если они не ...

import numpy as np

# Generate grids of coordinates from a min, max, and spacing
dx, dy = 0.5, 0.5

# For the larger grid...
Axmin, Axmax = -180, 180
Aymin, Aymax = -90, 90

# For the smaller grid...
Bxmin, Bxmax = -5, 10
Bymin, Bymax = 30, 40

# Generate the indicies on a 2D grid
Ax = np.arange(Axmin, Axmax+dx, dx)
Ay = np.arange(Aymin, Aymax+dy, dy)
Ax, Ay = np.meshgrid(Ax, Ay)

Bx = np.arange(Bxmin, Bxmax+dx, dx)
By = np.arange(Bymin, Bymax+dy, dy)
Bx, By = np.meshgrid(Bx, By)

# Find the corner of where the two grids overlap...
ix, iy = np.argwhere((Ax == Bxmin) & (Ay == Bymin))[0]

# Assert that the coordinates are identical.
assert np.all(Ax[ix:ix+Bx.shape[0], iy:iy+Bx.shape[1]] == Bx) 
assert np.all(Ay[ix:ix+Bx.shape[0], iy:iy+Bx.shape[1]] == By) 
0 голосов
/ 16 ноября 2010

Мне нужно найти индексы, при которых они начинают перекрываться

Так вы ищете индексы из А или В?И является ли B строго прямоугольным?

Найти ограничивающий прямоугольник или выпуклую оболочку B действительно дешево.

0 голосов
/ 11 ноября 2010

Можете ли вы сказать больше? Какую модель вы используете? Что ты моделируешь? Как это вычисляется?

Можете ли вы привести размеры в соответствие, чтобы избежать подгонки? (то есть, если B не зависит от всех A, подключите только ту часть A, которую моделирует B, или вычислите скучные значения для частей B, которые не будут перекрывать A, и отбросьте эти значения позже)

...