import numpy
import random

from matplotlib.pyplot import *

class Binomiale:
    def __init__(self,p):
        self.p = p
    def tirage(self):
        x = random.random()
        if x<self.p:
            return 1
        else:
            return 0
    def moyenne(self,N):
        somme = 0
        somme2 = 0
        for k in range(N):
            t = self.tirage()
            somme += t
            somme2 += t*t
        moyenne = somme*1.0/N
        variance = somme2*1.0/N-moyenne*moyenne
        ecart = numpy.sqrt(variance*1.0/N)
        return (moyenne,variance,ecart)
    def distribution_moyennes(self,N,N_tirages,M):
        hist = numpy.zeros(M)
        for t in range(N_tirages):
            (m,v,e) = self.moyenne(N)
            i = int(m*M)
            hist[i] += 1
        return (numpy.arange(M)*1.0/M,hist/N_tirages*M)
    def binomiale(self,N,N_tirages):
        hist = numpy.zeros(N+1)
        for t in range(N_tirages):
            s = 0
            for k in range(N):
                s += self.tirage()
            hist[int(s)] += 1
        return hist/N_tirages
            

p=0.5
binom = Binomiale(p)
N=1000
(m,v,e) = binom.moyenne(N)
            

M=N
Nt=100000
(moy,hist) = binom.distribution_moyennes(N,Nt,M)
x=numpy.linspace(0,1,1000)
sigma=numpy.sqrt(p*(1-p)/N)
normale=1.0/(numpy.sqrt(2*numpy.pi)*sigma)*numpy.exp(-(x-p)**2/(2*sigma**2))
figure()
plot(moy,hist,"ko")
plot(x,normale,"r")
xlabel("m")
ylabel("P")
grid()
             

p=0.5
binom = Binomiale(p)
hist = binom.binomiale(10,100000)
figure()
plot(hist,"o-")
xlabel("n")
ylabel("P")
grid()
title("Binomiale")
             

p=0.3
binom = Binomiale(p)
hist = binom.binomiale(10,100000)
figure()
plot(hist,"o-")
xlabel("n")
ylabel("P")
grid()
title("Binomiale")