IMinuit

Minuit - программа численной минимизации функций многих переменных, широко применяемая в физике элементарных частиц. Есть два питонских интерфейса, PyMinuit и IMinuit (он особенно удобен в ipython).

In [1]:
from iminuit import Minuit
%pylab inline
Populating the interactive namespace from numpy and matplotlib

Простой пример

Определим квадратичную функцию от двух параметров.

In [2]:
def f(a,b):
    return 10*a**2+10*b**2-16*a*b+12*a-24*b

Создадим объект класса Minuit. a и b - грубые догадки, около чего надо искать минимум; error_a и error_b - оценки точности этих догадок (в начале минимизации программа будет делать шаги порядка этих величин, потом они будут уменьшаться). Пределы изменения задавать не обязательно. Валичина errordef показывает, насколько функция должна быть выше своего минимума, чтобы это считалось отклонением на одну сигму; поскольку минимизируемая функция - это, как правило, $\chi^2$, значение 1 по умолчанию вполне годится.

In [3]:
m=Minuit(f,a=0,error_a=1,limit_a=(-10,10),
         b=0,error_b=1,limit_b=(-10,10))
/usr/lib64/python3.4/site-packages/ipykernel/__main__.py:2: InitialParamWarning: errordef is not given. Default to 1.
  from ipykernel import kernelapp as app

Наиболее популярный метод минимизации - migrad.

In [4]:
m.migrad()

FCN = -17.999997281186666 TOTAL NCALL = 36 NCALLS = 36
EDM = 2.7187527367825896e-06 GOAL EDM = 1e-05 UP = 1.0
Valid Valid Param Accurate Covar PosDef Made PosDef
True True True True False
Hesse Fail HasCov Above EDM Reach calllim
False True False False
+ Name Value Parab Error Minos Error- Minos Error+ Limit- Limit+ FIXED
1 a 0.999823 0.526788 0 0 -10.0 10.0
2 b 1.99935 0.526776 0 0 -10.0 10.0

Out[4]:
({'fval': -17.999997281186666, 'up': 1.0, 'is_above_max_edm': False, 'has_valid_parameters': True, 'edm': 2.7187527367825896e-06, 'has_made_posdef_covar': False, 'has_covariance': True, 'has_posdef_covar': True, 'nfcn': 36, 'has_accurate_covar': True, 'hesse_failed': False, 'is_valid': True, 'has_reached_call_limit': False},
 [{'is_const': False, 'number': 0, 'upper_limit': 10.0, 'lower_limit': -10.0, 'name': 'a', 'is_fixed': False, 'error': 0.5267876810263843, 'value': 0.9998227792225478, 'has_upper_limit': True, 'has_limits': True, 'has_lower_limit': True},
  {'is_const': False, 'number': 1, 'upper_limit': 10.0, 'lower_limit': -10.0, 'name': 'b', 'is_fixed': False, 'error': 0.5267762744877214, 'value': 1.9993477581584944, 'has_upper_limit': True, 'has_limits': True, 'has_lower_limit': True}])

Значения параметров.

In [5]:
m.values
Out[5]:
{'a': 0.9998227792225478, 'b': 1.9993477581584944}

Значение функции в точке минимума.

In [6]:
m.fval
Out[6]:
-17.999997281186666

Ошибки параметров.

In [7]:
m.errors
Out[7]:
{'a': 0.5267876810263843, 'b': 0.5267762744877214}

Если, скажем, a - наш окончательный физический результат, то мы напишем в статье $a=1\pm0.5$. На самом деле у нас есть больше информации, поскольку ошибки a и b сильно скоррелированы. Матрица корреляции ошибок:

In [8]:
m.matrix()
Out[8]:
((0.27776493841711364, 0.22220727646667665),
 (0.22220727646667665, 0.27776101882046583))

Минимизация квадратичной формы сводится к решению системы линейных уравнений, а матрица корреляции ошибок - обратная матрица этой системы. В таком простом случае не имеет смысла использовать инструмент минимизации произвольных функций, такой, как Minuit.

In [9]:
M=array([[10.,-8.],[-8.,10.]])
M=inv(M)
M
Out[9]:
array([[ 0.27777778,  0.22222222],
       [ 0.22222222,  0.27777778]])
