import numpy as np
import matplotlib.pyplot as plt

xs = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20]
fx = [5, 6, 8, 10, 12, 13, 12, 10, 8, 10, 8, 11, 7, 9, 11, 10, 9, 12, 11, 6]

degree = 7


def interpolate(x):
    index = 0
    assert (xs[0] <= x <= xs[-1])
    for i in range(len(xs)):
        if xs[i] > x:
            index = i
            break
    left = index - ((degree - 1) // 2)
    right = index + ((degree - 1) // 2)
    left = max(0, left)
    right = min(right, len(xs) - 1)
    result = 0
    for i in range(left, right + 1):
        cur = fx[i]
        for j in range(left, right + 1):
            if i == j:
                continue
            cur *= (x - xs[j])
            cur /= (xs[i] - xs[j])
        result += cur
    return result


x = xs[0]
dx = 0.01

gx = []
gy = []

while x < xs[-1]:
    gx.append(x)
    gy.append(interpolate(x))
    x += dx

gx = np.array(gx)
gy = np.array(gy)

fig, ax = plt.subplots()
ax.scatter(gx, gy)
ax.plot(np.array(xs), np.array(fx), color='red')
plt.show()

for i in range(len(gx)):
    for j in range(len(xs)):
        if abs(gx[i] - xs[j]) < 1e-5:
            print(xs[j], gy[i], fx[j])