Выполнение многих средств в NumPy - PullRequest
3 голосов
/ 03 июня 2011

Доброе утро, я использую фильтр Крессмана для выполнения взвешенных по расстоянию средних значений в Numpy. Я использую имплементацию Ball Tree (благодаря Джейку ВандерПласу), чтобы вернуть список местоположений для каждой точки в массиве запросов .. запросмассив (q) имеет форму [n, 3], и в каждой точке есть x, y, z, в точке. Я хочу сделать средневзвешенное значение точек, хранящихся в дереве .. код, обернутый вокруг дерева, возвращает точки в пределах определеннойрасстояние, поэтому я получаю массивы массивов переменной длины .. Я использую где найти непустые записи (то есть позиции, где было хотя бы несколько точек в радиусе влияния), создавая массив isgood ...

Затем я перебираю все точки запроса, чтобы получить средневзвешенное значение значений self.z (обратите внимание, что это может быть либо dims = 1, либо dims = 2, чтобы разрешить множественную совместную работу с сеткой)

, так чтокомпилирование с использованием карты или других более быстрых методов - это неоднородность длин массивов в пределах self.distances и self.locations ...до довольно зеленого до numpy / python, но я не могу придумать способ сделать этот массив мудрым (то есть не возвращаясь к циклам)

self.locations, self.distances = self.tree.query_radius( q, r, return_distance=True)
t2=time()
if debug: print "Removing voids"
isgood=np.where( np.array([len(x) for x in self.locations])!=0)[0]
interpol = np.zeros( (len(self.locations),) + np.shape(self.z[0]) )
interpol.fill(np.nan)
for dist, ix, posn, roi in zip(self.distances[isgood], self.locations[isgood], isgood, r[isgood]):
    interpol[isgood[jinterpol]] = np.average(self.z[ix], weights=(roi**2-dist**2) / (roi**2 + dist**2), axis=0)
    jinterpol += 1

так что ... Есть какие-нибудь намеки на то, как ускорить цикл?..

Для типичного отображения применительно к отображению данных метеорологического радара из диапазона, азимута, сетки высот и декартовой сетки, где у меня есть точки 240x240x34 и 4 переменные, для запроса дерева требуется 99 секунд (написанный Джейком вC и Cython ... это трудный шаг, так как вам нужно искать данные!) И 100 секунд, чтобы сделать расчет ... который, по моему мнению, медленный ??где мои накладные расходы?является ли np.mean эффективным или, как его называют миллионы раз, здесь можно добиться ускорения?я бы выиграл, используя float32, а не default64 ... или даже масштабируя до целых чисел (что было бы очень трудно избежать перестановки в весе ... любые советы с благодарностью получены!

1 Ответ

0 голосов
/ 03 июня 2011

Обсуждение относительных достоинств схемы Крессмана и использования весовой функции Гаусса можно найти по адресу:

http://www.flame.org/~cdoswell/publications/radar_oa_00.pdf

Ключ должен соответствовать параметру сглаживания данным (я рекомендую использовать значение, близкое к среднему интервалу между точками данных). Как только вы знаете параметр сглаживания, вы можете установить «радиус влияния», равный радиусу, в котором весовая функция падает до 0,01 (или что-то еще).

Насколько важна скорость? Если вы хотите вместо вызова экспоненциальной функции для определения веса, вы можете составить дискретную таблицу весов для некоторого фиксированного числа приращений радиуса, что значительно ускоряет вычисления. В идеале у вас должны быть данные вне границ сетки, которые можно использовать при отображении значений, окружающих точки сетки (даже в граничных точках сетки). Обратите внимание, что это НЕ истинная схема интерполяции - она ​​не будет возвращать наблюдаемые значения в точках данных точно. Как и схема Крессмана, это фильтр нижних частот.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...