Точки данных, выходящие за пределы интерполируемой сетки, в то время как сетка, безусловно, покрывает эти точки - PullRequest
Новая
       104

Точки данных, выходящие за пределы интерполируемой сетки, в то время как сетка, безусловно, покрывает эти точки

2 голосов
/ 16 июня 2020

Я пытаюсь интерполировать разреженные данные по сетке, но наблюдаю довольно странное поведение. Белые точки - это именно то место, где у меня есть значения, и я полагаюсь на алгоритм линейной интерполяции, чтобы заполнить другие сетки, где это возможно. Я понимаю, что этот тип интерполяции не идеален из-за очевидного отсутствия данных, но почему некоторые точки, в которых у меня есть данные, выходят за пределы сетки, по которой я интерполирую? Это обычное явление? Это не изменится, даже если я сделаю сетку грубее.

Я был бы признателен за понимание того, почему это происходит (возможно, как работает линейная интерполяция), или есть ли способы исправить это. См. Красные кружки на рисунке ниже, например: enter image description here
Точки данных, предоставленные для интерполяции, выходящие за пределы сетки, которая интерполирована по

Следующие - это некоторый код интерполяции, который сгенерировал данные с координатной сеткой. 

#mesh grid
xg = np.linspace(-130, -60, num=70)
yg = np.linspace(20,50,num=30)
Xg,Yg = np.meshgrid(xg,yg)

zg1 = griddata(points1, df2['tempratio'], (Xg, Yg), method = 'linear')

from mpl_toolkits.basemap import Basemap

lon_0 = xg.mean()
lat_0 = yg.mean()

m = Basemap(width=5000000, height=3500000,
           resolution='l', projection='stere',\
           lat_ts=40, lat_0=lat_0, lon_0=lon_0)

xm, ym = m(Xg, Yg)

cs = m.pcolormesh(xm,ym,zg1,shading='flat',cmap=plt.cm.Reds)

1 Ответ

2 голосов
/ 17 июня 2020

griddata присваивает значения вершинам сетки, то есть 70x30 точек. pcolormesh окрашивает не вершины, а прямоугольники между ними. Всего из заданных вершин образованы прямоугольники 69x29. Таким образом, одна строка и один столбец zg1 будут удалены. Чтобы противостоять этому, к координатам можно добавить дополнительную строку и дополнительный столбец и сдвинуть все на половину прямоугольника в каждом направлении.

Это по-прежнему не заставляет griddata включать все заданные точки, но делает шаг к желаемый результат. Также может помочь более плотная сетка. (Выбор «ближайшей» вместо «линейной» интерполяции заполнит всю сетку.)

Вот код, иллюстрирующий происходящее: 

import numpy as np
from scipy.interpolate import griddata
from matplotlib import pyplot as plt

def extend_range(x):
    dx = (x[1] - x[0]) / 2
    return np.append( x - dx, x[-1] + dx)

N = 10
points1 = np.vstack([np.random.randint(-130, -60, N), np.random.randint(20, 50, N)]).T
tempratio = np.random.randint(0, 20, N)

xg = np.linspace(-130, -60, num=15)
yg = np.linspace(20, 50, num=10)
Xg, Yg = np.meshgrid(xg, yg)

zg1 = griddata(points1, tempratio, (Xg, Yg), method='linear')

fig, axs = plt.subplots(ncols=2, figsize=(12, 4))
for ax in axs:
    ax.scatter(Xg, Yg, c=zg1, cmap='coolwarm', ec='g', s=80, zorder=2, label='griddata')
    ax.scatter(points1[:,0], points1[:,1], c=tempratio, cmap='coolwarm', ec='black', s=150, zorder=3, label='given data')

    if ax == axs[0]:
        ax.pcolormesh(xg, yg, zg1, shading='flat', cmap='coolwarm')
        ax.set_title('given x and y ranges')
    else:
        #todo: convert xg and yg to map coordinates
        ax.pcolormesh(extend_range(xg), extend_range(yg), zg1, shading='flat', cmap='coolwarm')
        ax.set_title('extended x and y ranges')
    ax.legend()
plt.show()

example plot

...