20081111

Python на службе науки (черновик)

Не знаю, соберусь ли это дело систематизировать, пока буду просто добавлять в этот пост короткие «рецепты» по применению Python в околонаучных целях (создание отчётов, построение графиков и тому подобное). Поживём-увидим.

Преобразование числа из экспоненциальной формы (например, –1.2e+003) в формулу LaTeX (например, -1.2 10³). Регулярное выражение:
>>> import re
>>> re.sub('([\d.+-]+)[Ee]\+?([-]?)0*([\d]+)', '\\1 \\cdot 10^{\\2\\3} ','-1.2e+003')
'-1.2 \\cdot 10^{3} '


Прозрачные/полупрозрачные легенды (подписи) в графиках pylab:
l=legend()
l.get_frame().set_alpha(0.2)


Формулы LaTeX в подписях и надписях в pylab: просто заключить нужный фрагмент надписи в $…$, например '$\rho$'.

Открытый вопрос — как делать русские подписи в pylab? Дополнение: вопрос решился сам. См. заметку Русские буквы в matplotlib/pylab.

Произвольная надпись, выровненная по центру — (xlim()[0]+xlim()[1])*0.5, ha="center" — нижнего края — lim()[0], va="bottom") — графика в pylab:
text((xlim()[0]+xlim()[1])*0.5,ylim()[0],"your message",ha="center",va="bottom",fontsize=10)


Дана пара массивов (списков) t и x с временем измерения и измеренной величиной, отобрать только те пары время-значение, которые попадают в заданный временной интервал [t0, t1] (или какому другому условию):
xt=[ pair for pair in zip(x,t) if t0 <= pair[0] <= t1 ]


Из полученных пар выделить опять отдельные массивы t и x:
t=[pair[0] for pair in tx]
x=[pair[1] for pair in tx]
или с использованием array из scipy:
t,x=array(tx)[:,0],array(tx)[:,1]


Кусочно-линейная интерполяция зависимости t→x на другом наборе точек t2 (тоже список, как и t):
from scipy import interp
interp(t2,t,x)


Транспонирование двумерного массива — проще всего array([[1,2,],[3,4]]).transpose().

Линейная регрессия (аналогично — фиттирование более сложных функций, модифицировать количество параметров и вычисление ошибки в residual()). Как и прежде, предполагаю, что точки для фиттирования заданы в массивах t и x:
from scipy import array
from scipy.optimize import leastsq
def residual(params,t,x):
    a,b=params
    x_ab=array([a*x=t_i+b for t_i in t])
    return x_ab-x
[a,b],ier=leastsq(residual,(1,1),args=(t,x))


Интегрирование системы обыкновенных дифференциальных уравнений (ОДУ):
from scipy import linspace, array
from scipy.integrate import odeint

def dotu(u,t0):
        "расчёт производных в момент времени t0"
        x,y,z=u # вектор переменных
        sigma,rho,beta=10,28,8.0/3 # параметры модели
        dotx=sigma*(y-x)
        doty=x*(rho-z)-y
        dotz=x*y-beta*z
        return [dotx,doty,dotz]

t=linspace(0,50,1e4) # точки интегрирования
u0=[-5,-5,1] # начальное условие
u=odeint(dotu,u0,t) # решение


Чтобы построить график полученного решения (2D), транспонируем u и строим как обычно:
from pylab import plot, show
u=array(u).transpose()
plot(u[1],u[2]) # строим фазовый портрет y-z
show()
Ярлыки: , ,