Newer
Older
R_phipi / fitter.ipynb
@Davide Lancierini Davide Lancierini on 23 Oct 2018 14 KB created fitter
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ROOT as r\n",
    "import pickle\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "l_index=0\n",
    "l_flv=['e','mu']\n",
    "\n",
    "data = np.load('/disk/lhcb_data/davide/Rphipi/selected_data/'+l_flv[l_index]+l_flv[l_index]+'/sel_data_NN_'+l_flv[l_index]+l_flv[l_index]+'.npy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "#cut on mass\n",
    "cut_indices=[]\n",
    "\n",
    "for i in range(len(data)):\n",
    "    if 1790<data[i]<2030:\n",
    "        cut_indices.append(i)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA6gAAAJCCAYAAADN6ep4AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAG/xJREFUeJzt3V+MpXd93/HPN15CJUyDkRfX8Z8MQm4TJ1INXRlUS5UrVDDeC5NKVHYlcC0q58JUIOXCE25IEyFN1Ya0KCmtIyxAIlBLYOFqrSSuBUK54M/adWyM47IyU7yxZW/qCIyQgmx+vZjHZLBnd2Z35pzznTOvl7Samec8Z+bn1bPH5z3f5zynxhgBAACARfu5RS8AAAAAEoEKAABAEwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoIVDi15Aklx44YVjZWVl0csAAABgBh544IG/HmMc3m6/FoG6srKS48ePL3oZAAAAzEBV/d+d7OcUXwAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0cWvQCAABYTiurx854+/ra0TmtBNgvTFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0MKhRS8AAID9aWX12KKXACwZE1QAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALSwbaBW1WVV9eWqeqyqHq2qD07bf7uq/qqqHpr+XL/pPr9VVSeq6vGqeucs/wMAAABYDod2sM8LSX5zjPFgVb02yQNVdd902++PMf7T5p2r6sokNyb51SS/mOR/VdU/HGO8uJcLBwAAYLlsO0EdYzw9xnhw+vz5JI8lueQMd7khyefHGH87xvhukhNJrt6LxQIAALC8zuo1qFW1kuTNSb4+bfpAVT1cVXdW1QXTtkuSPLnpbiezRdBW1a1Vdbyqjp86deqsFw4AAMBy2XGgVtX5Sb6Q5ENjjB8k+USSNyW5KsnTSX7vpV23uPt4xYYx7hhjHBljHDl8+PBZLxwAAIDlsqNArapXZSNOPzvG+GKSjDGeGWO8OMb4SZI/yt+dxnsyyWWb7n5pkqf2bskAAAAso51cxbeSfDLJY2OMj23afvGm3X49ybemz+9JcmNVvbqq3pjkiiTf2LslAwAAsIx2chXfa5K8N8kjVfXQtO3DSW6qqquycfruepLfSJIxxqNVdVeSb2fjCsC3uYIvAAAvt7J6bNt91teOzmElQBfbBuoY48+z9etK7z3DfT6a5KO7WBcAAAAHzE4mqAAAHEA7mXAC7KWzepsZAAAAmBUTVAAA2tpuius1qrBcTFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQwqFFLwAAgMVYWT226CUA/AwTVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKCFbQO1qi6rqi9X1WNV9WhVfXDa/vqquq+qvjN9vGDaXlX18ao6UVUPV9VbZv0fAQAAwP53aAf7vJDkN8cYD1bVa5M8UFX3Jfk3Se4fY6xV1WqS1SS3J3lXkiumP29N8onpIwAA7KmV1WNnvH197eicVgLshW0nqGOMp8cYD06fP5/ksSSXJLkhyaen3T6d5N3T5zck+czY8LUkr6uqi/d85QAAACyVs3oNalWtJHlzkq8nuWiM8XSyEbFJ3jDtdkmSJzfd7eS07eXf69aqOl5Vx0+dOnX2KwcAAGCp7DhQq+r8JF9I8qExxg/OtOsW28YrNoxxxxjjyBjjyOHDh3e6DAAAAJbUjgK1ql6VjTj97Bjji9PmZ146dXf6+Oy0/WSSyzbd/dIkT+3NcgEAAFhWO7mKbyX5ZJLHxhgf23TTPUlunj6/OcmXNm1/33Q137cl+f5LpwIDAADA6ezkKr7XJHlvkkeq6qFp24eTrCW5q6ren+R7Sd4z3XZvkuuTnEjyoyS37OmKAQAAWErbBuoY48+z9etKk+TtW+w/kty2y3UBAABwwJzVVXwBAABgVgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaOHQohcAAMBsrKweW/QSFm67v4P1taNzWgmwEyaoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWDi16AQCdrKwe23af9bWju/oe290fYKd28pgFsJ+YoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWtg2UKvqzqp6tqq+tWnbb1fVX1XVQ9Of6zfd9ltVdaKqHq+qd85q4QAAACyXnUxQP5Xkui22//4Y46rpz71JUlVXJrkxya9O9/mvVXXeXi0WAACA5bVtoI4xvprkuR1+vxuSfH6M8bdjjO8mOZHk6l2sDwAAgANiN69B/UBVPTydAnzBtO2SJE9u2ufktO0VqurWqjpeVcdPnTq1i2UAAACwDM41UD+R5E1JrkrydJLfm7bXFvuOrb7BGOOOMcaRMcaRw4cPn+MyAAAAWBbnFKhjjGfGGC+OMX6S5I/yd6fxnkxy2aZdL03y1O6WCAAAwEFwToFaVRdv+vLXk7x0hd97ktxYVa+uqjcmuSLJN3a3RAAAAA6CQ9vtUFWfS3Jtkgur6mSSjyS5tqquysbpu+tJfiNJxhiPVtVdSb6d5IUkt40xXpzN0gEAAFgm2wbqGOOmLTZ/8gz7fzTJR3ezKAAAAA6ebQMVgJ+1snps0UsAAFhKu3mbGQAAANgzAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaOHQohcAME8rq8cWvQQAAE7DBBUAAIAWTFCBpWJCCgCwf5mgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGjh0KIXAABwEK2sHtt2n/W1o3NYCUAfJqgAAAC0YIIKMGfbTU1MTACAg0qgAgBwYPmlIfTiFF8AAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFrwPKgBAU9u9RyfAsjFBBQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0cGjRCwAAWEYrq8cWvQSAfccEFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCt5kBaGa7t6ZYXzs6p5UAAMyXCSoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBYEKAABACwIVAACAFg4tegEAm62sHjvj7etrR+e0EgAA5s0EFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0cGjRCwAOlpXVY4teAgAATZmgAgAA0MK2gVpVd1bVs1X1rU3bXl9V91XVd6aPF0zbq6o+XlUnqurhqnrLLBcPAADA8tjJBPVTSa572bbVJPePMa5Icv/0dZK8K8kV059bk3xib5YJAADAstv2NahjjK9W1crLNt+Q5Nrp808n+UqS26ftnxljjCRfq6rXVdXFY4yn92rBAAAwL9tdO2F97eicVgIHw7m+BvWil6Jz+viGafslSZ7ctN/JaRsAAACc0V5fJKm22Da23LHq1qo6XlXHT506tcfLAAAAYL8510B9pqouTpLp47PT9pNJLtu036VJntrqG4wx7hhjHBljHDl8+PA5LgMAAIBlca6Bek+Sm6fPb07ypU3b3zddzfdtSb7v9acAAADsxLYXSaqqz2XjgkgXVtXJJB9Jspbkrqp6f5LvJXnPtPu9Sa5PciLJj5LcMoM1AwAs3HYXzwHg7O3kKr43neamt2+x70hy224XBQAAwMGzbaAC0Iu3PAAAltVeX8UXAAAAzolABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaOLToBQCcjZXVY4teAgAAM2KCCgAAQAsmqAAAW3DGBsD8maACAADQggkqAACco+0m7etrR+e0ElgOJqgAAAC0IFABAABowSm+AEBLTp1kGTiO4eyYoAIAANCCCSoAsJRMrgD2HxNUAAAAWjBBBQAOpO0mrADMnwkqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFo4tOgFALC3VlaPnfH29bWjc1oJB9l2x2HiWATglUxQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANDCoUUvAAAADqqV1WNnvH197eicVgI9mKACAADQgkAFAACgBYEKAABACwIVAACAFgQqAAAALbiKL8AB44qRAEBXAhXYU9vFDwAAnI5TfAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWXMUXAFiI3V7121XDAZaPCSoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQgkAFAACgBW8zAwB7bCdvf7K+dnQOKwGA/cUEFQAAgBZMUAE4K6aDAL1s97jsMZn9xAQVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IK3mQHOyk7eYgQA2D+8TQ2dmKACAADQggkqAPuO3/YDwHIyQQUAAKAFE1QA5s4EdPf8HQKwjExQAQAAaGFXE9SqWk/yfJIXk7wwxjhSVa9P8j+SrCRZT/Kvxhh/s7tlAsDOmS4C7B2PqczTXkxQ//kY46oxxpHp69Uk948xrkhy//Q1AAAAnNEsXoN6Q5Jrp88/neQrSW6fwc8BAE7DxAOWg/cf56DZ7QR1JPmzqnqgqm6dtl00xng6SaaPb9jlzwAAAOAA2O0E9ZoxxlNV9YYk91XVX+70jlPQ3pokl19++S6XAQCcDRNWADra1QR1jPHU9PHZJHcnuTrJM1V1cZJMH589zX3vGGMcGWMcOXz48G6WAQAAwBI45wlqVb0myc+NMZ6fPn9Hkt9Jck+Sm5OsTR+/tBcLBWD/8Jqp3fN3CMBBtJtTfC9KcndVvfR9/niM8SdV9c0kd1XV+5N8L8l7dr9MANg7Tm8FgJ7OOVDHGE8k+cdbbP9/Sd6+m0UBAABw8OzF+6ACAADArglUAAAAWtjt28wAsGRcnAdguXhcZz8xQQUAAKAFgQoAAEALAhUAAIAWvAYV9pD3VoS94fVSAHAwmaACAADQgkAFAACgBaf4QiNOEQb2C6dhA3vF8x82M0EFAACgBYEKAABACwIVAACAFrwGFfgprymDnfFvBaCPnTwmex3r/mGCCgAAQAsmqACwAN2nsN3XB8ByMkEFAACgBRNUAADgnHkfU38He8kEFQAAgBZMUGGTZf/tl9eUAQDz5vkHZ0OgAgAAB9qyDyn2E6f4AgAA0IIJKiwRp9AAALCfmaACAADQggkqzJEJJwAAnJ4JKgAAAC2YoALAyzjbAQAWwwQVAACAFgQqAAAALQhUAAAAWhCoAAAAtCBQAQAAaMFVfGEfcWVRAIDls91zvPW1o3NayeIJVOZm0f/wxB0AwP7jOdzB4hRfAAAAWjBBZWnM47drfoMHAMC8LfpMxHkyQQUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABacBVf2jhIVycDAABeyQQVAACAFkxQD4idvH+nCSUAAMvIe9nvHwJ1SczjH51TcAEAgFlyii8AAAAtCFQAAABaEKgAAAC04DWo/JQXjwMAwCt5njw/JqgAAAC0YIIKAAAwQyawO2eCCgAAQAsmqPuE37oAAADLzgQVAACAFkxQ2TdMkQEAYLkJVPaMgAQAAHbDKb4AAAC0IFABAABoQaACAADQgteg7tB2r69cXzs60/sDAAAsOxNUAAAAWjBB3SO7vYKtK+ACAAAHnQkqAAAALQhUAAAAWhCoAAAAtCBQAQAAaEGgAgAA0IJABQAAoAWBCgAAQAsCFQAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEALAhUAAIAWBCoAAAAtCFQAAABaEKgAAAC0IFABAABoQaACAADQwswCtaquq6rHq+pEVa3O6ucAAACwHGYSqFV1XpI/TPKuJFcmuamqrpzFzwIAAGA5zGqCenWSE2OMJ8YYP07y+SQ3zOhnAQAAsARmFaiXJHly09cnp20AAACwpUMz+r61xbbxMztU3Zrk1unLH1bV4zNay7K5MMlfL3oRcAaOUbpzjNKdY5TuHKPN1H9Y9Ap25Jd2stOsAvVkkss2fX1pkqc27zDGuCPJHTP6+Uurqo6PMY4seh1wOo5RunOM0p1jlO4co8zSrE7x/WaSK6rqjVX180luTHLPjH4WAAAAS2AmE9QxxgtV9YEkf5rkvCR3jjEencXPAgAAYDnM6hTfjDHuTXLvrL7/Aea0aLpzjNKdY5TuHKN05xhlZmqMsf1eAAAAMGOzeg0qAAAAnBWB2kRV3VlVz1bVt05z+y9U1f+sqr+oqker6pZNt91cVd+Z/tw8v1VzkOzyGH2xqh6a/rhgGjOxg2P0gqq6u6oerqpvVNWvbbrtuqp6vKpOVNXq/FbNQbLLY3S9qh6ZHkePz2/VHCRVdVlVfbmqHpv+X/7BLfapqvr49Hj5cFW9ZdNtnpOya07xbaKq/lmSHyb5zBjj17a4/cNJfmGMcXtVHU7yeJJ/kOT8JMeTHMnGe80+kOSfjDH+Zm6L50A412N0jPHjqvrhGOP8OS+ZA2YHx+h/TPLDMca/r6pfTvKHY4y3V9V5Sf5Pkn+RjbdJ+2aSm8YY357j8jkAzvUYnW5bT3JkjOG9J5mZqro4ycVjjAer6rXZeF757s2Ph1V1fZJ/l+T6JG9N8l/GGG+tqtfHc1L2gAlqE2OMryZ57ky7JHltVVU2ovS5JC8keWeS+8YYz00PAPcluW7W6+Xg2cUxCnOxg2P0yiT3T/v+ZZKVqrooydVJTowxnhhj/DjJ55PcMOv1cvDs4hiFuRhjPD3GeHD6/PkkjyW55GW73ZCNX7KMMcbXkrxuClvPSdkTAnX/+IMkv5LkqSSPJPngGOMn2XjQeHLTfifzygcSmIfTHaNJ8veq6nhVfa2q3r2wFXLQ/UWSf5kkVXV1kl9Kcmk8jtLH6Y7RZOOXgH9WVQ9U1a0LWh8HSFWtJHlzkq+/7KbTPWZ6LGVPCNT9451JHkryi0muSvIHVfX3k9QW+zpvm0U43TGaJJePMY4k+ddJ/nNVvWlBa+RgW0tyQVU9lI3T0/53Nqb8Hkfp4nTHaJJcM8Z4S5J3JbltOl0YZqKqzk/yhSQfGmP84OU3b3GXcYbtcFYE6v5xS5IvTqdTnEjy3SS/nI3fTl22ab9LszHBgnk73TGaMcZT08cnknwlG7+RhbkaY/xgjHHLGOOqJO9Lcjgbx6nHUVo4wzG6+XH02SR3Z+PUdNhzVfWqbMTpZ8cYX9xil9M9ZnosZU8I1P3je0leulDCRUn+UZInkvxpkndMV/67IMk7pm0wb1seo9Ox+epp+4VJrkni4jPMXVW9rqp+fvry3yb56jQZ+GaSK6rqjdPtNyZxtWnm7nTHaFW9ZrpgTarqNdn4f/2WVwKG3ZiuI/HJJI+NMT52mt3uSfK+6Wq+b0vy/THG0/GclD1yaNELYENVfS7JtUkurKqTST6S5FVJMsb4b0l+N8mnquqRbJxCcftLV/Krqt/NxhOsJPmdMcaZLsAA5+Rcj9Gq+qdJ/ntV/SQbvxRbc3VUZmEHx+ivJPlMVb2YjV+SvH+67YWq+kA2nkidl+TOMcaj8/8vYNmd6zGa5KIkd2+0Qw4l+eMxxp/Md/UcENckeW+SR6ZTzZPkw0kuT356nN6bjSv4nkjyo2ycQZUxxnOek7IXvM0MAAAALTjFFwAAgBYEKgAAAC0IVAAAAFoQqAAAALQgUAEAAGhBoAIAANCCQAUAAKAFgQoAAEAL/x+YSyWT5ZcpRgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "data_cut=data[cut_indices]/1000\n",
    "plt.hist(data_cut,range=(1.79,2.03),bins=100);\n",
    "fig=plt.gcf();\n",
    "fig.set_size_inches(16,10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "h1 = r.TH1D(\"h1\",\"data histogram\",100,1.79,2.03)\n",
    "for i in range(len(data_cut)): h1.Fill(data_cut[i])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = r.RooRealVar(\"m\",\"m\",1.8,2.5)\n",
    "l = r.RooArgList(m)\n",
    "data_hist = r.RooDataHist(\"data histogram\", \"data\", l, h1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Creating the signal PDF\n",
    "\n",
    "\n",
    "mean_dplus= r.RooRealVar(\"mean_D+\",\"mean_D+\",1.87,1.77,1.97)\n",
    "sigma_dplus = r.RooRealVar(\"sigma_D+\",\"sigma_D+\",0.020,0.,0.2)\n",
    "sig_dplus = r.RooGaussian(\"signal_D+\",\"signal_D+\", m, mean_dplus, sigma_dplus)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Creating the signal PDF\n",
    "\n",
    "mean_ds= r.RooRealVar(\"mean_Ds\",\"mean_Ds\",1.97,1.87,2.07)\n",
    "sigma_ds = r.RooRealVar(\"sigma_Ds\",\"sigma_Ds\",0.020,0.,0.2)\n",
    "sig_ds = r.RooGaussian(\"signal_Ds\",\"signal_Ds\", m, mean_ds, sigma_ds)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "coef0 = r.RooRealVar(\"c0\",\"coefficient #0\",1.0,-1.,1)\n",
    "coef1 = r.RooRealVar(\"c1\",\"coefficient #1\",0.1,-1.,1)\n",
    "coef2 = r.RooRealVar(\"c2\",\"coefficient #2\",-0.1,-1.,1)\n",
    "bkg = r.RooChebychev(\"bkg\",\"background p.d.f.\",m,r.RooArgList(coef0,coef1,coef2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Add 2 sources of signal\n",
    "\n",
    "# first add the sig and the peak together with fraction fpeak\n",
    "fpeak = r.RooRealVar(\"fpeak\",\"peaking signals fraction\",0.5,0.,1.)\n",
    "# sigpeak(x) = fpeak*sig(x) + (1-fpeak)*bkg(x)\n",
    "sigpeaks = r.RooAddPdf(\"sigpeak\",\"two mass peaks\",r.RooArgList(sig_ds,sig_dplus),r.RooArgList(fpeak))\n",
    "\n",
    "\n",
    "# now we can add (sig+peak) to bkg with the fraction fpeak\n",
    "fbkg = r.RooRealVar(\"fbkg\",\"background fraction\",0.5,0.,1.)\n",
    "# model2(x) = fbkg*bkg(x) + (1-fbkg)*sigpeak(x)\n",
    "model = r.RooAddPdf(\"model\",\"bkg+(sig_ds+sig_dplus)\",r.RooArgList(bkg,sigpeaks),r.RooArgList(fbkg))\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "##Plotting composite models\n",
    "#\n",
    "#xframe = m.frame()\n",
    "#sigpeaks.plotOn(xframe)\n",
    "#bkg_component = r.RooArgSet(bkg)\n",
    "#model.plotOn(xframe,r.RooFit.Components(bkg_component),r.RooFit.LineStyle(2))\n",
    "#\n",
    "#xframe.Draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.fitTo(data_hist,r.RooFit.Extended(),r.RooFit.PrintLevel(-1))\n",
    "\n",
    "\n",
    "xframe = m.frame(r.RooFit.Title(\"Fit to muon data\"))\n",
    "sigpeaks.plotOn(xframe)\n",
    "bkg_component = r.RooArgSet(bkg)\n",
    "model.plotOn(xframe,r.RooFit.Components(bkg_component),r.RooFit.LineStyle(2))\n",
    "data_hist.plotOn(xframe)\n",
    "model.plotOn(xframe)\n",
    "xframe.Draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}