sigmanc=input('Inserisci sigmanc')
#Devo simulare q volte 450 misure di pendolo
outfile_9=open('stimesim.dat', 'w')
outfile_9.write('numero media varianzanc\n')
outfile_10=open('datisim.dat', 'w')
q=10000
m=450
#CICLO DI SIMULAZIONE
for i in range(0,q):
 a=np.random.normal(muT,sigmanc,m)
 s=a
 
 outfile_9.write(('%f %f\n')%(i,j))
 musim=sum(s)/m
 #print i
 
 differenzesim=[(j-musim)**2 for j in s ]
 varncsim=sum(differenzesim)/m
 #print musim, varncsim
 outfile_9.write(('%f %f %f\n')%(i,musim,varncsim))
 
 
outfile_9.close() 
ni,musim,varncsim=np.loadtxt('stimesim.dat',skiprows=1,usecols=(0,1,2),unpack=True)
#print len(mu)
#print len(varnc)
mediamedia=sum(musim)/len(musim) #media delle medie
mediavar=sum(varncsim)/len(varncsim) #Media delle varianze

differenzesim2=[(k-mediamedia)**2 for k in musim]
varianzamedia=(sum(differenzesim2))/(len(musim)*(len(musim)-1)) #Varianza della media
sigmamedia=math.sqrt(varianzamedia)
print 'La media delle medie e', mediamedia
print 'La media delle varianze e', mediavar
print 'La varianza della media e', varianzamedia
print 'La deviazione standard della media e', sigmamedia
print 'Test di Pearson per i dati della simulazione',

#istogrammi
num_binsmu=70
ns, binssim, patches=plt.hist(musim,num_binsmu,normed=0) #bins=Estremo destro
plt.xlabel('Valore')
plt.ylabel('Frequenza')
plt.title("Medie")
plt.show()

num_binsvar=65
ns2, binssim2, patches=plt.hist(varncsim,num_binsvar,normed=0) #bins=Estremo destro
plt.xlabel('Valore')
plt.ylabel('Frequenza')
plt.title("Varianze")
plt.show()