matplotlib.mlab.griddata очень медленный и возвращает массив nan при вводе правильных данных

Question

Я пытаюсь отобразить набор данных с нерегулярной сеткой (необработанные спутниковые данные) с соответствующими широтами и долготами в набор широт и долгот с регулярной сеткой, заданный basemap.makegrid().Я использую matplotlib.mlab.griddata с mpl_toolkits.natgrid установленным.Ниже приведен список переменных, которые используются в качестве выходных данных whos в ipython и некоторые статистические данные по переменным:

Variable   Type       Data/Info
-------------------------------
datalat    ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)
datalon    ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)
gridlat    ndarray    1200x1000: 1200000 elems, type `float64`, 9600000 bytes (9 Mb)
gridlon    ndarray    1200x1000: 1200000 elems, type `float64`, 9600000 bytes (9 Mb)
var        ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)

In [11]: var.min()
Out[11]: -30.0

In [12]: var.max()
Out[12]: 30.0

In [13]: datalat.min()
Out[13]: 27.339874

In [14]: datalat.max()
Out[14]: 47.05302

In [15]: datalon.min()
Out[15]: -137.55658

In [16]: datalon.max()
Out[16]: -108.41629

In [17]: gridlat.min()
Out[17]: 30.394031556984299

In [18]: gridlat.max()
Out[18]: 44.237140350357713

In [19]: gridlon.min()
Out[19]: -136.17646180595321

In [20]: gridlon.max()
Out[20]: -113.82353819404671

datalat и datalon - координаты оригинальных данных

gridlat и gridlon - это координаты для интерполяции в

var содержит фактические данные

Использование этих переменных, когда я вызываю griddata(datalon, datalat, var, gridlon, gridlat), заняло целых 20минут до завершения и возвращает массив nan.При просмотре данных широта и долгота кажутся правильными: исходные координаты перекрывают часть новой области и несколько точек данных, лежащих за пределами новой области.У кого-нибудь есть предложения?Значения наночастиц предполагают, что я делаю что-то глупое ...

Answer 1

Похоже, что подпрограмма mlab.griddata может вводить дополнительные ограничения в ваши выходные данные, которые могут быть необязательными. Хотя входные местоположения могут быть любыми, выходные местоположения должны быть регулярной сеткой - так как ваш пример находится в широтном / долгом пространстве, ваш выбор проекции карты может нарушить это (то есть обычная сетка в х / у не регулярная сетка в широте / долготе).

Вы можете попробовать подпрограмму interpolate.griddata из SciPy в качестве альтернативы - вам нужно будет объединить ваши переменные местоположения в один массив, так как сигнатура вызова отличается: что-то вроде

import scipy.interpolate
data_locations = np.vstack(datalon.ravel(), datalat.ravel()).T
grid_locations = np.vstack(gridlon.ravel(), gridlat.ravel()).T
grid_data      = scipy.interpolate.griddata(data_locations, val.ravel(),
                                            grid_locations, method='nearest')

для интерполяции ближайшего соседа. Это получает местоположения в массиве с 2 столбцами, соответствующими вашим 2 измерениям. Вы также можете выполнить интерполяцию в преобразованном пространстве вашей картографической проекции.

Answer 2

Скорее всего, griddata слишком сложен. Он предназначен для работы со случайно выбранными данными. Ваши данные почти наверняка регулярно отбираются - но не в той же самой сетке, что и целевая выходная сетка.

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

Существуют некоторые готовые решения, которые могут помочь. GDAL является хорошим примером.

Кроме того, этот тип проблемы часто обсуждается в ГИС. См:

https://gis.stackexchange.com/questions/10430/changing-image-projection-using-python

Answer 3

Если вы используете pclormesh, вам не нужно выполнять какую-либо интерполяцию.pcolormesh с радостью примет структуру данных, как вы указали здесь:

from mpl_toolkits.basemap import Basemap
m = Basemap(-----)
x,y = m(datalon, datalat)
m.pcolormesh(x,y,var)
plt.show()

, пожалуйста, используйте это и скажите мне, работает ли это или нет.

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

Использование pcolormesh для построения данных орбиты

Answer 4

Если ваши данные находятся в сетке так, что точка данных в точке (datalon[i], datalat[j]) находится в data[i,j], тогда вы можете использовать scipy.interpolate.RectBivariateSpline вместо griddata.Однако некоторые географические библиотеки могут предлагать больше функциональных возможностей.