In [10]:
M.dot(array([[-6],[12]]))
Out[10]:
array([[ 1.],
       [ 2.]])

Нарисуем контуры, соответствующие отклонению на 1, 2 и 3 сигмы от оптимальной точки.

In [11]:
m.draw_mncontour('a','b',nsigma=3)
Out[11]:
(array([-1.03835166, -0.997266  , -0.95618035, -0.9150947 , -0.87400905,
        -0.8329234 , -0.79183775, -0.7507521 , -0.70966645, -0.6685808 ,
        -0.62749514, -0.58640949, -0.54532384, -0.50423819, -0.46315254,
        -0.42206689, -0.38098124, -0.33989559, -0.29880994, -0.25772428,
        -0.21663863, -0.17555298, -0.13446733, -0.09338168, -0.05229603,
        -0.01121038,  0.02987527,  0.07096092,  0.11204657,  0.15313223,
         0.19421788,  0.23530353,  0.27638918,  0.31747483,  0.35856048,
         0.39964613,  0.44073178,  0.48181743,  0.52290309,  0.56398874,
         0.60507439,  0.64616004,  0.68724569,  0.72833134,  0.76941699,
         0.81050264,  0.85158829,  0.89267395,  0.9337596 ,  0.97484525,
         1.0159309 ,  1.05701655,  1.0981022 ,  1.13918785,  1.1802735 ,
         1.22135915,  1.26244481,  1.30353046,  1.34461611,  1.38570176,
         1.42678741,  1.46787306,  1.50895871,  1.55004436,  1.59113001,
         1.63221567,  1.67330132,  1.71438697,  1.75547262,  1.79655827,
         1.83764392,  1.87872957,  1.91981522,  1.96090087,  2.00198652,
         2.04307218,  2.08415783,  2.12524348,  2.16632913,  2.20741478,
         2.24850043,  2.28958608,  2.33067173,  2.37175738,  2.41284304,
         2.45392869,  2.49501434,  2.53609999,  2.57718564,  2.61827129,
         2.65935694,  2.70044259,  2.74152824,  2.7826139 ,  2.82369955,
         2.8647852 ,  2.90587085,  2.9469565 ,  2.98804215,  3.0291278 ]),
 array([ -3.83516801e-02,   2.73397445e-03,   4.38196290e-02,
          8.49052835e-02,   1.25990938e-01,   1.67076592e-01,
          2.08162247e-01,   2.49247901e-01,   2.90333556e-01,
          3.31419210e-01,   3.72504865e-01,   4.13590519e-01,
          4.54676174e-01,   4.95761828e-01,   5.36847483e-01,
          5.77933138e-01,   6.19018792e-01,   6.60104447e-01,
          7.01190101e-01,   7.42275756e-01,   7.83361410e-01,
          8.24447065e-01,   8.65532719e-01,   9.06618374e-01,
          9.47704028e-01,   9.88789683e-01,   1.02987534e+00,
          1.07096099e+00,   1.11204665e+00,   1.15313230e+00,
          1.19421796e+00,   1.23530361e+00,   1.27638926e+00,
          1.31747492e+00,   1.35856057e+00,   1.39964623e+00,
          1.44073188e+00,   1.48181754e+00,   1.52290319e+00,
          1.56398885e+00,   1.60507450e+00,   1.64616015e+00,
          1.68724581e+00,   1.72833146e+00,   1.76941712e+00,
          1.81050277e+00,   1.85158843e+00,   1.89267408e+00,
          1.93375974e+00,   1.97484539e+00,   2.01593105e+00,
          2.05701670e+00,   2.09810235e+00,   2.13918801e+00,
          2.18027366e+00,   2.22135932e+00,   2.26244497e+00,
          2.30353063e+00,   2.34461628e+00,   2.38570194e+00,
          2.42678759e+00,   2.46787324e+00,   2.50895890e+00,
          2.55004455e+00,   2.59113021e+00,   2.63221586e+00,
          2.67330152e+00,   2.71438717e+00,   2.75547283e+00,
          2.79655848e+00,   2.83764414e+00,   2.87872979e+00,
          2.91981544e+00,   2.96090110e+00,   3.00198675e+00,
          3.04307241e+00,   3.08415806e+00,   3.12524372e+00,
          3.16632937e+00,   3.20741503e+00,   3.24850068e+00,
          3.28958633e+00,   3.33067199e+00,   3.37175764e+00,
          3.41284330e+00,   3.45392895e+00,   3.49501461e+00,
          3.53610026e+00,   3.57718592e+00,   3.61827157e+00,
          3.65935723e+00,   3.70044288e+00,   3.74152853e+00,
          3.78261419e+00,   3.82369984e+00,   3.86478550e+00,
          3.90587115e+00,   3.94695681e+00,   3.98804246e+00,
          4.02912812e+00]),
 masked_array(data =
  [[-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  ..., 
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]],
              mask =
  [[ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  ..., 
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]],
        fill_value = 1e+20),
 <matplotlib.contour.QuadContourSet at 0x7f98793e9780>)

