import random
from matplotlib import pyplot as plt
import numpy as np
from scipy.special import binom


def binomial_sample(n, p):
	"""
	Take a sample of size n from a binomial distribution with
	success probability equal p.
	Returns the number of sucesses.
	"""
	s = 0
	for i in xrange(n):
		x = random.random()
		if x < p:
			s += 1
	return s


def binomial(x, n, p):
	return binom(n, x) * p**x * (1-p)**(n-x)


def exp1():
	n = 15
	p = 0.75
	s = binomial_sample(n, p)
	hist_data = [binomial_sample(n, p) for x in xrange(1000)]
	theory_x = range(0, n + 1)
	theory_y = [binomial(x, n, p) for x in theory_x]
	plt.hist(hist_data, normed=1)
	plt.plot(theory_x, theory_y, '--')
	plt.xlabel('Fotones detectados')
	plt.ylabel('Tasa')
	plt.savefig('figb.jpg')
	plt.show()


def main():
	n = 10
	p = 0.7
	s = binomial_sample(n, p)
	print '-----Ejercicio 15-------'
	print '-----a)-------'
	print 'n = %s, p = %s, s= %s' % (n, p, s)
	print '-----b)-------'
	exp1()



if __name__ == '__main__':
	main()
@@ -95,16 +95,30 @@
 t binom%0A
+import inspect
 %0A%0Adef bi
@@ -644,16 +644,41 @@
 normed=1
+, label=r'$Simulaci%5C'on$'
 )%0A%09plt.p
@@ -702,18 +702,56 @@
 ory_y, '
---
+k-', label=r'$Distribuci%5C'on%5C Te%5C'orica$
 ')%0A%09plt.
@@ -799,16 +799,46 @@
 'Tasa')%0A
+%09plt.legend(loc='upper left')%0A
 %09plt.sav
@@ -950,32 +950,32 @@
 icio 15-------'%0A
-
 %09print '-----a)-
@@ -974,32 +974,74 @@
 -----a)-------'%0A
+%09print inspect.getsource(binomial_sample)%0A
 %09print 'n = %25s, 

add print function definition. add legend to plot