En observaciones astronómicas e imágenes en general, llamamos resolución espacial a la distancia angular minima a la que pueden estar dos fuentes puntuales de luz y aun poder ser reconocidas como objetos individuales.
En el caso de la astronomía, este efecto tiene que ver con la dispersión de la luz al atravezar la atmósfera, la cual hace que una estrella, que debería en principio aparecer como una fuente puntual (pues las estrellas están muy lejos), aparezca en cambio como una mancha. Así, si dos estrellas están demasiado cerca sus manchas se superpondrán hasta el punto en que sea imposible distinguirlas como fuentes individuales (Ver imágenes en este link)
Para modelar este efecto, típicamente consideramos la acción de la atmósferacomo la convolución de la imagen "perfecta" (como se vería desde el espacio) con un kernel gaussiano. El ancho de esa función gaussiana 2D caracteriza las condiciones de observación, varía con las condiciones climáticas y para cada sitio de la Tierra.
La resolución espacial normalmente se toma como el FWHM de la gaussiana caracteristica registrada durante una observación. Es decir, si dos estrellas están a una distancia aparente en el cielo menor que el FWHM del efecto atmosférico, la luz de ambas fuentes se mezclará después de la convolución hasta el punto de impedir reconocerlas de modo individual.
Además, la atmósfera puede interactuar de maneras distintas con la luz de distintas longitudes de onda, de manera que el ancho de la gaussiana puede ser distinto para observaciones con diferentes filtros.
El objetivo de esta tarea es medir de forma aproximada la resolución espacial en una noche de observación en Zapatoca, Santander (Colombia), a partir de una foto del cielo estrellado.
* Leer la imagen almacenada en la carpeta data como un array de numpy. Ese array estará compuesto de 3 matrices, una tras otra, correspondiente a los canales R,G,B * Combinar los 3 array para generar una versión blanco y negro de la imagen, en la cual ella consiste de una sola matriz 2D. Puede usar su intuición y prueba y error para combinar las 3 imágenes, explicando el procedimiento elegido. Esto será más interesante que usar un comando desconocido de una biblioteca sofisticada que haga las cosas como una caja negra (queremos practicar numpy) * Recorte un sector de la imagen conteniendo una estrella individual y ajuste una gaussiana 2D simétrica a la imagen de la estrella por mínimos cuadrados, incluyendo una constante aditiva (el cielo "vacio" brilla) * Repita este procedimiento para varias estrellas y presente alguna estadística sobre las medidas de la FWHM de las distintas gaussianas: histograma, media, mediana, desviación estándar * Repita el mismo ejercicio sobre cada una de las bandas R,G,B separadamente y comente si observa diferencias en los resultados
import numpy as np
import imageio as iio
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.patches as patches
from scipy.optimize import leastsq
# %matplotlib notebook
La imagen se lee con imageio y se convierte en array con np.asarray
img = np.asarray(iio.imread('data/zapatocaImage.jpeg'))
# -- visualizar la imagen
fig, ax = plt.subplots(figsize=(20, 7))
ax.imshow(img)
plt.show()
# -- shape de la imagen
# -- el ultimo numero indica el numero de canales, al ser una imagen hace referencia a RedGreenBlue
img.shape
(789, 1184, 3)
Para combinar los canales, vamos a promediar los tres canales para así "crear" la imagen en escala de grises.
# -- promedio de la imagen en axis = 2 (la dimension de los canales), el tercer numero de la tupla
# -- ahora, no hay tercera dimension
gr_img = img.mean(axis=2)
gr_img.shape
(789, 1184)
# -- visualizar la imagen en escala de grises
fig, ax = plt.subplots(figsize=(20, 7))
# -- usar un mapa de color en grises
ax.imshow(gr_img, cmap ='gist_gray')
plt.show()
fig, ax = plt.subplots(figsize=(20, 7))
ax.imshow(gr_img, cmap ='gist_gray')
rect1 = patches.Rectangle((889-5, 283-5), 15, 10, linewidth=1, edgecolor='r', facecolor='none')
ax.text(889-8,283-10,'1',fontsize=10, c='orange')
rect2 = patches.Rectangle((946-5, 250-5), 20, 20, linewidth=1, edgecolor='r', facecolor='none')
ax.text(946-8,250-10,'2',fontsize=10, c='orange')
rect3 = patches.Rectangle((617, 308-5), 20, 20, linewidth=1, edgecolor='r', facecolor='none')
ax.text(617,308-10,'3',fontsize=10, c='orange')
rect4 = patches.Rectangle((605, 16-5), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(605+18,16,'4',fontsize=10, c='orange')
rect5 = patches.Rectangle((542, 347-5), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(542,347-9,'5',fontsize=10, c='orange')
rect6 = patches.Rectangle((305, 295-5), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(305,295-9,'6',fontsize=10, c='orange')
rect7 = patches.Rectangle((386, 141-5), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(386,141-9,'7',fontsize=10, c='orange')
rect8 = patches.Rectangle((320, 458), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(320,458-9,'8',fontsize=10, c='orange')
rect9 = patches.Rectangle((1071-2, 277), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(1071-2,277-9,'9',fontsize=10, c='orange')
rect10 = patches.Rectangle((882, 42-9), 15, 15, linewidth=1, edgecolor='r', facecolor='none')
ax.text(882,42-9,'10',fontsize=10, c='orange')
ax.add_patch(rect1)
ax.add_patch(rect2)
ax.add_patch(rect3)
ax.add_patch(rect4)
ax.add_patch(rect5)
ax.add_patch(rect6)
ax.add_patch(rect7)
ax.add_patch(rect8)
ax.add_patch(rect9)
ax.add_patch(rect10)
plt.show()
Función gaussiana generalizada en dos dimensiones
$f(x,y) = A \exp\left(- \left(a(x - x_o)^2 + 2b(x-x_o)(y-y_o) + c(y-y_o)^2 \right)\right)$
${\displaystyle {\begin{aligned}a&={\frac {\cos ^{2}\theta }{2\sigma _{X}^{2}}}+{\frac {\sin ^{2}\theta }{2\sigma _{Y}^{2}}}\\[4pt]b&=-{\frac {\sin 2\theta }{4\sigma _{X}^{2}}}+{\frac {\sin 2\theta }{4\sigma _{Y}^{2}}}\\[4pt]c&={\frac {\sin ^{2}\theta }{2\sigma _{X}^{2}}}+{\frac {\cos ^{2}\theta }{2\sigma _{Y}^{2}}}\end{aligned}}}$
def function(params, x, y):
'''
Define función gaussiana generalizada en dos dimensiones.
'''
A = params[0]
sigma_x = params[1]
sigma_y = params[2]
x0 = params[3]
y0 = params[4]
theta = params[5]
K = params[6]
a = np.cos(theta)**2/(2*sigma_x**2) + np.sin(theta)**2/(2*sigma_y**2)
b = -np.sin(2*theta)/(4*sigma_x**2) + np.sin(2*theta)/(4*sigma_y**2)
c = np.sin(theta)**2/(2*sigma_x**2) + np.cos(theta)**2/(2*sigma_y**2)
z = A*np.exp( - (a*(x-x0)**2 + 2*b*(x-x0)*(y-y0) + c*(y-y0)**2)) + K
return z
def Error(tpl,x,y,z,incer):
'''
Calcula el error entre el valor z dado por la imagen y el determinado por el método de mínimos cuadrados
o cualquier otro. Al ser el resultado dos dimensional, fue necesario convertirlo en unidimensional.
'''
# x son las posiciones reportadas
# y son los valores medidos en las posiciones x
# tpl es una tupla con los parámetros para calcular el modelo
zmodel = function(tpl,x,y)
errors = (zmodel.flatten() - z.flatten()) / np.sqrt(abs(incer))
return errors
def calculator(star,p0,Error,incer=1):
'''
Extrae las coordenadas x & y de la imagen y calcula la mejor gaussiana mediante el método de mínimos
cuadrados
'''
x = np.arange(0,star.shape[0],1)
y = np.arange(0,star.shape[1],1)
x, y = np.meshgrid(y,x)
best,suss = leastsq(Error, p0, args=(x,y,star,incer))
return best,x,y
# -- listas para guardar los valore del Full Width at Half Maximum
FWHM_x = []
FWHM_y = []
Formúla para determinar el Full width at half maximum
${\mathrm {FWHM}}=2{\sqrt {2\ln 2}}\;\sigma$
def determinar(star,p0, FWHMx = FWHM_x, FWHMy = FWHM_y,incer=1):
'''
Determina la mejor gaussiana dado un punto inicial, usando la función calculator.
Determina el valor del FWHM en dirección x & y, guarda el valor en una lista.
Este es determinado usando la formula de https://en.wikipedia.org/wiki/Full_width_at_half_maximum
Genera una gráfica donde se muestra el contorno de la gaussiana determinada encima de la imagen.
'''
best_star,x,y = calculator(star, p0, Error,incer)
print(best_star)
FWHM_star_x = 2*np.sqrt(2*np.log(2))*best_star[1]
FWHM_star_y = 2*np.sqrt(2*np.log(2))*best_star[2]
FWHMx.append(FWHM_star_x)
FWHMy.append(FWHM_star_y)
print(FWHM_star_x,FWHM_star_y)
# -- plot el contorno de la gaussiana encontrada encima de la imagen de la estrella
fig, ax = plt.subplots(figsize=(20, 7))
ax.imshow(star,cmap='gist_gray')
ax.contour(x, y, function(best_star,x,y), 3, colors='r')
plt.show()
Para cada uno de los casos mostrados, los parámetros iniciales p0 son elegidos convenientemente para que la gaussiana determinada por el método sea la mejor posible.
str_1 = gr_img.copy()
str_1 = str_1[283:288,890:895]
plt.imshow(str_1, cmap ='gist_gray')
plt.show()
p0 = [70, 1.2, 11, 1, 1, 0, 0]
determinar(str_1,p0)
[186.07969768 0.92039049 1.01048277 2.00654537 2.0493886 -0.300277 69.41303048] 2.1673539689071366 2.3795050880398305
str_1
array([[ 76.33333333, 72.33333333, 92. , 89. , 76. ], [ 75.33333333, 123.33333333, 182. , 136. , 70. ], [ 84.33333333, 176.33333333, 247. , 176. , 78. ], [ 82.33333333, 143.33333333, 191. , 141. , 78. ], [ 84.33333333, 80.33333333, 84. , 86. , 83. ]])
# -- extraer valores de cada columna para formar un vector con cada columna
# -- el tercer vector es el que parte la estrella en dos
wa = list(zip(*str_1))
list(wa[2])
[92.0, 182.0, 247.0, 191.0, 84.0]
# -- grafiquemos esos puntos
# -- no se ve tan linda como la del profe
plt.scatter(np.arange(0,len(wa[2]),1),wa[2],c='green')
plt.ylabel('Intensidad del pixel')
plt.xlabel('No del pixel fila')
plt.show()
def function1(params, x):
'''
Define función gaussiana generalizada en una dimension.
'''
A = params[0]
sigma_x = params[1]
#sigma_y = params[2]
x0 = params[2]
#y0 = params[4]
theta = params[3]
K = params[4]
a = np.cos(theta)**2/(2*sigma_x**2) #+ np.sin(theta)**2/(2*sigma_y**2)
#b = -np.sin(2*theta)/(4*sigma_x**2) + np.sin(2*theta)/(4*sigma_y**2)
#c = np.sin(theta)**2/(2*sigma_x**2) + np.cos(theta)**2/(2*sigma_y**2)
z = A*np.exp( - (a*(x-x0)**2 )) + K
return z
def Error1(tpl,x,y):
'''
Calcula el error entre el valor y dado por la imagen y el determinado por el método de mínimos cuadrados
o cualquier otro.
'''
# x son las posiciones reportadas
# y son los valores medidos en las posiciones x
# tpl es una tupla con los parámetros para calcular el modelo
ymodel = function1(tpl,x)
errors = (ymodel - y)
#print(errors)
return errors
p0 = [10, 10, 50, 4,0]
x = np.arange(0,len(wa[2]),1)
y = np.array(wa[2])
best,suss = leastsq(Error1, p0, args=(x,y))
plt.scatter(x,function1(best,x),c='orange',label='Fit')
plt.scatter(x,wa[2], c = 'green',label='Original')
plt.legend()
<matplotlib.legend.Legend at 0x7f63d48ae080>
str_2 = gr_img.copy()
str_2 = str_2[252:257,949:954]
plt.imshow(str_2, cmap='gist_gray')
plt.show()
p0 = [50, 1.3, 2.5, 4, 1, 0.5, 0]
#p0 = [60, 30.4, 11, 0.01, 1, 0, 0]
determinar(str_2,p0)
[199.17204186 1.08795697 1.01652482 2.3507139 2.25572616 -1.63067875 71.03541102] 2.5619428907503394 2.393733033767761
str_3 = gr_img.copy()
str_3 = str_3[309:324,619:636]
plt.imshow(str_3, cmap='gist_gray')
plt.show()
p0 = [50, 15, 10, 4, 1, 0, 0]
determinar(str_3,p0)
[167.1691 2.54202832 2.43955445 8.26104046 6.7936942 -1.57329251 97.72363193] 5.986019246498141 5.744711731315182
str_4 = gr_img.copy()
str_4 = str_4[18:23,607:612]
plt.imshow(str_4, cmap='gist_gray')
plt.show()
p0 = [10, 3, 6, 0.9, 1.4, 0, 0]
determinar(str_4,p0)
[171.67801287 0.72367521 0.95397152 1.9468729 2.26661796 -3.22904682 85.47495094] 1.704124894388063 2.2464312597874616
str_5 = gr_img.copy()
str_5 = str_5[347:356,542:552]
plt.imshow(str_5, cmap='gist_gray')
plt.show()
p0 = [10, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_5,p0)
[164.57072408 1.12967336 1.30573147 4.56858863 4.23103384 6.39463879 115.86139744] 2.660177475704662 3.0747626286144585
str_6 = gr_img.copy()
str_6 = str_6[295:300,307:312]
plt.imshow(str_6, cmap='gist_gray')
plt.show()
p0 = [10, 3.3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6,p0)
[154.74939288 1.06584842 0.86557986 2.09066242 1.67091097 -1.79900438 106.39500057] 2.509881229793308 2.038284809401981
str_7 = gr_img.copy()
str_7 = str_7[142:154,387:399]
plt.imshow(str_7, cmap='gist_gray')
plt.show()
p0 = [10, 2, 5, 3, 5, 0, 0]
determinar(str_7,p0)
[140.94787706 1.79851031 2.03013514 6.17006805 5.63716518 3.68873831 108.41886214] 4.235168131706809 4.780602919231583
str_8 = gr_img.copy()
str_8 = str_8[459:465,322:327]
plt.imshow(str_8, cmap='gist_gray')
plt.show()
p0 = [10, 3.43, 6.5, 0.01, 0.5, 0, 0]
determinar(str_8,p0)
[133.76503279 1.29368943 0.97282715 2.12474538 2.37096544 10.63544949 129.69305459] 3.0464058103768084 2.2908328681370387
str_9 = gr_img.copy()
str_9 = str_9[277:284,1071:1076]
plt.imshow(str_9, cmap='gist_gray')
plt.show()
p0 = [10, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9,p0)
[188.16928656 0.91808841 1.25701937 2.10259749 3.14037841 3.29910297 61.1276509 ] 2.161933002277811 2.96005440206742
str_10 = gr_img.copy()
str_10 = str_10[43:49,883:889]
plt.imshow(str_10, cmap='gist_gray')
plt.show()
p0 = [700, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_10,p0)
[2.12780258e+02 9.98780282e-01 1.13708459e+00 2.28269387e+00 2.89718094e+00 1.60170734e-03 6.12198320e+01] 2.3519478278094055 2.677629596603565
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_x = abs(np.array(FWHM_x))
print(FWHM_x)
stdx = FWHM_x.std()
print(stdx)
meanx = FWHM_x.mean()
print(meanx)
medianx = np.median(FWHM_x)
print(medianx)
plt.hist(FWHM_x.flatten())
plt.title('Valores FHWM en la dirección x para escala de grises')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.16735397 2.56194289 5.98601925 1.70412489 2.66017748 2.50988123 4.23516813 3.04640581 2.161933 2.35194783] 1.2020310334918178 2.9384954478212486 2.5359120602718237
FWHM_y = abs(np.array(FWHM_y))
print(FWHM_y)
stdy = FWHM_y.std()
print(stdy)
meany = FWHM_y.mean()
print(meany)
mediany = np.median(FWHM_y)
print(mediany)
plt.hist(FWHM_y.flatten())
plt.title('Valores FHWM en la dirección y para escala de grises')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.37950509 2.39373303 5.74471173 2.24643126 3.07476263 2.03828481 4.78060292 2.29083287 2.9600544 2.6776296 ] 1.1631400032813368 3.0586548336966284 2.5356813151856628
Para extraer el canal rojo, usamos todas las filas, todas las columnas, pero solo el primero de los tres canales
$img[:,:,0]$
FWHM_xr = []
FWHM_yr = []
str_1r = img[:,:,0].copy()
str_1r = str_1r[283:288,890:895]
plt.imshow(str_1r, cmap ='Reds')
plt.show()
p0 = [70, 1.2, 11, 1, 1, 0, 0]
determinar(str_1r,p0,FWHM_xr,FWHM_yr)
[178.76060595 0.9511361 1.04253055 2.00150507 2.05170422 -0.30877778 88.35069654] 2.2397543569080045 2.454971833929507
str_2r = img[:,:,0].copy()
str_2r = str_2r[252:257,949:955]
plt.imshow(str_2r, cmap='Reds')
plt.show()
p0 = [50, 1.3, 2.5, 4, 1, 0, 0]
determinar(str_2r,p0,FWHM_xr,FWHM_yr)
[190.6122258 1.09953422 0.9819177 2.34511946 2.27447217 -1.62016404 85.47014377] 2.589205231289726 2.312239480960313
str_3r = img[:,:,0].copy()
str_3r = str_3r[309:324,619:636]
plt.imshow(str_3r, cmap='Reds')
plt.show()
p0 = [50, 15, 10, 4, 1, 0, 0]
determinar(str_3r,p0,FWHM_xr,FWHM_yr)
[143.67174102 2.32328961 2.16090506 8.28633162 6.8475628 -1.5779969 114.72814395] 5.470928949222188 5.0885425478166635
str_4r = img[:,:,0].copy()
str_4r = str_4r[18:23,607:612]
plt.imshow(str_4r, cmap='Reds')
plt.show()
p0 = [10, 3, 6, 0.9, 1.4, 0, 0]
determinar(str_4r,p0,FWHM_xr,FWHM_yr)
[170.51074692 0.72404326 0.94899811 1.94940078 2.28406376 6.190567 100.11931426] 1.704991580249038 2.2347197813928976
str_5r = img[:,:,0].copy()
str_5r = str_5r[347:356,542:552]
plt.imshow(str_5r, cmap='Reds')
plt.show()
p0 = [10, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_5r,p0,FWHM_xr,FWHM_yr)
[159.16045587 1.07357801 1.27234105 4.5756869 4.21360054 6.28738116 128.40525837] 2.528083008621747 2.9961342032205103
str_6r = img[:,:,0].copy()
str_6r = str_6r[295:300,307:312]
plt.imshow(str_6r, cmap='Reds')
plt.show()
p0 = [10, 3.3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6r,p0,FWHM_xr,FWHM_yr)
[149.35942882 1.0738311 0.87541102 2.08507806 1.68940448 -1.81607981 124.90198335] 2.528679002803907 2.061435424259867
str_7r = img[:,:,0].copy()
str_7r = str_7r[142:154,387:399]
plt.imshow(str_7r, cmap='Reds')
plt.show()
p0 = [10, 2, 5, 3, 5, 0, 0]
determinar(str_7r,p0,FWHM_xr,FWHM_yr)
[138.7387436 1.9386357 -2.22391119 6.19193906 5.64616817 0.47421667 129.67231319] 4.565138195369105 -5.236910658145347
str_8r = img[:,:,0].copy()
str_8r = str_8r[459:465,322:327]
plt.imshow(str_8r, cmap='Reds')
plt.show()
p0 = [10, 3.43, 6.4, 0.01, 0.5, 0, 0]
determinar(str_8r,p0,FWHM_xr,FWHM_yr)
[117.96132322 0.8918261 1.23916011 2.140237 2.33663264 5.9503522 150.18634797] 2.1000899865145097 2.9179990671423646
str_9r = img[:,:,0].copy()
str_9r = str_9r[277:284,1071:1076]
plt.imshow(str_9r, cmap='Reds')
plt.show()
p0 = [10, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9r,p0,FWHM_xr,FWHM_yr)
[182.9258356 0.94142815 -1.29220082 2.10748647 3.15148273 3.294021 80.34963419] 2.2168938787115526 -3.042900394987331
str_10r = img[:,:,0].copy()
str_10r = str_10r[43:49,883:889]
plt.imshow(str_10r, cmap='Reds')
plt.show()
p0 = [700, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_10r,p0,FWHM_xr,FWHM_yr)
[2.09299250e+02 1.01842846e+00 1.15336547e+00 2.28048409e+00 2.89212867e+00 1.67741306e-02 7.77637647e+01] 2.398215753935267 2.7159681216478253
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xr = abs(np.array(FWHM_xr))
print(FWHM_xr)
stdxr = FWHM_xr.std()
print(stdxr)
meanxr = FWHM_xr.mean()
print(meanxr)
medianxr = np.median(FWHM_xr)
print(medianxr)
plt.hist(FWHM_xr.flatten())
plt.title('Valores FHWM en la dirección x para la image roja')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.23975436 2.58920523 5.47092895 1.70499158 2.52808301 2.528679 4.5651382 2.10008999 2.21689388 2.39821575] 1.1372820062135576 2.8341979943625044 2.463149381278507
FWHM_yr = abs(np.array(FWHM_yr))
print(FWHM_yr)
stdyr = FWHM_yr.std()
print(stdyr)
meanyr = FWHM_yr.mean()
print(meanyr)
medianyr = np.median(FWHM_yr)
print(medianyr)
plt.hist(FWHM_yr.flatten())
plt.title('Valores FHWM en la dirección y para la image roja')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.45497183 2.31223948 5.08854255 2.23471978 2.9961342 2.06143542 5.23691066 2.91799907 3.04290039 2.71596812] 1.0760536200525121 3.1061821513502625 2.816983594395095
Para extraer el canal verde, usamos todas las filas, todas las columnas, pero solo el segundo de los tres canales
$img[:,:,1]$
FWHM_xg = []
FWHM_yg = []
str_1g = img[:,:,1].copy()
str_1g = str_1g[283:288,890:895]
plt.imshow(str_1g, cmap ='Greens')
plt.show()
p0 = [70, 1.5, 11, 1, 1, 0, 0]
determinar(str_1g,p0,FWHM_xg,FWHM_yg)
[189.66980515 0.90810494 0.99594246 2.00139712 2.04842509 -0.31654792 67.1450173 ] 2.138423711514054 2.3452652626833705
str_2g = img[:,:,1].copy()
str_2g = str_2g[252:257,949:955]
plt.imshow(str_2g, cmap='Greens')
plt.show()
p0 = [50, 1.3, 2.5, 4, 1, 0, 0]
determinar(str_2g,p0,FWHM_xg,FWHM_yg)
[200.49323003 1.0742782 1.00605482 2.34595298 2.25284053 -1.63474596 71.37063765] 2.5297318424914206 2.3690780599175865
str_3g = img[:,:,1].copy()
str_3g = str_3g[309:324,619:636]
plt.imshow(str_3g, cmap='Greens')
plt.show()
p0 = [50, 15, 10, 4, 1, 0, 0]
determinar(str_3g,p0,FWHM_xg,FWHM_yg)
[170.3427937 2.50270952 2.42460793 8.20606974 6.75458075 -1.59836573 96.7865555 ] 5.893430538924626 5.709515350529837
str_4g = img[:,:,1].copy()
str_4g = str_4g[18:23,607:612]
plt.imshow(str_4g, cmap='Greens')
plt.show()
p0 = [10, 4, 6, 0.9, 1.4, 0, 0]
determinar(str_4g,p0,FWHM_xg,FWHM_yg)
[171.70162698 0.72207629 0.95112594 1.94675696 2.26657667 9.33559274 84.76461655] 1.7003597200006297 2.2397304322878533
str_5g = img[:,:,1].copy()
str_5g = str_5g[347:356,542:552]
plt.imshow(str_5g, cmap='Greens')
plt.show()
p0 = [10, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_5g,p0,FWHM_xg,FWHM_yg)
[168.26372326 1.10033428 1.29450467 4.5902933 4.22695556 9.45794583 113.51970684] 2.5910892262579335 3.048325549295308
str_6g = img[:,:,1].copy()
str_6g = str_6g[295:300,307:312]
plt.imshow(str_6g, cmap='Greens')
plt.show()
p0 = [10, 3.5, 6, 0.78, 1.7, 0, 0]
determinar(str_6g,p0,FWHM_xg,FWHM_yg)
[156.84218534 1.06030104 0.85529269 2.10415204 1.66655508 -1.79293087 104.3264753 ] 2.4968181499130977 2.0140603781124353
str_7g = img[:,:,1].copy()
str_7g = str_7g[142:154,387:399]
plt.imshow(str_7g, cmap='Greens')
plt.show()
p0 = [15, 2, 5, 3, 5, 0, 0]
determinar(str_7g,p0,FWHM_xg,FWHM_yg)
[145.65432612 1.76651186 1.99543641 6.18715013 5.63414873 0.57532504 107.4753585 ] 4.159817530728484 4.6988936558118475
str_8g = img[:,:,1].copy()
str_8g = str_8g[459:465,322:327]
plt.imshow(str_8g, cmap='Greens')
plt.show()
p0 = [10, 3.43, 6.4, 1, 0.5, 0, 0]
determinar(str_8g,p0,FWHM_xg,FWHM_yg)
[132.62532774 0.95533301 1.27035965 2.15830574 2.3619008 5.91555374 128.72993202] 2.249637328620205 2.9914683602476924
str_9g = img[:,:,1].copy()
str_9g = str_9g[277:284,1071:1076]
plt.imshow(str_9g, cmap='Greens')
plt.show()
p0 = [10, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9g,p0,FWHM_xg,FWHM_yg)
[191.01057596 0.90731452 1.24170476 2.10161186 3.14651353 3.29961415 61.2009776 ] 2.1365624274379273 2.923991268052027
str_10g = img[:,:,1].copy()
str_10g = str_10g[43:49,883:889]
plt.imshow(str_10g, cmap='Greens')
plt.show()
p0 = [700, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_10g,p0,FWHM_xg,FWHM_yg)
[ 2.16319665e+02 9.99190581e-01 1.12653091e+00 2.27946151e+00 2.90122310e+00 -1.93991955e-02 6.06937104e+01] 2.3529140096920815 2.6527775769756388
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xg = abs(np.array(FWHM_xg))
print(FWHM_xg)
stdxg = FWHM_xg.std()
print(stdxr)
meanxg = FWHM_xg.mean()
print(meanxg)
medianxg = np.median(FWHM_xg)
print(medianxg)
plt.hist(FWHM_xg.flatten())
plt.title('Valores FHWM en la dirección x para la imagen verde')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.13842371 2.52973184 5.89343054 1.70035972 2.59108923 2.49681815 4.15981753 2.24963733 2.13656243 2.35291401] 1.1372820062135576 2.824878448558046 2.4248660798025896
FWHM_yg = abs(np.array(FWHM_yg))
print(FWHM_yg)
stdyg = FWHM_yg.std()
print(stdyg)
meanyg = FWHM_yg.mean()
print(meanyg)
medianyg = np.median(FWHM_yg)
print(medianyg)
plt.hist(FWHM_yg.flatten())
plt.title('Valores FHWM en la dirección y para la imagen verde')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.34526526 2.36907806 5.70951535 2.23973043 3.04832555 2.01406038 4.69889366 2.99146836 2.92399127 2.65277758] 1.124040162412368 3.0993105893913597 2.7883844225138326
Para extraer el canal azul, usamos todas las filas, todas las columnas, pero solo el último de los tres canales
$img[:,:,2]$
FWHM_xb = []
FWHM_yb = []
str_1b = img[:,:,2].copy()
str_1b = str_1b[283:288,890:895]
plt.imshow(str_1b, cmap ='Blues')
plt.show()
p0 = [70, 1.7, 11, 1, 1, 0, 0]
determinar(str_1b,p0,FWHM_xb,FWHM_yb)
[190.2169832 0.90335067 0.99435099 2.01610059 2.04808311 -0.27729555 52.66186116] 2.1272282649036094 2.3415176386132233
str_2b = img[:,:,2].copy()
str_2b = str_2b[252:257,949:955]
plt.imshow(str_2b, cmap='Blues')
plt.show()
p0 = [50, 2.9, 2, 4, 1, 0, 0]
determinar(str_2b,p0,FWHM_xb,FWHM_yb)
[205.9355405 1.07443705 1.03024455 2.35673966 2.24362524 1.41576639 58.80752854] 2.5301059117924476 2.42604050765156
str_3b = img[:,:,2].copy()
str_3b = str_3b[309:324,619:636]
plt.imshow(str_3b, cmap='Blues')
plt.show()
p0 = [50, 15, 10, 4, 1, 0, 0]
determinar(str_3b,p0,FWHM_xb,FWHM_yb)
[191.02580433 2.73160061 2.64058181 8.29780833 6.77398083 -1.59662047 81.16874448] 6.432427878391856 6.218094970738713
str_4b = img[:,:,2].copy()
str_4b = str_4b[18:23,607:612]
plt.imshow(str_4b, cmap='Blues')
plt.show()
p0 = [9, 3.4, 6, 1.23, 1.46, 0, 0]
determinar(str_4b,p0,FWHM_xb,FWHM_yb)
[172.85830022 0.7251488 0.96189654 1.94448225 2.24924495 -6.36452675 71.51882385] 1.707594921415306 2.265093259490133
str_5b = img[:,:,2].copy()
str_5b = str_5b[347:356,542:552]
plt.imshow(str_5b, cmap='Blues')
plt.show()
p0 = [10, 12, 10, 2, 1, 0, 0]
determinar(str_5b,p0,FWHM_xb,FWHM_yb)
[165.9928307 1.21438496 1.37811936 4.53438415 4.24417785 3.5362753 105.2004878 ] 2.8596580347727287 3.2452230957557635
str_6b = img[:,:,2].copy()
str_6b = str_6b[295:300,307:312]
plt.imshow(str_6b, cmap='Blues')
plt.show()
p0 = [10, 3.3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6b,p0,FWHM_xb,FWHM_yb)
[158.14516944 1.06375236 0.86631248 2.08220315 1.65755899 -1.78824612 89.94146728] 2.5049453718991064 2.0400099880751132
str_7b = img[:,:,2].copy()
str_7b = str_7b[142:154,387:399]
plt.imshow(str_7b, cmap='Blues')
plt.show()
p0 = [11, 2, 5, 3, 5, 0, 0]
determinar(str_7b,p0,FWHM_xb,FWHM_yb)
[140.2375672 1.69057103 1.87880031 6.13700838 5.63296036 0.59694216 87.61531461] 3.9809905380035873 4.424236622195258
str_8b = img[:,:,2].copy()
str_8b = str_8b[459:465,322:327]
plt.imshow(str_8b, cmap='Blues')
plt.show()
p0 = [10, 3.43, 6.4, 0.01, 0.5, 0, 0]
determinar(str_8b,p0,FWHM_xb,FWHM_yb)
[152.46960483 1.37726272 1.08272843 2.08071456 2.39950791 -1.95926817 107.90558619] 3.2432058673937094 2.5496306060688854
str_9b = img[:,:,2].copy()
str_9b = str_9b[277:284,1071:1076]
plt.imshow(str_9b, cmap='Blues')
plt.show()
p0 = [13, 3, 6.4, 1, 1.1, 0, 0]
determinar(str_9b,p0,FWHM_xb,FWHM_yb)
[190.8831359 0.90622759 1.23752934 2.09884658 3.12378005 3.30383272 41.79158734] 2.1340028974084926 2.9141588888311842
str_10b = img[:,:,2].copy()
str_10b = str_10b[43:49,883:889]
plt.imshow(str_10b, cmap='Blues')
plt.show()
p0 = [50, 2, 4, 1, 1, 0, 0]
determinar(str_10b,p0,FWHM_xb,FWHM_yb)
[2.12869194e+02 9.79245215e-01 1.13164829e+00 2.28753325e+00 2.89802157e+00 6.97495764e-03 4.51627862e+01] 2.305946260209862 2.6648280839211322
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xb = abs(np.array(FWHM_xb))
print(FWHM_xb)
stdxb = FWHM_xb.std()
print(stdxb)
meanxb = FWHM_xb.mean()
print(meanxb)
medianxb = np.median(FWHM_xb)
print(medianxb)
plt.hist(FWHM_xb.flatten())
plt.title('Valores FHWM en la dirección x para la image azul')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.12722826 2.53010591 6.43242788 1.70759492 2.85965803 2.50494537 3.98099054 3.24320587 2.1340029 2.30594626] 1.3018950466931847 2.9826105946190706 2.517525641845777
FWHM_yb = abs(np.array(FWHM_yb))
print(FWHM_yb)
stdyb = FWHM_yb.std()
print(stdyb)
meanyb = FWHM_yb.mean()
print(meanyb)
medianyb = np.median(FWHM_yb)
print(medianyb)
plt.hist(FWHM_yb.flatten())
plt.title('Valores FHWM en la dirección y para la image azul')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.34151764 2.42604051 6.21809497 2.26509326 3.2452231 2.04000999 4.42423662 2.54963061 2.91415889 2.66482808] 1.2195655968098804 3.1088833661340964 2.607229344995009
En dirección x
print(stdx,stdxr,stdxg,stdxb)
print(meanx,meanxr,meanxg,meanxb)
print(medianx,medianxr,medianxg,medianxb)
1.2020310334918178 1.1372820062135576 1.1920847689204435 1.3018950466931847 2.9384954478212486 2.8341979943625044 2.824878448558046 2.9826105946190706 2.5359120602718237 2.463149381278507 2.4248660798025896 2.517525641845777
En dirección y
print(stdy,stdyr,stdyg,stdyb)
print(meany,meanyr,meanyg,meanyb)
print(mediany,medianyr,medianyg,medianyb)
1.1631400032813368 1.0760536200525121 1.124040162412368 1.2195655968098804 3.0586548336966284 3.1061821513502625 3.0993105893913597 3.1088833661340964 2.5356813151856628 2.816983594395095 2.7883844225138326 2.607229344995009
El valor el FWHM es mas o menos el mismo para todas las estrellas, entre 2 y 3, en algunos casos es mayor y en otros menor.
Para la escala de grises y para los canales, debí cambiar el punto inicial p0 para que la gaussiana fuera mas apropiada, e.g. para la Estrella n, debí usar cuatro valores de p0, uno para cada canal y para la gris. También elegí los valores de p0 de tal manera que los valores de sigma_x y sigma_y fueran siempre positivos.
El valor para el FWHM si cambia dependiendo del canal del que se parta para hacer el ajuste. Mirando los valores de la desviación estándar, el promedio y la mediana se ve que hay diferencias, sin embargo, los valores podrian considerarse como iguales para los cuatro casos. Igual teniendo en cuenta que el valor inicial p0 para cada estrella debio ser modificado para cada canal.
# -- listas para guardar los valore del Full Width at Half Maximum
FWHM_xi = []
FWHM_yi = []
str_1 = gr_img.copy()
str_1 = str_1[283:288,890:895]
plt.imshow(str_1, cmap ='gist_gray')
plt.show()
p0 = [71, 1.2, 11, 1, 1, 0, 0]
determinar(str_1,p0,FWHM_xi,FWHM_yi,incer=str_1.flatten())
[191.03306608 0.88416723 0.97198929 2.00370142 2.04447666 -0.30033141 70.9737103 ] 2.0820547263026454 2.28885986684923
str_2 = gr_img.copy()
str_2 = str_2[252:257,949:954]
plt.imshow(str_2, cmap='gist_gray')
plt.show()
p0 = [57, 1.3, 2.5, 4, 1, 0.5, 0]
#p0 = [60, 30.4, 11, 0.01, 1, 0, 0]
determinar(str_2,p0,FWHM_xi,FWHM_yi,incer=str_2.flatten())
[200.04931224 1.06771002 1.01254502 2.37072665 2.22891783 -1.63770305 71.38762112] 2.514264960806871 2.3843613096910805
str_3 = gr_img.copy()
str_3 = str_3[309:324,619:636]
plt.imshow(str_3, cmap='gist_gray')
plt.show()
p0 = [70, 7, 10, 4, 1, 0, 0]
determinar(str_3,p0,FWHM_xi,FWHM_yi,incer=str_3.flatten())
[165.09149595 2.56740779 2.48379149 8.25953555 6.79491264 1.52624294 97.1682973 ] 6.045783318082368 5.848881983108898
str_4 = gr_img.copy()
str_4 = str_4[18:23,607:612]
plt.imshow(str_4, cmap='gist_gray')
plt.show()
p0 = [11, 3, 6, 0.9, 1.4, 0, 0]
determinar(str_4,p0,FWHM_xi,FWHM_yi,incer=str_4.flatten())
[181.39720161 0.73469636 0.98245806 1.94925035 2.25033059 3.05248272 78.45962331] 1.7300777056892003 2.3135119445798464
str_5 = gr_img.copy()
str_5 = str_5[347:356,542:552]
plt.imshow(str_5, cmap='gist_gray')
plt.show()
p0 = [12, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_5,p0,FWHM_xi,FWHM_yi,incer=str_5.flatten())
[166.83763898 1.12850905 1.28552333 4.56532092 4.21399213 3.32990317 115.39600089] 2.657435737167189 3.0271761137992788
str_6 = gr_img.copy()
str_6 = str_6[295:300,307:312]
plt.imshow(str_6, cmap='gist_gray')
plt.show()
p0 = [70, 2.3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6,p0,FWHM_xi,FWHM_yi,incer=str_6.flatten())
[155.79025228 0.85378903 1.05568553 2.08202815 1.70316379 -0.22008333 106.48341521] 2.0105195108028022 2.4859494493339502
str_7 = gr_img.copy()
str_7 = str_7[142:154,387:399]
plt.imshow(str_7, cmap='gist_gray')
plt.show()
p0 = [10, 2, 5, 3, 5, 0, 0]
determinar(str_7,p0,FWHM_xi,FWHM_yi,incer=str_7.flatten())
[138.57467342 1.84344675 2.08148352 6.17117586 5.62605279 0.63536992 107.33799354] 4.340985359062775 4.901519126606775
str_8 = gr_img.copy()
str_8 = str_8[459:465,322:327]
plt.imshow(str_8, cmap='gist_gray')
plt.show()
p0 = [10, 3.43, 6.5, 0.01, 0.5, 0, 0]
determinar(str_8,p0,FWHM_xi,FWHM_yi,incer=str_8.flatten())
[133.55992457 0.98677464 1.29549702 2.12033093 2.36053394 2.75037941 129.10974029] 2.3236767077234046 3.050662357379447
str_9 = gr_img.copy()
str_9 = str_9[277:284,1071:1076]
plt.imshow(str_9, cmap='gist_gray')
plt.show()
p0 = [10, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9,p0,FWHM_xi,FWHM_yi,incer=str_9.flatten())
[190.85299257 1.24748929 0.89798364 2.10782787 3.12535911 1.72494314 61.26764694] 2.937612779019202 2.1145898720577216
str_10 = gr_img.copy()
str_10 = str_10[43:49,883:889]
plt.imshow(str_10, cmap='gist_gray')
plt.show()
p0 = [43, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_10,p0,FWHM_xi,FWHM_yi,incer=str_10.flatten())
[219.64233685 1.0717617 0.97269196 2.28430393 2.86515034 -1.59071136 62.34114315] 2.523805945849836 2.2905145186106575
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xi = abs(np.array(FWHM_xi))
print(FWHM_xi)
stdxi = FWHM_xi.std()
print(stdxi)
meanxi = FWHM_xi.mean()
print(meanxi)
medianxi = np.median(FWHM_xi)
print(medianxi)
plt.hist(FWHM_xi.flatten())
plt.title('Valores FHWM en la dirección x para escala de grises considerando incertidumbre')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.08205473 2.51426496 6.04578332 1.73007771 2.65743574 2.01051951 4.34098536 2.32367671 2.93761278 2.52380595] 1.2442375593713781 2.916621675050629 2.519035453328353
FWHM_yi = abs(np.array(FWHM_yi))
print(FWHM_yi)
stdyi = FWHM_yi.std()
print(stdyi)
meanyi = FWHM_yi.mean()
print(meanyi)
medianyi = np.median(FWHM_yi)
print(medianyi)
plt.hist(FWHM_yi.flatten())
plt.title('Valores FHWM en la dirección y para escala de grises considerando incertidumbre')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.28885987 2.38436131 5.84888198 2.31351194 3.02717611 2.48594945 4.90151913 3.05066236 2.11458987 2.29051452] 1.2080118897663938 3.0706026542016884 2.4351553795125156
FWHM_xir = []
FWHM_yir = []
str_1r = img[:,:,0].copy()
str_1r = str_1r[283:288,890:895]
plt.imshow(str_1r, cmap ='Reds')
plt.show()
p0 = [70, 1.2, 11, 1, 1, 0, 0]
determinar(str_1r,p0,FWHM_xir,FWHM_yir, incer = str_1r.flatten())
[182.91431936 0.91235058 1.00060875 1.99849318 2.04670046 -0.31143976 90.21305084] 2.1484214221733677 2.3562535363489063
str_2r = img[:,:,0].copy()
str_2r = str_2r[252:257,949:955]
plt.imshow(str_2r, cmap='Reds')
plt.show()
p0 = [24, 1.3, 2.5, 4, 1, 0, 0]
determinar(str_2r,p0,FWHM_xir,FWHM_yir, incer = str_2r.flatten())
[192.45744813 0.97887011 1.08961655 2.35993751 2.25608062 -6.34642439 84.92551405] 2.3050629462193077 2.5658509037264716
str_3r = img[:,:,0].copy()
str_3r = str_3r[309:324,619:636]
plt.imshow(str_3r, cmap='Reds')
plt.show()
p0 = [11, 15, 9, 4, 1, 0, 0]
determinar(str_3r,p0,FWHM_xir,FWHM_yir, incer = str_3r.flatten())
[141.9874119 2.34167953 2.19328524 8.28237019 6.86111281 20.41751836 114.33748466] 5.514233897956943 5.16479203779562
str_4r = img[:,:,0].copy()
str_4r = str_4r[18:23,607:612]
plt.imshow(str_4r, cmap='Reds')
plt.show()
p0 = [12, 3, 6, 0.9, 1.4, 0, 0]
determinar(str_4r,p0,FWHM_xir,FWHM_yir, incer = str_4r.flatten())
[179.47671572 0.73445207 0.9714805 1.95111908 2.27151171 3.04648514 93.82977231] 1.729502460722425 2.287661764280814
str_5r = img[:,:,0].copy()
str_5r = str_5r[347:356,542:552]
plt.imshow(str_5r, cmap='Reds')
plt.show()
p0 = [15, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_5r,p0,FWHM_xir,FWHM_yir, incer = str_5r.flatten())
[162.50745202 1.0652327 1.2441567 4.57620325 4.20175448 3.18754763 128.05035619] 2.5084313108756726 2.929765136330512
str_6r = img[:,:,0].copy()
str_6r = str_6r[295:300,307:312]
plt.imshow(str_6r, cmap='Reds')
plt.show()
p0 = [10, 3.3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6r,p0,FWHM_xir,FWHM_yir, incer = str_6r.flatten())
[150.54027437 1.06288385 0.86245768 2.07589231 1.71659829 -1.80244502 125.04225717] 2.5029001928382875 2.0309326365113836
str_7r = img[:,:,0].copy()
str_7r = str_7r[142:154,387:399]
plt.imshow(str_7r, cmap='Reds')
plt.show()
p0 = [14, 10, 15, 1, 5, 0, 0]
determinar(str_7r,p0,FWHM_xir,FWHM_yir, incer = str_7r.flatten())
[137.08273591 1.97954232 2.26715815 6.20129932 5.63488281 0.52724304 128.6534564 ] 4.66146594576438 5.338749467405802
str_8r = img[:,:,0].copy()
str_8r = str_8r[459:465,322:327]
plt.imshow(str_8r, cmap='Reds')
plt.show()
p0 = [11, 3.43, 6.4, 0.01, 0.5, 0, 0]
determinar(str_8r,p0,FWHM_xir,FWHM_yir, incer = str_8r.flatten())
[118.37738938 1.2301783 0.88804317 2.13659917 2.33042575 4.36580063 150.26652019] 2.896848514525203 2.091181859056181
str_9r = img[:,:,0].copy()
str_9r = str_9r[277:284,1071:1076]
plt.imshow(str_9r, cmap='Reds')
plt.show()
p0 = [13, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9r,p0,FWHM_xir,FWHM_yir, incer = str_9r.flatten())
[185.59903905 1.27508233 0.91857923 2.11100885 3.13854179 1.72086826 80.86325465] 3.0025894335728154 2.163088772629732
str_10r = img[:,:,0].copy()
str_10r = str_10r[43:49,883:889]
plt.imshow(str_10r, cmap='Reds')
plt.show()
p0 = [700, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_10r,p0,FWHM_xir,FWHM_yir, incer = str_10r.flatten())
[214.94414166 1.08939095 0.9932327 2.27855954 2.86276987 1.57713906 79.07029388] 2.565319646114155 2.338884281994375
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xir = abs(np.array(FWHM_xir))
print(FWHM_xir)
stdxir = FWHM_xir.std()
print(stdxir)
meanxir = FWHM_xir.mean()
print(meanxir)
medianxir = np.median(FWHM_xir)
print(medianxir)
plt.hist(FWHM_xir.flatten())
plt.title('Valores FHWM en la dirección x para la image roja')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.14842142 2.30506295 5.5142339 1.72950246 2.50843131 2.50290019 4.66146595 2.89684851 3.00258943 2.56531965] 1.1220401417419346 2.9834775770762554 2.536875478494914
FWHM_yir = abs(np.array(FWHM_yir))
print(FWHM_yir)
stdyir = FWHM_yir.std()
print(stdyir)
meanyir = FWHM_yir.mean()
print(meanyir)
medianyir = np.median(FWHM_yir)
print(medianyir)
plt.hist(FWHM_yir.flatten())
plt.title('Valores FHWM en la dirección y para la image roja')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.35625354 2.5658509 5.16479204 2.28766176 2.92976514 2.03093264 5.33874947 2.09118186 2.16308877 2.33888428] 1.1882827213448535 2.9267160396079794 2.347568909171641
FWHM_xig = []
FWHM_yig = []
str_1g = img[:,:,1].copy()
str_1g = str_1g[283:288,890:895]
plt.imshow(str_1g, cmap ='Greens')
plt.show()
p0 = [71, 1.5, 11, 1, 1, 0, 0]
determinar(str_1g,p0,FWHM_xig,FWHM_yig, incer = str_1g.flatten())
[194.58735935 0.87499245 0.96168531 1.99756044 2.04355755 -0.315852 68.48717297] 2.060449762326214 2.2645958531909964
str_2g = img[:,:,1].copy()
str_2g = str_2g[252:257,949:955]
plt.imshow(str_2g, cmap='Greens')
plt.show()
p0 = [22, 1.3, 2.5, 4, 1, 0, 0]
determinar(str_2g,p0,FWHM_xig,FWHM_yig, incer = str_2g.flatten())
[202.16491005 1.00204564 1.05684175 2.36233933 2.22523495 -3.24379486 71.25427427] 2.35963715103674 2.4886721282798
str_3g = img[:,:,1].copy()
str_3g = str_3g[309:324,619:636]
plt.imshow(str_3g, cmap='Greens')
plt.show()
p0 = [50, 15, 10, 4, 1, 0, 0]
determinar(str_3g,p0,FWHM_xig,FWHM_yig, incer = str_3g.flatten())
[168.35593595 2.52518256 2.46583528 8.20589042 6.75978717 -1.67205716 96.27653287] 5.946350509252893 5.806598337482345
str_4g = img[:,:,1].copy()
str_4g = str_4g[18:23,607:612]
plt.imshow(str_4g, cmap='Greens')
plt.show()
p0 = [14, 3.4, 6, 0.9, 1.4, 0, 0]
determinar(str_4g,p0,FWHM_xig,FWHM_yig, incer = str_4g.flatten())
[181.55452689 0.73372479 0.98023727 1.94887108 2.25037477 3.05065648 77.61482048] 1.7277898473442026 2.3082823783862074
str_5g = img[:,:,1].copy()
str_5g = str_5g[347:356,542:552]
plt.imshow(str_5g, cmap='Greens')
plt.show()
p0 = [17, 3.4, 6, 0.8, 1.4, 0, 0]
determinar(str_5g,p0,FWHM_xig,FWHM_yig, incer = str_5g.flatten())
[171.08439474 1.09751253 1.27162261 4.58951748 4.21422199 3.22642819 113.06046966] 2.5844445024818627 2.994442416716548
str_6g = img[:,:,1].copy()
str_6g = str_6g[295:300,307:312]
plt.imshow(str_6g, cmap='Greens')
plt.show()
p0 = [32, 3.5, 6, 0.78, 1.7, 0, 0]
determinar(str_6g,p0,FWHM_xig,FWHM_yig, incer = str_6g.flatten())
[157.82945535 0.84466944 1.05013524 2.09746462 1.69957051 -0.21827805 104.40660141] 1.9890445292987982 2.4728795039574236
str_7g = img[:,:,1].copy()
str_7g = str_7g[142:154,387:399]
plt.imshow(str_7g, cmap='Greens')
plt.show()
p0 = [11, 2, 5, 3, 5, 0, 0]
determinar(str_7g,p0,FWHM_xig,FWHM_yig, incer = str_7g.flatten())
[142.44945667 1.81945909 2.06026126 6.18751642 5.62400824 0.67907496 106.22773724] 4.284498726614627 4.851544509614211
str_8g = img[:,:,1].copy()
str_8g = str_8g[459:465,322:327]
plt.imshow(str_8g, cmap='Greens')
plt.show()
p0 = [16, 3.43, 6.4, 1, 0.5, 0, 0]
determinar(str_8g,p0,FWHM_xig,FWHM_yig, incer = str_8g.flatten())
[132.27002448 0.97039729 1.27318881 2.15941699 2.35249593 -0.40139803 128.14877035] 2.2851110001414745 2.9981305351356577
str_9g = img[:,:,1].copy()
str_9g = str_9g[277:284,1071:1076]
plt.imshow(str_9g, cmap='Greens')
plt.show()
p0 = [10, 3, 6, 0.8, 1.4, 0, 0]
determinar(str_9g,p0,FWHM_xig,FWHM_yig, incer = str_9g.flatten())
[193.35932279 1.23687204 0.89029534 2.10755642 3.13488928 1.72376695 61.15734226] 2.9126110808108847 2.0964853040996143
str_10g = img[:,:,1].copy()
str_10g = str_10g[43:49,883:889]
plt.imshow(str_10g, cmap='Greens')
plt.show()
p0 = [10, 3.3, 3, 0.9, 1.4, 0, 0]
determinar(str_10g,p0,FWHM_xig,FWHM_yig, incer = str_10g.flatten())
[223.03027334 0.97687939 1.06282319 2.28196065 2.86747905 -15.75893458 61.71980515] 2.300375167201353 2.5027573548739817
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xig = abs(np.array(FWHM_xig))
print(FWHM_xig)
stdxig = FWHM_xig.std()
print(stdxr)
meanxig = FWHM_xig.mean()
print(meanxig)
medianxig = np.median(FWHM_xig)
print(medianxig)
plt.hist(FWHM_xig.flatten())
plt.title('Valores FHWM en la dirección x para la imagen verde considerando incertidumbre')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.06044976 2.35963715 5.94635051 1.72778985 2.5844445 1.98904453 4.28449873 2.285111 2.91261108 2.30037517] 1.1372820062135576 2.845031227650905 2.330006159119047
FWHM_yig = abs(np.array(FWHM_yig))
print(FWHM_yig)
stdyig = FWHM_yig.std()
print(stdyig)
meanyig = FWHM_yig.mean()
print(meanyig)
medianyig = np.median(FWHM_yig)
print(medianyig)
plt.hist(FWHM_yig.flatten())
plt.title('Valores FHWM en la dirección y para la imagen verde considerando incertidumbre')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.26459585 2.48867213 5.80659834 2.30828238 2.99444242 2.4728795 4.85154451 2.99813054 2.0964853 2.50275735] 1.1775361127169395 3.0784388321736786 2.495714741576891
FWHM_xib = []
FWHM_yib = []
str_1b = img[:,:,2].copy()
str_1b = str_1b[283:288,890:895]
plt.imshow(str_1b, cmap ='Blues')
plt.show()
p0 = [71, 1.7, 11, 1, 1, 0, 0]
determinar(str_1b,p0,FWHM_xib,FWHM_yib, incer = str_1b.flatten())
[196.08431561 0.86622721 0.95511986 2.01460987 2.04323191 -0.27430574 54.08557709] 2.039809207469812 2.249135380109503
str_2b = img[:,:,2].copy()
str_2b = str_2b[252:257,949:955]
plt.imshow(str_2b, cmap='Blues')
plt.show()
p0 = [50, 2.9, 2, 4, 1, 0, 0]
determinar(str_2b,p0,FWHM_xib,FWHM_yib, incer = str_2b.flatten())
[207.39523454 1.05340756 1.02290869 2.3785664 2.2034983 1.28470635 58.97679989] 2.480585247068367 2.4087658908107965
str_3b = img[:,:,2].copy()
str_3b = str_3b[309:324,619:636]
plt.imshow(str_3b, cmap='Blues')
plt.show()
p0 = [14, 15, 10, 4, 1, 0, 0]
determinar(str_3b,p0,FWHM_xib,FWHM_yib, incer = str_3b.flatten())
[189.3566925 2.75095638 2.67312878 8.30245528 6.75026981 -4.8848669 80.55741834] 6.47800722275355 6.294737245460569
str_4b = img[:,:,2].copy()
str_4b = str_4b[18:23,607:612]
plt.imshow(str_4b, cmap='Blues')
plt.show()
p0 = [9, 3, 6, 1.23, 1.46, 0, 0]
determinar(str_4b,p0,FWHM_xib,FWHM_yib, incer = str_4b.flatten())
[183.38176381 0.73671432 0.99916499 1.94780626 2.22736734 3.06024167 63.6038059 ] 1.7348296501661935 2.3528537457692438
str_5b = img[:,:,2].copy()
str_5b = str_5b[347:356,542:552]
plt.imshow(str_5b, cmap='Blues')
plt.show()
p0 = [8, 1.7, 1, 3.5, 1, 0, 0]
determinar(str_5b,p0,FWHM_xib,FWHM_yib, incer = str_5b.flatten())
[166.04481043 -1.21990199 1.38650382 4.52448629 4.2147847 -2.61241384 104.34156413] -2.8726496561210997 3.2649669805187034
str_6b = img[:,:,2].copy()
str_6b = str_6b[295:300,307:312]
plt.imshow(str_6b, cmap='Blues')
plt.show()
p0 = [123, 3, 5.9, 0.78, 1.7, 0, 0]
determinar(str_6b,p0,FWHM_xib,FWHM_yib, incer = str_6b.flatten())
[159.01793668 0.8552177 1.05546949 2.07193554 1.69511554 -0.21091069 89.9142516 ] 2.0138837865360038 2.4854407042295907
str_7b = img[:,:,2].copy()
str_7b = str_7b[142:154,387:399]
plt.imshow(str_7b, cmap='Blues')
plt.show()
p0 = [100, 2, 5, 3, 5, 0, 0]
determinar(str_7b,p0,FWHM_xib,FWHM_yib, incer = str_7b.flatten())
[137.6552333 1.92930042 1.73162599 6.12640336 5.62271455 -0.84345339 86.564405 ] 4.543155312707578 4.077667599505995
str_8b = img[:,:,2].copy()
str_8b = str_8b[459:465,322:327]
plt.imshow(str_8b, cmap='Blues')
plt.show()
p0 = [30, 3.43, 6.4, 0.01, 0.5, 0, 0]
determinar(str_8b,p0,FWHM_xib,FWHM_yib, incer = str_8b.flatten())
[152.36430635 1.12628322 1.39881391 2.06937349 2.37438637 -0.44230575 105.37899408] 2.6521942943381274 3.293955036969304
str_9b = img[:,:,2].copy()
str_9b = str_9b[277:284,1071:1076]
plt.imshow(str_9b, cmap='Blues')
plt.show()
p0 = [10, 3, 6.4, 7.3, 1.1, 0, 0]
determinar(str_9b,p0,FWHM_xib,FWHM_yib, incer = str_9b.flatten())
[193.86710559 0.88686474 1.23340715 2.10631332 3.10028245 -15.54702541 41.51349323] 2.088406856621961 2.9044518727402315
str_10b = img[:,:,2].copy()
str_10b = str_10b[43:49,883:889]
plt.imshow(str_10b, cmap='Blues')
plt.show()
p0 = [41, 12, 4, 1, 1, 0, 0]
determinar(str_10b,p0,FWHM_xib,FWHM_yib, incer = str_10b.flatten())
[ 2.21501512e+02 9.47560514e-01 1.06364334e+00 2.29364326e+00 2.86413861e+00 -1.70292345e-02 4.60375723e+01] 2.231334493416455 2.5046886503061434
Considere el valor en dirección x & en dirección y. Saque la desviación estándar, el promedio, la mediana y el histograma
FWHM_xib = abs(np.array(FWHM_xib))
print(FWHM_xib)
stdxib = FWHM_xib.std()
print(stdxib)
meanxib = FWHM_xib.mean()
print(meanxib)
medianxib = np.median(FWHM_xib)
print(medianxib)
plt.hist(FWHM_xib.flatten())
plt.title('Valores FHWM en la dirección x para la image azul')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.03980921 2.48058525 6.47800722 1.73482965 2.87264966 2.01388379 4.54315531 2.65219429 2.08840686 2.23133449] 1.4045244740696985 2.9134855727199147 2.355959870242411
FWHM_yib = abs(np.array(FWHM_yib))
print(FWHM_yib)
stdyib = FWHM_yib.std()
print(stdyib)
meanyib = FWHM_yib.mean()
print(meanyib)
medianyib = np.median(FWHM_yib)
print(medianyib)
plt.hist(FWHM_yib.flatten())
plt.title('Valores FHWM en la dirección y para la image azul')
plt.xlabel('Valor FHWM')
plt.ylabel('No de imagénes con ese valor')
plt.show()
[2.24913538 2.40876589 6.29473725 2.35285375 3.26496698 2.4854407 4.0776676 3.29395504 2.90445187 2.50468865] 1.1699287444022204 3.183666310642008 2.7045702615231875
En dirección x
print(stdxi,stdxir,stdxig,stdxib)
print(meanxi,meanxir,meanxig,meanxib)
print(medianxi,medianxir,medianxig,medianxib)
1.2442375593713781 1.1220401417419346 1.233085493598224 1.4045244740696985 2.916621675050629 2.9834775770762554 2.845031227650905 2.9134855727199147 2.519035453328353 2.536875478494914 2.330006159119047 2.355959870242411
En dirección y
print(stdyi,stdyir,stdyig,stdyib)
print(meanyi,meanyir,meanyig,meanyib)
print(medianyi,medianyir,medianyig,medianyib)
1.2080118897663938 1.1882827213448535 1.1775361127169395 1.1699287444022204 3.0706026542016884 2.9267160396079794 3.0784388321736786 3.183666310642008 2.4351553795125156 2.347568909171641 2.495714741576891 2.7045702615231875
Para x
# -- crear un array con todas las cantidades calculadas anteriormente para el FWHM
valuesx = np.array([[stdx,stdxr,stdxg,stdxb, meanx,meanxr,meanxg,meanxb, medianx,medianxr,medianxg,medianxb],
[stdxi,stdxir,stdxig,stdxib, meanxi,meanxir,meanxig,meanxib, medianxi,medianxir,medianxig,medianxib]])
valuesx.shape
(2, 12)
# -- crear un data frame con esa información
FWHMxx = pd.DataFrame(data=valuesx,
columns=['stdxL','stdxR','stdxG','stdxB','meanxL','meanxR','meanxG','meanxB','medianxL','medianxR','medianxG','medianxB'],
index=['Sin-incertidumbre','Con-incertidumbre'])
# -- crear un data frame con la diferencia entre valores con y sin incertidumbre
resta = pd.DataFrame(data=abs(FWHMxx.loc['Sin-incertidumbre'] - FWHMxx.loc['Con-incertidumbre']).values.reshape(1,12),
columns=['stdxL','stdxR','stdxG','stdxB','meanxL','meanxR','meanxG','meanxB','medianxL','medianxR','medianxG','medianxB'],
index=['diferencia'])
# -- combinar ambos data frame
FWHMxx = FWHMxx.append(resta)
Para y
# -- lo mismo pero en dirección y
valuesy = np.array([[stdy,stdyr,stdyg,stdyb, meany,meanyr,meanyg,meanyb, mediany,medianyr,medianyg,medianyb],
[stdyi,stdyir,stdyig,stdyib, meanyi,meanyir,meanyig,meanyib, medianyi,medianyir,medianyig,medianyib]])
FWHMyy = pd.DataFrame(data=valuesy,
columns=['stdyL','stdyR','stdyG','stdyB','meanyL','meanyR','meanyG','meanyB','medianyL','medianyR','medianyG','medianyB'],
index=['Sin-incertidumbre','Con-incertidumbre'])
resta = pd.DataFrame(data=abs(FWHMyy.loc['Sin-incertidumbre'] - FWHMyy.loc['Con-incertidumbre']).values.reshape(1,12),
columns=['stdyL','stdyR','stdyG','stdyB','meanyL','meanyR','meanyG','meanyB','medianyL','medianyR','medianyG','medianyB'],
index=['diferencia'])
FWHMyy = FWHMyy.append(resta)
FWHMxx
stdxL | stdxR | stdxG | stdxB | meanxL | meanxR | meanxG | meanxB | medianxL | medianxR | medianxG | medianxB | |
---|---|---|---|---|---|---|---|---|---|---|---|---|
Sin-incertidumbre | 1.202031 | 1.137282 | 1.192085 | 1.301895 | 2.938495 | 2.834198 | 2.824878 | 2.982611 | 2.535912 | 2.463149 | 2.424866 | 2.517526 |
Con-incertidumbre | 1.244238 | 1.122040 | 1.233085 | 1.404524 | 2.916622 | 2.983478 | 2.845031 | 2.913486 | 2.519035 | 2.536875 | 2.330006 | 2.355960 |
diferencia | 0.042207 | 0.015242 | 0.041001 | 0.102629 | 0.021874 | 0.149280 | 0.020153 | 0.069125 | 0.016877 | 0.073726 | 0.094860 | 0.161566 |
FWHMyy
stdyL | stdyR | stdyG | stdyB | meanyL | meanyR | meanyG | meanyB | medianyL | medianyR | medianyG | medianyB | |
---|---|---|---|---|---|---|---|---|---|---|---|---|
Sin-incertidumbre | 1.163140 | 1.076054 | 1.124040 | 1.219566 | 3.058655 | 3.106182 | 3.099311 | 3.108883 | 2.535681 | 2.816984 | 2.788384 | 2.607229 |
Con-incertidumbre | 1.208012 | 1.188283 | 1.177536 | 1.169929 | 3.070603 | 2.926716 | 3.078439 | 3.183666 | 2.435155 | 2.347569 | 2.495715 | 2.704570 |
diferencia | 0.044872 | 0.112229 | 0.053496 | 0.049637 | 0.011948 | 0.179466 | 0.020872 | 0.074783 | 0.100526 | 0.469415 | 0.292670 | 0.097341 |
Al incluir la incertidumbre los valores de FWHM si cambian. No se puede afirmar que el valor sea mas o menos preciso, sin embargo al darle un peso a cada pixel que conforma la imagen de la estrella, es mas facil manejar imagenes con mayor ruido como la Estrella 3 y 7.
El módelo es extremadamente sensible a cualquier cambio como el cambio de los parametros iniciales p0, cambio en los pixeles alrededor de la estrella.