То же в виде цветов.

In [12]:
a,b,g,r=m.mncontour_grid('a','b',nsigma=3)
pcolormesh(a,b,g)
colorbar()
Out[12]:
<matplotlib.colorbar.Colorbar at 0x7f98792e7cf8>

Дайте мне 3 параметра, и я профитирую слона. С 4 параметрами он будет махать хоботом.

Пусть у нас есть экспериментальные данные, и мы хотим профитировать их прямой.

In [13]:
def fit(a,b,x):
    return a*x+b

Данные не настоящие, а сгенерированные. Все имеют ошибки 0.1.

In [14]:
x=linspace(0,1,11)
dy=0.1*ones(11)
y=x+dy*normal(size=11)

Функция $\chi^2$.

In [15]:
def chi2(a,b):
    global x,y,dy
    return (((y-fit(a,b,x))/dy)**2).sum()

Минимизируем.

In [16]:
m=Minuit(chi2,a=0,b=0,error_a=1,error_b=1)
/usr/lib64/python3.4/site-packages/ipykernel/__main__.py:1: InitialParamWarning: errordef is not given. Default to 1.
  if __name__ == '__main__':
In [17]:
m.migrad()

FCN = 7.341081582910074 TOTAL NCALL = 32 NCALLS = 32
EDM = 6.20027421166483e-24 GOAL EDM = 1e-05 UP = 1.0
Valid Valid Param Accurate Covar PosDef Made PosDef
True True True True False
Hesse Fail HasCov Above EDM Reach calllim
False True False False
+ Name Value Parab Error Minos Error- Minos Error+ Limit- Limit+ FIXED
1 a 0.952919 0.0953463 0 0
2 b -0.0578366 0.0564076 0 0

Out[17]:
({'fval': 7.341081582910074, 'up': 1.0, 'is_above_max_edm': False, 'has_valid_parameters': True, 'edm': 6.20027421166483e-24, 'has_made_posdef_covar': False, 'has_covariance': True, 'has_posdef_covar': True, 'nfcn': 32, 'has_accurate_covar': True, 'hesse_failed': False, 'is_valid': True, 'has_reached_call_limit': False},
 [{'is_const': False, 'number': 0, 'upper_limit': 0.0, 'lower_limit': 0.0, 'name': 'a', 'is_fixed': False, 'error': 0.0953462587249722, 'value': 0.9529192501482413, 'has_upper_limit': False, 'has_limits': False, 'has_lower_limit': False},
  {'is_const': False, 'number': 1, 'upper_limit': 0.0, 'lower_limit': 0.0, 'name': 'b', 'is_fixed': False, 'error': 0.05640760738050268, 'value': -0.0578365609047049, 'has_upper_limit': False, 'has_limits': False, 'has_lower_limit': False}])
In [18]:
m.values
Out[18]:
{'a': 0.9529192501482413, 'b': -0.0578365609047049}
In [19]:
m.fval
Out[19]:
7.341081582910074
In [22]:
m.matrix()
Out[22]:
((0.009090909052849335, -0.004545454524510073),
 (-0.004545454524510073, 0.003181818170392941))
