Tarea 1, parte 2

Pregunta 1

Para convertir del sistema de filtros de la HRC al sistema UBVRI, se realiza el ajuste siguiendo la pauta del ejemplo (que transforma del sistema HRC a WFC) de Sirianni et al. (2005): $$(1): V = F555W + c_{1,1}(V-I)+c_{2,1}(V-I)^2$$ $$(2): I = F814W + c_{2,1}(V-I)+c_{2,2}(V-I)^2$$ Los valores de las constantes se encuentran tabuladas para cada filtro, y dependen del valor $TCOL = F555W - F814W$, el cual también se asumirá como el valor $(V-I)$ para una primera aproximación. Luego se aplica la ecuación (1) y (2) para calcular V e I. Se calcula $TCOL$ con los valores actualizados de $(V-I)$ y se vuelven a aplicar las ecuaciones (1) y (2). Repitiendo el procedimiento, los valores $(V-I)$ cambian cada vez menos. Se comprueba que realizar 6 repeticiones basta para obtener un resultado con mejor precisión que la magnitud de entrada, pero se realizarán un par más en esta ocasión.

Se escribe un catálogo y un CMD con esta nueva información.



In [1]:

    
from astropy.io import ascii
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline

tbl = ascii.read('n121_match.cat')

def transform(v,i):
	c1f555 = [-0.09,-0.124]
	c2f555 = [0.034,0.018]
	c1f814 = [0.06,0.001]
	c2f814 = [-0.099,0.013]
	for j in range(8):
		tcol = v-i
		v = np.where(tcol<0.6, (tbl['f555wMAG']+c1f555[0]*tcol+c2f555[0]*tcol*tcol), (tbl['f555wMAG']+c1f555[1]*tcol+c2f555[1]*tcol*tcol))
		i = np.where(tcol<0.1, (tbl['f814wMAG']+c1f814[0]*tcol+c2f814[0]*tcol*tcol), (tbl['f814wMAG']+c1f814[1]*tcol+c2f814[1]*tcol*tcol))
	return v,i

tbl['V'], tbl['I'] = transform(tbl['f555wMAG'],tbl['f814wMAG'])

# Genera catalogo
data = [tbl['V'],tbl['I'],tbl['f814wALPHA'],tbl['f814wDELTA']]
ascii.write(data,'n121_ubvri.cat',delimiter='\t')

# Genera plots
fig, ax = plt.subplots(2,1,figsize=(13,20))
for axis in ax:
    axis.set_xlim(-0.5,1.5)
    axis.set_ylim(24,16)
    axis.set_xlabel('$V-I$',fontsize=20)
    axis.set_ylabel('$V$',fontsize=20)

ax[0].plot(tbl['V']-tbl['I'],tbl['V'],'k.',alpha=0.4,ms=3)
hess = ax[1].hexbin(tbl['V']-tbl['I'],tbl['V'],bins='log',gridsize=220,mincnt=1)
leg = fig.colorbar(hess,ax=ax[1])
leg.set_label('Cuentas (log)')

plt.savefig('cmd_n121_ubvri', dpi=300)
plt.show()
plt.close()

Pregunta 2

Usando como parámetros de entrada:

Edad = 10.5 Gyr
Fe/H = -1.56
Distancia = 61KPc
$A_V = 0.1272$



In [2]:

    
from isochrones.dartmouth import Dartmouth_Isochrone
iso = Dartmouth_Isochrone(bands=['V','I'])
model = iso.isochrone(age=10.021189299,feh=-1.56,distance=61000,AV=0.1272,dm=1e-5)

fig, ax = plt.subplots(figsize=(10,10))
ax.set_xlim(-0.5,1.5)
ax.set_ylim(24,16)
ax.set_xlabel('$V-I$',fontsize=20)
ax.set_ylabel('$V$',fontsize=20)

ax.plot(tbl['V']-tbl['I'],tbl['V'],'k.',alpha=0.4,ms=3)
ax.plot(model.V_mag - model.I_mag,model.V_mag,'g',lw=2)

