Хабрахабр

Контроллер Arduino с датчиком температуры и Python интерфейсом для динамической идентификации объектов управления

# -*- coding: utf8 -*-
import matplotlib.pyplot as plt
import time
start = time.time()
from scipy.interpolate import splev, splrep
import scipy.integrate as spint
import numpy as np
from scipy.integrate import quad
xx =np.array(np.arange(0,230,1))
yy1 =np.array([ 24.87, 25.06, 25.31, 25.5, 25.81, 26.06, 26.37, 26.62, 26.87, 27.12, 27.44, 27.69, 27.87, 28.12, 28.31, 28.56, 28.75, 28.94, 29.12, 29.31, 29.5, 29.69, 29.81, 30.0, 30.12, 30.25, 30.37, 30.56, 30.69, 30.81, 30.87, 31.0, 31.12, 31.25, 31.31, 31.44, 31.5, 31.62, 31.69, 31.81, 31.87, 31.94, 32.06, 32.13, 32.19, 32.25, 32.31, 32.38, 32.44, 32.5, 32.56, 32.63, 32.69, 32.75, 32.75, 32.81, 32.88, 32.94, 32.94, 33.0, 33.06, 33.13, 33.13, 33.19, 33.25, 33.25, 33.31, 33.31, 33.38, 33.38, 33.44, 33.44, 33.5, 33.5, 33.56, 33.56, 33.56, 33.63, 33.63, 33.69, 33.69, 33.69, 33.75, 33.75, 33.81, 33.81, 33.81, 33.88, 33.88, 33.88, 33.88, 33.94, 33.94, 33.94, 34.0, 34.0, 34.0, 34.0, 34.06, 34.06, 34.06, 34.06, 34.06, 34.13, 34.13, 34.13, 34.13, 34.13, 34.19, 34.19, 34.19, 34.19, 34.19, 34.25, 34.25, 34.25, 34.25, 34.25, 34.25, 34.31, 34.31, 34.31, 34.31, 34.31, 34.31, 34.31, 34.38, 34.38, 34.38, 34.38, 34.38, 34.38, 34.38, 34.38, 34.38, 34.44, 34.44, 34.44, 34.44, 34.44, 34.44, 34.44, 34.44, 34.44, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.5, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.56, 34.63, 34.56, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.63, 34.69, 34.63, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.69, 34.75, 34.69, 34.75, 34.69, 34.69, 34.75, 34.75, 34.75, 34.75, 34.75])
yy2=yy1-24.87#компенсация смещения нуля
yy=yy2/max(yy2)#нормирование """ Интерполяция переходной характеристики при помощи сплайнов"""
def h(x): spl = splrep(xx , yy ) return splev(x, spl) """ Численное интегрирование без смены координаты времени в соответствии с (6)"""
S1=(spint.quad(lambda x:1-h(x),xx[0],xx[len(xx)-1])[0])
S2=(spint.quad(lambda x:(1-h(x))*(S1-x),xx[0],xx[len(xx)-1])[0])
S3=(spint.quad(lambda x:(1-h(x))*(S2-S1*x+(1/2)*x**2),xx[0],xx[len(xx)-1])[0])
S4=(spint.quad(lambda x:(1-h(x))*(S3-S2*x+S1*(1/2)*x**2-(1/6)*x**3),xx[0],xx[len(xx)-1])[0]) """ Определение коэффициентов передаточной функции"""
b1=-S4/S3
a1=b1+S1
a2=b1*S1+S2
a3=b1*S2+S3 """ Возврат во временную область"""
def ff(x,t): j=(-1)**0.5 return (2/np.pi)*( ((b1*x*j+1)*np.e**(-25*x)/(a3*(x*j)**3+a2*(x*j)**2+a1*x*j+1)).real)*(np.sin(x*t)/x)
y=np.array([round(quad(lambda x: ff(x,t),0, 0.6)[0],2) for t in xx]) """ Определение критерия адекватности модели """
k=round(1-sum([(yy[i]-y[i])**2 for i in np.arange(0,len(yy)-1,1)])/sum([(yy[i])**2 for i in np.arange(0,len(yy)-1,1)]),5)
stop = time.time()
print ("Время работы программы :",round(stop-start,3))
plt.title('Идентификация модифицированным методом площадей.\n Адекватность модели: %s'%k)
plt.plot(xx, yy,label='W=(%s*p+1)/(%s*p**3+%s*p**2+%s*p+1)'%(round(b1,1),round(a3,1),round(a2,1),round(a1,1)))
plt.legend(loc='best')
plt.grid(True)
plt.show()

Теги
Показать больше

Похожие статьи

Кнопка «Наверх»
Закрыть