import numpy as np
import matplotlib.pyplot as plt
from iminuit import Minuit

evt = np.load('dimuon.npy')[:200]
evt = evt[abs(evt-3.1)<0.5]

xmin, xmax, xbinwidth = 2.6, 3.6, 0.01
vy,edges = np.histogram(evt, bins=100, range=(xmin,xmax))
vx = 0.5*(edges[1:]+edges[:-1])
vyerr = vy**0.5

def model(x, norm, mean, sigma, c1):
    
    linear = (1. + c1*x)/((xmax-xmin) + c1*(xmax**2-xmin**2)/2.)
    
    gaussian = 1./(2.*np.pi)**0.5/sigma * \
        np.exp(-0.5*((x-mean)/sigma)**2)
    
    fs = norm/len(evt)
    return fs*gaussian + (1.-fs)*linear

def fcn(norm, mean, sigma, c1):
    
    L = model(evt, norm, mean, sigma, c1)
    if np.any(L<=0.): return 1E100

    return -2.*np.log(L).sum()

m = Minuit(fcn, norm=40., mean=3.09, sigma=0.04, c1=0.)
m.migrad()
m.minos()
m.print_param()

fig = plt.figure(figsize=(6,6), dpi=80)
m.draw_mnprofile('norm', bins=1000, subtract_min=True)
plt.show()