import numpy as np
from scipy.optimize import curve_fit

x = np.arange(100, dtype=np.float)
y = np.zeros_like(x)
yn = y + np.random.normal(size=len(y))


def const(x, a):
    return a * np.ones_like(x)

popt, pcov = curve_fit(const, x, yn)
print 'Best fit constant =', popt[0]
print 'Uncertainty on constant =', np.sqrt(pcov[0, 0])
plot(x, yn, '.')
plot(x, const(x, popt[0]))
grid()

sigma = 100
popt, pcov = curve_fit(const, x, yn, sigma=sigma)
print 'Best fit constant =', popt[0]
print 'Uncertainty on constant =', np.sqrt(pcov[0, 0])

chi = (yn - const(x, *popt)) / sigma
chi2 = (chi ** 2).sum()
dof = len(x) - len(popt)
factor = (chi2 / dof)
pcov_sigma = pcov / factor
print chi2
print dof
print factor
print np.sqrt(pcov_sigma[0, 0])

import sherpa.ui as ui  # sherpa user interface

ui.load_arrays(1, x, yn, 100.0 * np.ones_like(yn))

ui.set_model(1, ui.const1d.c1)
c1.c0.min = -100  # set the parameter bounds to -100 to +100
c1.c0.max = 100

ui.fit()

ui.plot_fit()  # this shows the large error bars

ui.covar()