plt.show()
plt.close()

Se obtuvo la isócrona de la imagen (línea verde). Al menos pasa cerca de los puntos obtenidos por fotometría. Sin embargo, para una mejro aproximación sería necesaria una magnitud V mayor, de forma que todo el CMD se mueva a la derecha y hacia abajo.

Pregunta 3

Para calcular el ridge line del set de datos, primero se dividió el eje y en rangos de 0.22 magnitudes, que es lo mínimo posible para que no aparezcan secotres sin puntos. En cada rango se calculó la mediana de V-I y se ensambló un set de coordenadas que contiene (mediana,V). El segmento en el eje y que contiene la rama horizontal (entre 19.4 y 20) fue ignorado para este cálculo. Posteriormente, se ajustó un polinomio de grado 6.

Se grafica en azul la línea ajustada. Se ajustarán la magnitudes en la última pregunta con tal de poder comparar esta isócrona con el ridge line y un mejor fit.



In [3]:

    
s = 0.22
x = []
y1 = np.arange(16.5,19.4,step=s)
y2 = np.arange(19.5,23.5,step=s)
y = np.append(y1,y2)
for i in y:
	a = np.where((i<tbl['V']) & (tbl['V']<i+s) & (tbl['V']-tbl['I']>0.4))
	x.append(np.median(tbl['V'][a]-tbl['I'][a]))

    
fig, ax = plt.subplots(figsize=(10,10))
ax.set_xlim(-0.5,1.5)
ax.set_ylim(24,16)
ax.set_xlabel('$V-I$',fontsize=20)
ax.set_ylabel('$V$',fontsize=20)

ax.plot(tbl['V']-tbl['I'],tbl['V'],'k.',alpha=0.4,ms=3)

# ajuste de polinomio a la ridge line
p = np.poly1d(np.polyfit(y,x,7))
ax.plot(p(y),y,'b',lw=2)

#nueva isocrona que se acerque a nuestra ridge line
iso = Dartmouth_Isochrone(bands=['V','I'])
model = iso.isochrone(age=9.89,feh=-1.55,distance=61000,dm=1e-5)

ax.plot(model.V_mag - model.I_mag,model.V_mag,'r',lw=2)

plt.show()
plt.close()

pregunta 4

Ahora que se tiene la ridge line, se busca una isócrona que se acerque de la mejor manera. Suponiendo que la distancia al cúmulo es el dato mejor conocido, se cambian los valores de edad y metalicidad. Se encuentra que el mejor fit se logra con

Edad: 9.89Gyr
Fe/H: -1.55
distancia: 61KPc

El fit se traza en rojo. Si bien no es el más cercano a la línea azul, adopta una forma cercana. El trabajo de desplazar el gráfico del CMD puede hacerse al cambiar los valores V e I. Haciendo $V = V+0.15$ Se obtiene el siguiente fit, que corresponde a la primera isócrona, con parámetros de literatura. Hecho así, podría intentarse refinar más calibración en vista de que la forma de la curva es la del CMD.



In [4]:

    
fig, ax = plt.subplots(figsize=(10,10))
ax.set_xlim(-0.5,1.5)
ax.set_ylim(24,16)
ax.set_xlabel('$V-I$',fontsize=20)
ax.set_ylabel('$V$',fontsize=20)

ax.plot(tbl['V']-tbl['I']+0.15,tbl['V']+0.15,'k.',alpha=0.4,ms=3)

#nueva isocrona que se acerque a nuestra ridge line
iso = Dartmouth_Isochrone(bands=['V','I'])
model = iso.isochrone(age=10.021189299,feh=-1.56,distance=61000,AV=0.1272,dm=1e-5)

ax.plot(model.V_mag - model.I_mag,model.V_mag,'r',lw=2)

plt.show()
plt.close()

Otras isócronas que podrían haberse utilizado son las de Teramo y las de Padova. En la literatura, ambas entregan edades mayores de este cúmulo (11.8 Gyr y 11.2 Gyr respectivamente, según katharina Glatt et al. 2008)