In [18]:
from IPython.display import display
from IPython.display import HTML
import IPython.core.display as di # Example: di.display_html('<h3>%s:</h3>' % str, raw=True)
# This line will hide code by default when the notebook is exported as HTML
di.display_html('<script>jQuery(function() {if (jQuery("body.notebook_app").length == 0) { jQuery(".input_area").toggle(); jQuery(".prompt").toggle();}});</script>', raw=True)
# This line will add a button to toggle visibility of code blocks, for use with the HTML export version
di.display_html('''<button onclick="jQuery('.input_area').toggle(); jQuery('.prompt').toggle();">Pulse para ver el codigo</button>''', raw=True)
In [19]:
#!/usr/bin/env python
# -*- coding: utf-8 -*-
# Librerías
import matplotlib
from scipy import misc
from scipy import stats
from scipy import special
import pylab as plt
import numpy as np
import math
%matplotlib inline
font = {'weight' : 'bold',
'size' : 16}
matplotlib.rc('font', **font)
figura 2 del artículo: la tercera columna de la figura corresponde a la probabilidad de detectar un número n de fotones. Suponiendo que la distribución es de Poisson para cada gráfica, extraiga los puntos coordenados, determine el valor promedio, y grafique en una misma figura los puntos coordenados y la correspondiente curva de Poisson. Reporte el valor promedio de fotones que encontró.
para realizar este punto se usan los datos de S. Haroche, et. al. del articulo: "Quantum Rabi Oscillation: A Direct Test of Field Quantization in a Cavity"
se usa el proceso aprendido en clase para graficar los datos
In [56]:
##################################
## grafica 1
##################################
#graficas de probabilidades del paper:
hAlfa=np.array([0.95,0.05,0.0,0.0,0.0,0.0])
hBeta=np.array([0.65,0.3,0.05,0.0,0.0,0.0])
hGamma=np.array([0.45,0.35,0.15,0.08,0.05,0.0])
hDelta=np.array([0.2,0.3,0.25,0.15,0.08,0.05])
n=[0,1,2,3,4,5] #X axis
h=np.array([hAlfa,hBeta,hGamma,hDelta]) #todos los datos unidos
#print ( sum(h[:,0])/len(n))
promedios=np.array([sum(h[:,0])/len(n),sum(h[:,1])/4.,sum(h[:,2])/4.,sum(h[:,3])/4.,sum(h[:,4])/4.,sum(h[:,5])/4.])
lamb=sum(promedios)/len(promedios)
plt.figure(figsize=(12,7))
plt.plot(n,promedios,'*',ms=10, label='distribucion con media = %.1f' % lamb)
plt.legend()
plt.xlabel('n')
plt.ylabel('probabilidad media')
plt.xlim(0,5)
plt.grid()
plt.show()
In [58]:
num = 0.1 # probabilidad esperada
prob1 = (lamb**num*np.exp(-lamb)/special.gamma(num))*100 # Probabilidad de encontrar dicho número
x = np.arange(0,6) # rango de datos
histP1 = stats.poisson.pmf(x, lamb) # función de probabilidad de Poisson
ProbP1 = (np.sum(histP1[0:int(num)+1]))*100 # Probabilidad acumulada
print('Promedio de probabilidades de obtener n fotones = %.1f' % lamb)
print('La probabilidad de que se obtenga un valor de %.1f es = %.1f%% ' % (num,prob1))
print('Probabilidad de obtener hasta %.1f = %.1f%%' %(num,ProbP1))
In [59]:
plt.figure(figsize=(16,9))
plt.plot(x, histP1, 'go--', ms=8, label='Distribucion de Poisson con $\lambda=$ %.1f' % lamb)
plt.plot(n,promedios,'k*',ms=10,label='Conjunto de datos con promedio %.1f' % lamb)
plt.xlabel('n')
plt.ylabel('probabilidad media')
plt.legend()
plt.grid()
ya que los datos muestran los rangos de probabilidad de obtener n fotones emitidos en el experimento mencionado en el paper. se puede observar como estas probabilidades se comportan como una Poissoniana con un $\lambda = 0.1$ esto ya que la distribución aparenta no tener una simetría además de que el comportamiento de los datos y la distribución teórica se comportan practicamente igual. El pequeño margen de error entre los puntos y la teorica se puede deber a la forma de obtener los datos, la cual fue simplemente de observar las graficas de la figura 2 del paper y aproximar los valores con el eje coordenado de probabilidades P(n)
para crear la imagen de intensidad por cuestiones de tiempo, solo se recurrió a la librería PIL con la cual se invirtió el color de la imagen y luego se cargó nuevamente al programa en falso color en la escala de grises. finalmente se procede igual que en clase para hallar las distribuciones
In [2]:
from PIL import Image
import PIL.ImageOps
image = Image.open('stars.jpg')
inverted_image = PIL.ImageOps.invert(image)
inverted_image.save('stars_inv.png')
Ima = misc.imread('stars_inv.png') # Se lee la imagen como matriz en escala de 8 bit
plt.rcParams['figure.figsize'] = 20, 6 # para modificar el tamaño de la figura
Imab = Ima[100:500,100:700,1] # La imagen original tenía tres canales (RGB); se elige un canal y se recorta
plt.figure(2)
plt.imshow(Imab, cmap='gray')
Out[2]:
In [3]:
plt.rcParams['figure.figsize'] = 18, 15 # para modificar el tamaño de la figura
fil, col = Imab.shape # número de filas y columnas de la imagen
numlado = 7# Número de imágenes por lado
contar = 1
plt.figure(5)
for enfil in range(1,numlado+1):
for encol in range(1,numlado+1):
plt.subplot(numlado,numlado,contar)
plt.imshow(Imab[(enfil-1)*np.int(fil/numlado):enfil*np.int(fil/numlado), \
(encol-1)*np.int(col/numlado):encol*np.int(col/numlado)],cmap='gray')
frame1 = plt.gca()
frame1.axes.get_yaxis().set_visible(False)
frame1.axes.get_xaxis().set_visible(False)
contar = contar + 1
In [54]:
# Para el caso de 10x10 imágenes en gal se presentan el número de estrellas contadas
est2 = np.array([14., 13., 9., 12., 13., 15., 17., 14., 14., 10., \
7., 13., 15., 16., 10., 11., 15., 13., 13., 12., \
16., 15., 13., 13., 14., 13., 12., 10., 9., 12., \
11., 14., 14., 12., 12., 13., 10., 13., 12., 12., \
17., 9., 15., 13., 15., 12., 14., 10., 9., 16., \
8., 9., 12., 11., 8., 11., 11., 13., 13., 12., \
18., 11., 11., 7., 12., 14., 17., 10., 12., 9., \
15., 7., 18., 9., 7., 13., 9., 10., 17., 13., \
16., 13., 10., 12., 7., 7., 11., 15., 14., 12., \
11., 13., 11., 11., 10., 11., 11., 8., 11., 14.])
la2 = np.mean(est) # Valor promedio del conjunto de datos
# Distribución del conjunto de datos. La primera fila es el número de estrellas, la segunda es el número de veces que
# se repite dicho número de estrellas
distriEst2 = np.array([[0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28., 29., 30.],
[0., 0., 0., 0., 0., 0., 0., 6., 3., 9., 9., 14., 15., 17., 10., 8., 4., 4., 2., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])
print('Valor promedio del conjunto de datos = %.2f' % la2)
plt.figure(figsize=(16,9))
plt.plot(distriEst2[0,:],distriEst2[1,:]/est2.size,'r*',ms=10,label='Distribucion datos con promedio %.1f' % la2)
plt.legend()
plt.xlabel('Numero de estrellas en el intervalo')
plt.ylabel('Rata de ocurrencia')
plt.grid()
In [47]:
num = 12. # Número de estrellas que se espera encontrar
prob2 = (la2**num*np.exp(-la2)/math.factorial(num))*100 # Probabilidad de encontrar dicho número de estrellas
x = np.arange(0,30) # rango de datos: número de estrellas
histP2 = stats.poisson.pmf(x, la2) # función de probabilidad de Poisson
ProbP2 = (np.sum(histP2[0:int(num)+1]))*100 # Probabilidad acumulada
print('Promedio de estrellas en el área estudiada = %.2f , caso 2' % la2)
print('La probabilidad de que se observe en la imagen del espacio profundo %d estrellas es = %.1f%% caso2' % (num,prob2))
print('Probabilidad de observar hasta %d estrellas = %.1f%%, caso 2' %(num,ProbP2))
In [48]:
plt.figure(figsize=(16,9))
plt.plot(x, histP2, 'go--', ms=8, label='Distribucion de Poisson con $\lambda=$ %.1f c2' % la2)
plt.plot(distriEst2[0,:],distriEst2[1,:]/est2.size,'k*',ms=10,label='Conjunto de datos con promedio %.1f c2' % la2)
plt.xlabel('Numero de estrellas (sucesos)')
plt.ylabel('Rata de ocurrencia')
plt.legend()
plt.grid()
en este caso se puede observar como los datos siguen aproximadamente la distribución de poisson, pero se observa que los datos llevan a un pico mayor en el centro de la distribución, lo cual podría mostrar un comportamiento más gaussiano que poissoniano, ya que la primera para $\sigma$s pequeños suele tener una campana más angosta, además de que la segunda mientras mayor sea $\lambda$ más se aproximará esta a una distribussión normal o gaussiana.