Ajustement de la surface polynomiale en Python 3D, en fonction de l'ordre

Je suis actuellement en train de travailler avec les données astronomiques parmi lesquels j'ai des images de la comète. Je voudrais enlever le fond de ciel dégradé dans ces images, le temps de capture (crépuscule). Le premier programme que j'ai développé pour le faire a pris l'utilisateur a sélectionné des points de Matplotlib "ginput" (x,y) tiré les données pour chaque coordonnée (z), puis quadriller les données dans un nouveau tableau avec SciPy "griddata."

Depuis l'arrière-plan est supposé ne varient que très peu, je voudrais un ajustement 3d faible polynôme à cet ensemble de (x,y,z) des points. Cependant, le "griddata" ne permet pas une entrée de commande:

griddata(points,values, (dimension_x,dimension_y), method='nearest/linear/cubic')

Toutes les idées sur une autre fonction qui peut être utilisée ou une méthode pour l'élaboration d'un lea-places ajustement qui va me permettre de contrôler l'ordre?

source d'informationauteur Sebastian De Pascuale

  1. 31

    Griddata utilise une définition des splines. Une 3ème commande spline n'est pas la même chose qu'un 3ème ordre polynomial (au lieu de cela, c'est une autre 3ème ordre polynôme à chaque point).

    Si vous voulez juste pour s'adapter à une 2D, 3e ordre polynôme à vos données, puis faire quelque chose comme ce qui suit pour estimer le 16 coefficients à l'aide de tous de vos points de données.

    import itertools
import numpy as np
import matplotlib.pyplot as plt
def main():
# Generate Data...
numdata = 100
x = np.random.random(numdata)
y = np.random.random(numdata)
z = x**2 + y**2 + 3*x**3 + y + np.random.random(numdata)
# Fit a 3rd order, 2d polynomial
m = polyfit2d(x,y,z)
# Evaluate it on a grid...
nx, ny = 20, 20
xx, yy = np.meshgrid(np.linspace(x.min(), x.max(), nx), 
np.linspace(y.min(), y.max(), ny))
zz = polyval2d(xx, yy, m)
# Plot
plt.imshow(zz, extent=(x.min(), y.max(), x.max(), y.min()))
plt.scatter(x, y, c=z)
plt.show()
def polyfit2d(x, y, z, order=3):
ncols = (order + 1)**2
G = np.zeros((x.size, ncols))
ij = itertools.product(range(order+1), range(order+1))
for k, (i,j) in enumerate(ij):
G[:,k] = x**i * y**j
m, _, _, _ = np.linalg.lstsq(G, z)
return m
def polyval2d(x, y, m):
order = int(np.sqrt(len(m))) - 1
ij = itertools.product(range(order+1), range(order+1))
z = np.zeros_like(x)
for a, (i,j) in zip(m, ij):
z += a * x**i * y**j
return z
main()

    Ajustement de la surface polynomiale en Python 3D, en fonction de l'ordre

  2. 10

    La suite de la mise en œuvre de polyfit2d utilise le disponible numpy méthodes numpy.polynomial.polynomial.polyvander2d et numpy.polynomial.polynomial.polyval2d

    #!/usr/bin/env python3
import unittest
def polyfit2d(x, y, f, deg):
from numpy.polynomial import polynomial
import numpy as np
x = np.asarray(x)
y = np.asarray(y)
f = np.asarray(f)
deg = np.asarray(deg)
vander = polynomial.polyvander2d(x, y, deg)
vander = vander.reshape((-1,vander.shape[-1]))
f = f.reshape((vander.shape[0],))
c = np.linalg.lstsq(vander, f)[0]
return c.reshape(deg+1)
class MyTest(unittest.TestCase):
def setUp(self):
return self
def test_1(self):
self._test_fit(
[-1,2,3],
[ 4,5,6],
[[1,2,3],[4,5,6],[7,8,9]],
[2,2])
def test_2(self):
self._test_fit(
[-1,2],
[ 4,5],
[[1,2],[4,5]],
[1,1])
def test_3(self):
self._test_fit(
[-1,2,3],
[ 4,5],
[[1,2],[4,5],[7,8]],
[2,1])
def test_4(self):
self._test_fit(
[-1,2,3],
[ 4,5],
[[1,2],[4,5],[0,0]],
[2,1])
def test_5(self):
self._test_fit(
[-1,2,3],
[ 4,5],
[[1,2],[4,5],[0,0]],
[1,1])
def _test_fit(self, x, y, c, deg):
from numpy.polynomial import polynomial
import numpy as np
X = np.array(np.meshgrid(x,y))
f = polynomial.polyval2d(X[0], X[1], c)
c1 = polyfit2d(X[0], X[1], f, deg)
np.testing.assert_allclose(c1,
np.asarray(c)[:deg[0]+1,:deg[1]+1],
atol=1e-12)
unittest.main()