cubic_spĺines.py 2.05 KB
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Wed Jul 14 16:01:06 2021

@author: liaf-ankylosaurus-admin
"""

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

xs = np.arange(-0.5, 100, 0.1)

y = np.sin(xs)

cs = CubicSpline(xs, y)

xs = np.arange(-0.5, 100, 0.1)

fig, ax = plt.subplots(figsize=(6.5, 4))

ax.plot(xs, y, 'o', label='data')

ax.plot(xs, np.sin(xs), label='true')

ax.plot(xs, cs(xs), label="S")

#ax.plot(xs, cs(xs, 1), label="S'")

"""
ax.plot(xs, cs(xs, 2), label="S''")

ax.plot(xs, cs(xs, 3), label="S'''")

ax.set_xlim(-0.5, 9.5)
"""
ax.legend(loc='lower left', ncol=2)

plt.show()

#%%

from scipy import interpolate, optimize

freqs = [ 0, 1, 2, 3, 4, 5]
pd = [12,14,22,39,58,77]

def f(x, x_points, y_points):
    tck = interpolate.splrep(x_points, y_points)
    return interpolate.splev(x, tck)

print(f(1.25, freqs, pd))

#%%
from scipy import interpolate, optimize

amps = np.arange(0.05, 0.3, 0.01)
pd = [0, 0.01, 0.03, 0.032, 0.035, 0.039, 0.45, 0.48, 0.53, 0.55, 0.56, 0.59, 0.65, 0.69, 0.73, 0.81, 0.9, 0.93, 0.96, 0.98, 1.15, 1.23, 1.40, 1.41, 1.45]

pdvalue = 0.1

interp_fn = interpolate.interp1d(amps,pd,kind='cubic')
interp_fn2 = lambda x: interp_fn(x)-pdvalue

try:
    a = optimize.newton(interp_fn2, 0.5*np.mean(amps))
except:
    try:    
        a = optimize.newton(interp_fn2, 1*np.mean(amps))
    except:
        try:
            a = optimize.newton(interp_fn2, 1.5*np.mean(amps))
        except:
            print('bue')
            
print(a)

#%%
x = np.arange(0.05, 0.3, 0.01)
y = [0, 0.01, 0.03, 0.032, 0.035, 0.039, 0.45, 0.48, 0.53, 0.55, 0.56, 0.59, 0.65, 0.69, 0.73, 0.81, 0.9, 0.93, 0.96, 0.98, 1.15, 1.23, 1.40, 1.41, 1.45]
f = interpolate.interp1d(x, y, 'quadratic')

xs = np.arange(0.05, 0.29, 0.0001)
plt.plot(x,y,'x',xs,f(xs))
#xnew = np.arange(70,111,1)

#%%
yToFind = 0.5
yreduced = np.array(y) - yToFind
freduced = interpolate.UnivariateSpline(x, yreduced, s=0)
a = freduced.roots()


for i in range(len(a)):
    plt.axvline(a[i])
plt.axhline(yToFind)