In [20]:
m.draw_mncontour('a','b',nsigma=3)
Out[20]:
(array([ 0.58416759,  0.59160027,  0.59903294,  0.60646562,  0.61389829,
         0.62133097,  0.62876364,  0.63619631,  0.64362899,  0.65106166,
         0.65849434,  0.66592701,  0.67335969,  0.68079236,  0.68822503,
         0.69565771,  0.70309038,  0.71052306,  0.71795573,  0.7253884 ,
         0.73282108,  0.74025375,  0.74768643,  0.7551191 ,  0.76255178,
         0.76998445,  0.77741712,  0.7848498 ,  0.79228247,  0.79971515,
         0.80714782,  0.8145805 ,  0.82201317,  0.82944584,  0.83687852,
         0.84431119,  0.85174387,  0.85917654,  0.86660922,  0.87404189,
         0.88147456,  0.88890724,  0.89633991,  0.90377259,  0.91120526,
         0.91863794,  0.92607061,  0.93350328,  0.94093596,  0.94836863,
         0.95580131,  0.96323398,  0.97066666,  0.97809933,  0.985532  ,
         0.99296468,  1.00039735,  1.00783003,  1.0152627 ,  1.02269538,
         1.03012805,  1.03756072,  1.0449934 ,  1.05242607,  1.05985875,
         1.06729142,  1.0747241 ,  1.08215677,  1.08958944,  1.09702212,
         1.10445479,  1.11188747,  1.11932014,  1.12675282,  1.13418549,
         1.14161816,  1.14905084,  1.15648351,  1.16391619,  1.17134886,
         1.17878154,  1.18621421,  1.19364688,  1.20107956,  1.20851223,
         1.21594491,  1.22337758,  1.23081026,  1.23824293,  1.2456756 ,
         1.25310828,  1.26054095,  1.26797363,  1.2754063 ,  1.28283898,
         1.29027165,  1.29770432,  1.305137  ,  1.31256967,  1.32000235]),
 array([-0.27599298, -0.27159575, -0.26719852, -0.26280129, -0.25840407,
        -0.25400684, -0.24960961, -0.24521238, -0.24081515, -0.23641792,
        -0.23202069, -0.22762346, -0.22322623, -0.218829  , -0.21443177,
        -0.21003454, -0.20563731, -0.20124008, -0.19684285, -0.19244562,
        -0.18804839, -0.18365117, -0.17925394, -0.17485671, -0.17045948,
        -0.16606225, -0.16166502, -0.15726779, -0.15287056, -0.14847333,
        -0.1440761 , -0.13967887, -0.13528164, -0.13088441, -0.12648718,
        -0.12208995, -0.11769272, -0.11329549, -0.10889827, -0.10450104,
        -0.10010381, -0.09570658, -0.09130935, -0.08691212, -0.08251489,
        -0.07811766, -0.07372043, -0.0693232 , -0.06492597, -0.06052874,
        -0.05613151, -0.05173428, -0.04733705, -0.04293982, -0.0385426 ,
        -0.03414537, -0.02974814, -0.02535091, -0.02095368, -0.01655645,
        -0.01215922, -0.00776199, -0.00336476,  0.00103247,  0.0054297 ,
         0.00982693,  0.01422416,  0.01862139,  0.02301862,  0.02741585,
         0.03181308,  0.0362103 ,  0.04060753,  0.04500476,  0.04940199,
         0.05379922,  0.05819645,  0.06259368,  0.06699091,  0.07138814,
         0.07578537,  0.0801826 ,  0.08457983,  0.08897706,  0.09337429,
         0.09777152,  0.10216875,  0.10656598,  0.1109632 ,  0.11536043,
         0.11975766,  0.12415489,  0.12855212,  0.13294935,  0.13734658,
         0.14174381,  0.14614104,  0.15053827,  0.1549355 ,  0.15933273]),
 masked_array(data =
  [[-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  ..., 
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]],
              mask =
  [[ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  ..., 
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]],
        fill_value = 1e+20),
 <matplotlib.contour.QuadContourSet at 0x7f987926d208>)
In [23]:
a,b,g,r=m.mncontour_grid('a','b',nsigma=3)
pcolormesh(a,b,g)
colorbar()
Out[23]:
<matplotlib.colorbar.Colorbar at 0x7f98791a6b00>

Нарисуем на одном графике экспериментальные точки, наш фит (сплошная линия) и истинную теоретическую кривую (пунктир).

In [24]:
errorbar(x,y,dy,fmt='ro')
xt=linspace(0,1,101)
plot(xt,fit(m.values['a'],m.values['b'],xt),'b-')
plot(xt,fit(1,0,xt),'g--')
Out[24]:
[<matplotlib.lines.Line2D at 0x7f98791be080>]

