Wydaje się to proste, ale nie mogę tego dokładnie rozgryźć. Mam krzywą obliczoną na podstawie danych x, y. Potem mam linię. Chcę znaleźć wartości x, y dla miejsca, w którym te dwie przecinają się.Znaleźć przecięcie krzywej z polyfit
Oto, co mam do tej pory. Jest to bardzo mylące i nie daje prawidłowego wyniku. Mogę spojrzeć na wykres i znaleźć wartość przecięcia x oraz obliczyć prawidłową wartość y. Chciałbym usunąć ten ludzki krok.
import numpy as np
import matplotlib.pyplot as plt
from pylab import *
from scipy import linalg
import sys
import scipy.interpolate as interpolate
import scipy.optimize as optimize
w = np.array([0.0, 11.11111111111111, 22.22222222222222, 33.333333333333336, 44.44444444444444, 55.55555555555556, 66.66666666666667, 77.77777777777777, 88.88888888888889, 100.0])
v = np.array([0.0, 8.333333333333332, 16.666666666666664, 25.0, 36.11111111111111, 47.22222222222222, 58.333333333333336, 72.22222222222221, 86.11111111111111, 100.0])
z = np.polyfit(w, v, 2)
print (z)
p=np.poly1d(z)
g = np.polyval(z,w)
print (g)
N=100
a=arange(N)
b=(w,v)
b=np.array(b)
c=(w,g)
c=np.array(c)
print(c)
d=-a+99
e=(a,d)
print (e)
p1=interpolate.PiecewisePolynomial(w,v[:,np.newaxis])
p2=interpolate.PiecewisePolynomial(w,d[:,np.newaxis])
def pdiff(x):
return p1(x)-p2(x)
xs=np.r_[w,w]
xs.sort()
x_min=xs.min()
x_max=xs.max()
x_mid=xs[:-1]+np.diff(xs)/2
roots=set()
for val in x_mid:
root,infodict,ier,mesg = optimize.fsolve(pdiff,val,full_output=True)
# ier==1 indicates a root has been found
if ier==1 and x_min<root<x_max:
roots.add(root[0])
roots=list(roots)
print(np.column_stack((roots,p1(roots),p2(roots))))
plt.plot(w,v, 'r', a, -a+99, 'b-')
plt.show()
q=input("what is the intersection value? ")
print (p(q))
Jakieś pomysły, aby to zadziałało?
Dzięki
Fair ostrzegawczego , 'np.polynomial.polynomial.polyfit' zwraca współczynniki' [A, B, C] 'na' A + Bx + Cx^2 + ... ', które jest odwrotną kolejnością od tego, co' np.polyfit' (co pierwotnie użyłeś, @ user2843767) zwraca: '... + Ax^2 + Bx + C'. Nie jestem pewien, kto podjął taką decyzję, po prostu nie przyjmuj pierwszego wyniku i używaj go w 'np.poly1d' lub np.polyval, chyba że użyjesz' np.polyfit'. – askewchan
Rzeczywiście, uczciwe ostrzeżenie. Nie ma ostrzeżenia o wycofaniu i może nigdy nie być, ale dokumenty [są jasne] (http://docs.scipy.org/doc/numpy/reference/routines.polynomials.html), że sposób, aby przejść do nowego kodu jest pakietem wielomianowym, a nie starszym poly1d. – Jaime
Tak, i na szczęście nowy (er) pakiet ma również bardziej standardowe zamówienie. Dzięki za wskazanie tego linku, na pewno doradzę tylko pakiet wielomianowy. – askewchan