Создание сетки и интерполяция (x, y, z) для контурной графики sagemath

!У меня есть значения в виде (x,y,z). Создав график list_plot3d, я ясно вижу, что они не совсем равномерно расположены. Обычно они образуют маленькие «капли» от 3 до 5 точек на плоскости xy. Итак, чтобы интерполяция и окончательный «контурный» график были лучше, или, я должен сказать, более гладким (?), мне нужно создать прямоугольную сетку (например, квадраты на шахматной доске), чтобы капли данных каким-то образом " сглаженный"? Я понимаю, что это может быть тривиально для некоторых людей, но я пробую это впервые и немного борюсь. Я смотрел на пакеты scipy, такие как scipy.interplate.interp2d, но графики, полученные в конце, действительно плохие. Может быть, краткое руководство по 2D-интерполяции в sagemath для таких любителей, как я? Какой-то совет? Спасибо.

ИЗМЕНИТЬ:

https://docs.google.com/file/d/0Bxv8ab9PeMQVUFhBYWlldU9ib0E/edit?pli=1

В основном это такие графики, которые он создает вместе с этим сообщением:

Warning:     No more knots can be added because the number of B-spline
coefficients
    already exceeds the number of data points m. Probably causes:
either
    s or m too small. (fp>s)
    kx,ky=3,3 nx,ny=17,20 m=200 fp=4696.972223 s=0.000000

Чтобы получить этот график, я просто запускаю эту команду:

f_interpolation = scipy.interpolate.interp2d(*zip(*matrix(C)),kind='cubic')
               plot_interpolation = contour_plot(lambda x,y:
                   f_interpolation(x,y)[0], (22.419,22.439),(37.06,37.08) ,cmap='jet', contours=numpy.arange(0,1400,100), colorbar=True)

               plot_all = plot_interpolation

               plot_all.show(axes_labels=["m", "m"])

Где матрица (с) может быть огромной матрицей, например, 10000 х 3 или даже намного больше, например 1000000 х 3. Проблема плохих графиков сохраняется даже при меньшем количестве данных, таких как изображение, которое я прикрепил сейчас, где матрица (С) была только 200 х 3 , Вот почему я начинаю думать, что, кроме возможного сбоя в программе, мой подход к использованию этой команды может быть совершенно неправильным, поэтому я прошу совета по поводу использования сетки, а не просто " бросая" мои данные в команду.


person CosmoSurreal    schedule 12.07.2012    source источник
comment
Можете ли вы опубликовать сюжет и описать, что с ним не так?   -  person user545424    schedule 13.07.2012


Ответы (1)


У меня была аналогичная проблема с использованием функции scipy.interpolate.interp2d. Насколько я понимаю, проблема возникает из-за того, что функции interp1d/interp2d и связанные с ними функции используют более старую оболочку FITPACK для базовых вычислений. Мне удалось решить проблему, похожую на вашу, для работы с использованием функций сплайна, которые основаны на более новой оболочке FITPACK. Сплайн-функции можно идентифицировать, потому что все они имеют заглавные буквы в своих именах здесь http://docs.scipy.org/doc/scipy/reference/interpolate.html. В рамках установки scipy эти новые функции находятся в scipy/interpolate/fitpack2.py, а функции, использующие более старые оболочки, находятся в fitpack.py.

Думаю, вам нужно RectBivariateSpline для ваших целей. Вот пример кода для реализации RectBivariateSpline:

import numpy as np
from scipy import interpolate

# Generate unevenly spaced x/y data for axes
npoints = 25
maxaxis = 100
x = (np.random.rand(npoints)*maxaxis) - maxaxis/2.
y = (np.random.rand(npoints)*maxaxis) - maxaxis/2.
xsort = np.sort(x)
ysort = np.sort(y)

# Generate the z-data, which first requires converting
# x/y data into grids
xg, yg = np.meshgrid(xsort,ysort)
z = xg**2 - yg**2

# Generate the interpolated, evenly spaced data
# Note that the min/max of x/y isn't necessarily 0 and 100 since
# randomly chosen points were used. If we want to avoid extrapolation,
# the explicit min/max must be found
interppoints = 100
xinterp = np.linspace(xsort[0],xsort[-1],interppoints)
yinterp = np.linspace(ysort[0],ysort[-1],interppoints)

# Generate the kernel that will be used for interpolation
# Note that the default version uses three coefficients for
# interpolation (i.e. parabolic, a*x**2 + b*x +c). Higher order
# interpolation can be used by setting kx and ky to larger 
# integers, i.e. interpolate.RectBivariateSpline(xsort,ysort,z,kx=5,ky=5)
kernel = interpolate.RectBivariateSpline(xsort,ysort,z)

# Now calculate the linear, interpolated data
zinterp = kernel(xinterp, yinterp)
person Michelle Lynn Gill    schedule 16.07.2012
comment
Похоже, что этот вопрос также можно решить с помощью связанной функции scipy.interpolate.griddata: stackoverflow.com/questions/11348708/ - person Michelle Lynn Gill; 16.07.2012
comment
Спасибо за ответ. Я до сих пор этого не понимаю. Может быть, лучше задать более простой вопрос, чтобы прояснить ситуацию. Если у меня есть значения, скажем, (долгота, широта, высота), и мне просто нужен контурный график с одинаковыми контурами высоты, что мне делать в sagemath? Буду ли я использовать interp2d, RectBivariateSpline или что-то еще? - person CosmoSurreal; 23.07.2012
comment
Я никогда не использовал sagemath, но если предположить, что все вызовы функций одинаковы (я понятия не имею), последние две строки кода — это то, что вам нужно для интерполяции данных. В этом и был смысл этого сценария. Для контурных графиков требуется, чтобы данные долготы и широты представляли собой сетку того же размера, что и высота, поэтому вам придется использовать команду numpy meshgrid, а затем отображать данные с помощью контурного графика matplotlib. Построение контурных диаграмм было рассмотрено на этом веб-сайте в ряде вопросов. - person Michelle Lynn Gill; 23.07.2012
comment
В поваренной книге scipy есть пример, который использует данные сетки и отображает здесь два разных типа контурных диаграмм: org/Cookbook/Matplotlib/Gridding_irregularly_spaced_data, который использует упомянутую выше команду griddata. Интересно, что команда meshgrid не требуется. Клянусь, мне всегда приходилось его использовать, но, возможно, я ошибаюсь или это требование изменилось с тех пор, как я последний раз проверял. Опять же, я не касался sagemath, поэтому вам придется выяснить, как эти команды отображаются в пространстве имен. - person Michelle Lynn Gill; 23.07.2012
comment
После запуска приведенного выше примера поваренной книги griddata я понял, почему нет команды meshgrid — это потому, что данные z также являются одномерным массивом. Таким образом, я лишь частично ошибался в использовании meshgrid. Данные x/y/z, которые входят в контур, должны быть одного размера. Если данные z линейны, то x и y должны быть линейными. Но если z представляет собой матрицу M x N, а x является линейным по длине N, а y является линейным по длине M, то необходимо использовать meshgrid. - person Michelle Lynn Gill; 23.07.2012