np.where () для исключения данных, где координаты находятся слишком близко друг к другу - PullRequest
Новая
       40

np.where () для исключения данных, где координаты находятся слишком близко друг к другу

3 голосов
/ 26 сентября 2019

Я делаю апертурную фотометрию на скоплении звезд, и чтобы облегчить обнаружение фонового сигнала, я хочу смотреть только на звезды дальше, чем n пикселей (n = 16 в моем случае).У меня есть 2 массива, xs и ys, с x- и y-значениями всех координат звезд: используя np. Где я должен найти индексы всех звезд, где расстояние до всех других звезд>> =n

До сих пор мой метод представлял собой цикл for 

import numpy as np

# Lists of coordinates w. values between 0 and 2000 for 5000 stars

xs = np.random.rand(5000)*2000
ys = np.random.rand(5000)*2000

# for-loop, wherein the np.where statement in  question is situated
n = 16
for i in range(len(xs)):
    index = np.where( np.sqrt( pow(xs[i] - xs,2) + pow(ys[i] - ys,2)) >= n)

Из-за скопления звезд довольно близко друг к другу я ожидал серьезного сокращения данных, хотя даже когда пыталсяn = 1000 У меня осталось около 4000 точек данных

Ответы [ 3 ]

2 голосов
/ 26 сентября 2019

Использование только numpy (и часть ответа здесь ) 

X = np.random.rand(5000,2) * 2000
XX = np.einsum('ij, ij ->i', X, X)
D_squared = XX[:, None] + XX - 2 * X.dot(X.T)
out = np.where(D_squared.min(axis = 0) > n**2)

Использование scipy.spatial.pdist

from scipy.spatial import pdist, squareform
D_squared = squareform(pdist(x, metric = 'sqeuclidean'))
out = np.where(D_squared.min(axis = 0) > n**2)

Использование KDTree для максимально быстрого: 

from scipy.spatial import KDTree

X_tree = KDTree(X)
in_radius = np.array(list(X_tree.query_pairs(n))).flatten()
out = np.where(~np.in1d(np.arange(X.shape[0]), in_radius))
0 голосов
/ 26 сентября 2019

Вы можете использовать np.subtract.outer для создания парных сравнений.Затем вы проверяете для каждой строки, находится ли расстояние ниже 16 для ровно одного элемента (что является сравнением с самим конкретным началом): 

distances = np.sqrt(
    np.subtract.outer(xs, xs)**2
    + np.subtract.outer(ys, ys)**2
)
indices = np.nonzero(np.sum(distances < 16, axis=1) == 1)
0 голосов
/ 26 сентября 2019 
np.random.seed(seed=1)
xs = np.random.rand(5000,1)*2000
ys = np.random.rand(5000,1)*2000

n = 16

mask = (xs>=0)

for i in range(len(xs)):
    if mask[i]:
        index = np.where( np.sqrt( pow(xs[i] - x,2) + pow(ys[i] - y,2)) <= n)
        mask[index] = False
        mask[i] = True

x = xs[mask]
y = ys[mask]

print(len(x))

4220

...