Когда фитирующая функция есть линейная комбинация каких-то фиксированных функций с неизвестными коэффициентами, минимизация $\chi^2$ сводится к решению системы линейных уравнений. Нет надобности использовать Minuit.

Резонанс без фона

Пусть теперь наша фитирующая функция - Брейт-Вигнеровский резонанс (без фона), с двумя параметрами - положением и шириной (лучше бы ввести третий - высоту, но я не стал этого делать для простоты). Теперь $\chi^2$ - сложная нелинейная функция параметров.

In [25]:
def fit(x0,Gamma,x):
    return 1/((x-x0)**2+Gamma**2)

Вот наши экспериментальные данные (с ошибками 0.1).

In [26]:
x=linspace(-3,3,21)
dy=0.1*ones(21)
y=fit(0,1,x)+dy*normal(size=21)

Минимизируем $\chi^2$.

In [27]:
def chi2(x0,Gamma):
    global x,y,dy
    return (((y-fit(x0,Gamma,x))/dy)**2).sum()
In [28]:
m=Minuit(chi2,x0=0,error_x0=1,Gamma=1,error_Gamma=1)
/usr/lib64/python3.4/site-packages/ipykernel/__main__.py:1: InitialParamWarning: errordef is not given. Default to 1.
  if __name__ == '__main__':
In [29]:
m.migrad()

FCN = 18.097376960150946 TOTAL NCALL = 33 NCALLS = 33
EDM = 1.7129910452155893e-07 GOAL EDM = 1e-05 UP = 1.0
Valid Valid Param Accurate Covar PosDef Made PosDef
True True True True False
Hesse Fail HasCov Above EDM Reach calllim
False True False False
+ Name Value Parab Error Minos Error- Minos Error+ Limit- Limit+ FIXED
1 x0 -0.0351926 0.0673422 0 0
2 Gamma 1.03425 0.0300784 0 0

Out[29]:
({'fval': 18.097376960150946, 'up': 1.0, 'is_above_max_edm': False, 'has_valid_parameters': True, 'edm': 1.7129910452155893e-07, 'has_made_posdef_covar': False, 'has_covariance': True, 'has_posdef_covar': True, 'nfcn': 33, 'has_accurate_covar': True, 'hesse_failed': False, 'is_valid': True, 'has_reached_call_limit': False},
 [{'is_const': False, 'number': 0, 'upper_limit': 0.0, 'lower_limit': 0.0, 'name': 'x0', 'is_fixed': False, 'error': 0.06734219710972826, 'value': -0.0351926437845268, 'has_upper_limit': False, 'has_limits': False, 'has_lower_limit': False},
  {'is_const': False, 'number': 1, 'upper_limit': 0.0, 'lower_limit': 0.0, 'name': 'Gamma', 'is_fixed': False, 'error': 0.030078423287488033, 'value': 1.034253164644145, 'has_upper_limit': False, 'has_limits': False, 'has_lower_limit': False}])
In [30]:
m.values
Out[30]:
{'Gamma': 1.034253164644145, 'x0': -0.0351926437845268}
In [31]:
m.fval
Out[31]:
18.097376960150946
In [32]:
m.errors
Out[32]:
{'Gamma': 0.030078423287488033, 'x0': 0.06734219710972826}
In [33]:
m.matrix()
Out[33]:
((0.004534971511565494, -2.272188453114644e-05),
 (-2.272188453114644e-05, 0.0009047115474613024))
