TDE工作（光变曲线拟合）

对所得到的TDE光变数据，想要得到其光变的特征，包括峰值星等，爆发时间，上升和下降持续时间等信息，因此对其光变的拟合是极为重要的。这里采用简单的自定义函数拟合，上升阶段用高斯函数，下降用power-law函数。值得注意，做拟合时，采用取对数方法，从而让拟合函数更为简单，目的在于提高拟合精度。

import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize as op
from scipy.optimize import minimize
import os
from scipy.optimize import curve_fit

plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文
plt.rcParams['axes.unicode_minus'] = False    # 用来正常显示负号

name='AT2018hco'
n=52
quiet=19.5637726

data = np.genfromtxt(name+ ".csv", delimiter=",")
mjd=data[:,0]
mag=-np.log(data[:,1])
mag_true=-(data[:,1]) #方便画图显示
magerr=data[:,2]

x1 = mjd[n:n+400] #定义拟合区间
y1 = mag[n:n+400]
x2=mjd[n-9:n+2]
y2=mag[n-9:n+2]

#phase=data[:,3]
###### PL
def f_1(x, A, B, C):
    return  A - C*np.log(x-B)


A, B, C = op.curve_fit(f_1, x1, y1, maxfev = 1000000)[0]
print("PL:"'  %-10.4f' %np.exp(-A))
print("PL:"'  %-10.4f' %B)
print("PL:"'  %-10.4f' %np.exp(C))

###### Gauss

def f_2(x, a, u, sig):
    return   a-((x - u) ** 2 / (2 * sig ** 2))
popt, pcov = curve_fit(f_2, x2, y2, maxfev=8000000)
a = popt[0]
u = popt[1]
sig = popt[2]
print("Gauss:"'  %-10.4f' %np.exp(-a))
print("Gauss:"'  %-10.4f' %u)
#print(sig)
#########一半处画线

delta1=quiet-np.exp(-A) 
delta2=quiet+mag_true[n]
delta3=quiet-np.exp(-a)
yline1=(-quiet+0.5*delta1)
yline2=(-quiet+0.5*delta2)
yline3=(-quiet+0.5*delta3)
good=yline1-yline2
#print(good)
######得到一半处time
x_1 = np.arange(np.min(x1), np.max(x1), 0.01)
x_2 = np.arange(np.min(x2), np.max(x2), 0.01)
yfit1 = A - C*np.log(x_1-B)
yfit2 = f_2(x_2, *popt)
yli=-np.log(-yline1)
yli2=-np.log(-yline2)
tdecay=np.exp((A-yli)/C)+B
t_decay=tdecay-B
trise=-((a-yli2)*2*sig**2)**0.5 +u
t_rise=u-trise
print("1/2decay:"'  %-10.4f' %t_decay)  
print("1/2rise:"'  %-10.4f' %t_rise)
########作图
plt.xlabel('mjd')
plt.ylabel('-mag')
plt.axhline(yline1,color="grey",linewidth=1,linestyle='--',label='1/2 下降（拟合）')
plt.axhline(yline2,color="grey",linewidth=1,linestyle='-.',label='1/2 （数据）')
plt.axhline(yline3,color="grey",linewidth=1,linestyle=':',label='1/2 上升（拟合）')
plt.title(name)
plt.errorbar(mjd[n:-1], mag_true[n:-1], yerr=0, fmt=".b", capsize=0)
mag_true=-data[:,1]
plt.errorbar(mjd[0:n], mag_true[0:n], yerr=0, fmt=".b", capsize=0)
plt.plot(x_1, -np.exp(-yfit1),color='red',label='PL曲线')
plt.plot(x_2, -np.exp(-yfit2),color='y',label='高斯曲线')
plt.legend() #指定legend的位置右下角
plt.show()

结果显示：

这个方法能够得到较好的近似参数，但是对拟合区间和拟合参数调整有着一定挑战性。