forked from scikit-hep/probfit
-
Notifications
You must be signed in to change notification settings - Fork 0
/
testfit.py
95 lines (79 loc) · 3.1 KB
/
testfit.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
import unittest
from probfit.costfunc import UnbinnedLH, BinnedLH, BinnedChi2, Chi2Regression,\
SimultaneousFit
from probfit.pdf import gaussian, linear
from numpy.random import randn, seed
from math import log,pi,sqrt
import warnings
from iminuit.iminuit_warnings import InitialParamWarning
from probfit.util import describe
from nose.tools import *
from probfit.funcutil import rename
import numpy as np
from iminuit import Minuit
def assert_almost_equal(x,y,delta=1e-7):
if y-delta < x < y+delta:
pass
else:
raise AssertionError('x = %f and y = %f differs more than %g'%(x,y,delta))
class TestFit(unittest.TestCase):
def setUp(self):
warnings.simplefilter("ignore", InitialParamWarning)
seed(0)
self.ndata = 20000
self.data = randn(self.ndata)
self.analytic = self.ndata*0.5*(log(2*pi)+1)
def test_UnbinnedLH(self):
f = gaussian
assert_equal(list(describe(f)), ['x','mean','sigma'])
lh = UnbinnedLH(gaussian, self.data,)
assert_equal(list(describe(lh)), ['mean','sigma'])
assert_almost_equal(lh(0,1), 28188.201229348757)
m = Minuit(lh)
assert_equal(m.errordef, 0.5)
def test_BinnedLH(self):
#write a better test... this depends on subtraction
f = gaussian
assert_equal(list(describe(f)), ['x','mean','sigma'])
lh = BinnedLH(gaussian, self.data, bound=[-3,3])
assert_equal(list(describe(lh)), ['mean','sigma'])
assert_almost_equal(lh(0,1), 20.446130781601543, 1)
m = Minuit(lh)
assert_equal(m.errordef, 0.5)
def test_BinnedChi2(self):
f = gaussian
assert_equal(list(describe(f)), ['x','mean','sigma'])
lh = BinnedChi2(gaussian, self.data, bound=[-3,3])
assert_equal(list(describe(lh)), ['mean','sigma'])
assert_almost_equal(lh(0,1), 19951.005399882044, 1)
m = Minuit(lh)
assert_equal(m.errordef, 1.0)
def test_Chi2Regression(self):
x = np.linspace(1, 10, 10)
y = 10*x+1
f = linear
assert_equal(list(describe(f)), ['x','m','c'])
lh = Chi2Regression(f, x, y)
assert_equal(list(describe(lh)), ['m','c'])
assert_almost_equal(lh(10, 1), 0)
assert_almost_equal(lh(10, 0), 10.)
m = Minuit(lh)
assert_equal(m.errordef, 1.0)
def test_simultaneous(self):
seed(0)
data = randn(10000)
shifted = data+3.
g1 = rename(gaussian, ['x', 'lmu', 'sigma'])
g2 = rename(gaussian, ['x', 'rmu', 'sigma'])
ulh1 = UnbinnedLH(g1,data)
ulh2 = UnbinnedLH(g2,shifted)
sim = SimultaneousFit(ulh1, ulh2)
assert_equal(describe(sim),['lmu', 'sigma', 'rmu'])
m = Minuit(sim,sigma=1.2, pedantic=False, print_level=0)
m.migrad()
assert(m.migrad_ok())
assert_almost_equal(m.values['lmu'], 0., delta=2*m.errors['lmu'])
assert_almost_equal(m.values['rmu'], 3., delta=2*m.errors['rmu'])
assert_almost_equal(m.values['sigma'], 1., delta=2*m.errors['sigma'])
if __name__ == '__main__':
unittest.main()