In [34]:
m.draw_mncontour('x0','Gamma',nsigma=3)
Out[34]:
(array([-0.30116278, -0.29582729, -0.29049181, -0.28515632, -0.27982084,
        -0.27448535, -0.26914986, -0.26381438, -0.25847889, -0.25314341,
        -0.24780792, -0.24247243, -0.23713695, -0.23180146, -0.22646597,
        -0.22113049, -0.215795  , -0.21045952, -0.20512403, -0.19978854,
        -0.19445306, -0.18911757, -0.18378208, -0.1784466 , -0.17311111,
        -0.16777563, -0.16244014, -0.15710465, -0.15176917, -0.14643368,
        -0.14109819, -0.13576271, -0.13042722, -0.12509174, -0.11975625,
        -0.11442076, -0.10908528, -0.10374979, -0.0984143 , -0.09307882,
        -0.08774333, -0.08240785, -0.07707236, -0.07173687, -0.06640139,
        -0.0610659 , -0.05573041, -0.05039493, -0.04505944, -0.03972396,
        -0.03438847, -0.02905298, -0.0237175 , -0.01838201, -0.01304653,
        -0.00771104, -0.00237555,  0.00295993,  0.00829542,  0.01363091,
         0.01896639,  0.02430188,  0.02963736,  0.03497285,  0.04030834,
         0.04564382,  0.05097931,  0.0563148 ,  0.06165028,  0.06698577,
         0.07232125,  0.07765674,  0.08299223,  0.08832771,  0.0936632 ,
         0.09899869,  0.10433417,  0.10966966,  0.11500514,  0.12034063,
         0.12567612,  0.1310116 ,  0.13634709,  0.14168258,  0.14701806,
         0.15235355,  0.15768903,  0.16302452,  0.16836001,  0.17369549,
         0.17903098,  0.18436647,  0.18970195,  0.19503744,  0.20037292,
         0.20570841,  0.2110439 ,  0.21637938,  0.22171487,  0.22705036]),
 array([ 0.92950681,  0.93189365,  0.9342805 ,  0.93666734,  0.93905418,
         0.94144103,  0.94382787,  0.94621471,  0.94860156,  0.9509884 ,
         0.95337524,  0.95576209,  0.95814893,  0.96053577,  0.96292262,
         0.96530946,  0.96769631,  0.97008315,  0.97246999,  0.97485684,
         0.97724368,  0.97963052,  0.98201737,  0.98440421,  0.98679105,
         0.9891779 ,  0.99156474,  0.99395159,  0.99633843,  0.99872527,
         1.00111212,  1.00349896,  1.0058858 ,  1.00827265,  1.01065949,
         1.01304633,  1.01543318,  1.01782002,  1.02020686,  1.02259371,
         1.02498055,  1.0273674 ,  1.02975424,  1.03214108,  1.03452793,
         1.03691477,  1.03930161,  1.04168846,  1.0440753 ,  1.04646214,
         1.04884899,  1.05123583,  1.05362268,  1.05600952,  1.05839636,
         1.06078321,  1.06317005,  1.06555689,  1.06794374,  1.07033058,
         1.07271742,  1.07510427,  1.07749111,  1.07987795,  1.0822648 ,
         1.08465164,  1.08703849,  1.08942533,  1.09181217,  1.09419902,
         1.09658586,  1.0989727 ,  1.10135955,  1.10374639,  1.10613323,
         1.10852008,  1.11090692,  1.11329376,  1.11568061,  1.11806745,
         1.1204543 ,  1.12284114,  1.12522798,  1.12761483,  1.13000167,
         1.13238851,  1.13477536,  1.1371622 ,  1.13954904,  1.14193589,
         1.14432273,  1.14670958,  1.14909642,  1.15148326,  1.15387011,
         1.15625695,  1.15864379,  1.16103064,  1.16341748,  1.16580432]),
 masked_array(data =
  [[-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  ..., 
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]
  [-- -- -- ..., -- -- --]],
              mask =
  [[ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  ..., 
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]
  [ True  True  True ...,  True  True  True]],
        fill_value = 1e+20),
 <matplotlib.contour.QuadContourSet at 0x7f98790ff6d8>)
In [35]:
x0,Gamma,g,r=m.mncontour_grid('x0','Gamma',nsigma=3)
pcolormesh(x0,Gamma,g)
colorbar()
Out[35]:
<matplotlib.colorbar.Colorbar at 0x7f9879040278>

Теперь контуры постоянной высоты $\chi^2$ - уже не симметричные эллипсы с центром в оптимальной точке, а какие-то сложные кривые. Ошибки положения и ширины резонанса довольно-таки независимы.

Нарисуем на одном графике экспериментальные точки, наш фит (сплошная линия) и истинную теоретическую кривую (пунктир).

In [36]:
errorbar(x,y,dy,fmt='ro')
xt=linspace(-3.5,3.5,101)
plot(xt,fit(m.values['x0'],m.values['Gamma'],xt),'b-')
plot(xt,fit(0,1,xt),'g--')
Out[36]:
[<matplotlib.lines.Line2D at 0x7f987914b278>]