import numpy as np
import matplotlib.pyplot as plt

def confint95(data):
    """returns 95% confidence interval for data"""
    n = len(data)
    mean = np.mean(data)
    standard_error = np.std(fid)/np.sqrt(n)
    return mean-1.96*standard_error, mean+1.96*standard_error

def histogram(data, numbins):
    """returns centred bins"""
    # range = [np.min(data), np.max(data)]
    range = [0.977,1]
    n, bin_edges = np.histogram(data, bins=numbins, range=range)
    bin_size = bin_edges[1] - bin_edges[0]
    bin_centers = bin_edges[:-1] + bin_size/2.
    return n, bin_centers

fid = np.loadtxt('fidelities_1s_measurement.dat')
fid = fid**2 #forgot to squre...to
np.savetxt('corrected_fidelities_1s_measurement.dat', fid)

#CONFIDENCE INTERVAL#
print(confint95(fid))

#PLOT#
plt.figure()
n, b = histogram(fid, 10)
plt.errorbar(b, n, np.sqrt(n))
np.savetxt('corrected_fidelities_1s_measurement_histogram.dat', np.vstack([b, n, np.sqrt(n)]).T)
plt.savefig('corrected_fidelities_1s_measurement.pdf')
plt.show()

# print(np.array([np.percentile(fid, q) for q in [5,50,95]]))
# plt.figure()
# plt.hist(fid);
# plt.xlim(0.97,1)
# plt.axvline(np.mean(fid), color='black')
# plt.axvline(np.mean(fid)-np.std(fid), color='black')
# plt.axvline(np.mean(fid)+np.std(fid), color='black')
# plt.savefig('fidelities.pdf')
# plt.show()