import numpy import matplotlib.pyplot as plt r1, r2, r3, r4, r5 = numpy.load('random.npy') N = r1.size # PDF Histo print "PDF histogram" plt.figure() plt.subplot(221) plt.hist(r1, 100, histtype='step', normed=True) plt.hist(r2, 100, histtype='step', normed=True) plt.hist(r3, 100, histtype='step', normed=True) plt.hist(r4, 100, histtype='step', normed=True) plt.hist(r5, 100, histtype='step', normed=True) # Chi square print "chi square" def chisquare(rs, sample_size=10): N = rs.size n = N/sample_size means = numpy.empty(n) for i in range(n): means[i] = rs[i*sample_size:(i+1)*sample_size].mean() return means plt.subplot(222) def gaussian(x, mu, sigma): return 1./(sigma * numpy.sqrt(2.0*numpy.pi)) * numpy.exp(-(x-mu)**2 / (2.0*sigma**2)) x = numpy.linspace(0.0, 1.0, 100) plt.plot(x, gaussian(x, 0.5, 1.0/numpy.sqrt(120.0)), linewidth=2, color='r') plt.hist(chisquare(r1), 100, histtype='step', normed=True) plt.hist(chisquare(r2), 100, histtype='step', normed=True) plt.hist(chisquare(r3), 100, histtype='step', normed=True) plt.hist(chisquare(r4), 100, histtype='step', normed=True) plt.hist(chisquare(r5), 100, histtype='step', normed=True) plt.legend(("expected", "r1", "r2", "r3", "r4", "r5")) # Fourier analysis print "Fourier analysis" def fourier(rs): return numpy.abs(numpy.fft.rfft(rs))/rs.size plt.subplot(223) plt.semilogy(fourier(r1)) plt.semilogy(fourier(r2)) plt.semilogy(fourier(r3)) plt.semilogy(fourier(r4)) plt.semilogy(fourier(r5)) # Triplet correlation print "triplet correlation" def triplet_corr_error(rs, p=250): corr = numpy.empty(p-1) N = rs.size for k in range(1,p): c = rs[p:]*rs[:-p]*rs[p-k:-k] corr[k-1] = c.sum() / (N-p) error = abs(corr - 0.5**3) return error plt.subplot(224) plt.semilogy(triplet_corr_error(r1), '-') plt.semilogy(triplet_corr_error(r2), '-') plt.semilogy(triplet_corr_error(r3), '-') plt.semilogy(triplet_corr_error(r4), '-') plt.semilogy(triplet_corr_error(r5), '-') plt.show()