Newer
Older
Master_thesis / raremodel-nb.ipynb
{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Import"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\util\\execution.py:57: UserWarning: Not running on Linux. Determining available cpus for thread can failand be overestimated. Workaround (only if too many cpus are used):`zfit.run.set_n_cpu(your_cpu_number)`\n",
      "  warnings.warn(\"Not running on Linux. Determining available cpus for thread can fail\"\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "WARNING: The TensorFlow contrib module will not be included in TensorFlow 2.0.\n",
      "For more information, please see:\n",
      "  * https://github.com/tensorflow/community/blob/master/rfcs/20180907-contrib-sunset.md\n",
      "  * https://github.com/tensorflow/addons\n",
      "If you depend on functionality not listed there, please file an issue.\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import os\n",
    "\n",
    "# os.environ[\"CUDA_VISIBLE_DEVICES\"] = \"-1\"\n",
    "\n",
    "import numpy as np\n",
    "from pdg_const import pdg\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "import pickle as pkl\n",
    "import sys\n",
    "import time\n",
    "from helperfunctions import display_time, prepare_plot\n",
    "import cmath as c\n",
    "import scipy.integrate as integrate\n",
    "from scipy.optimize import fminbound\n",
    "from array import array as arr\n",
    "import collections\n",
    "from itertools import compress\n",
    "import tensorflow as tf\n",
    "import zfit\n",
    "from zfit import ztf\n",
    "# from IPython.display import clear_output\n",
    "import os\n",
    "import tensorflow_probability as tfp\n",
    "tfd = tfp.distributions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# chunksize = 10000\n",
    "# zfit.run.chunking.active = True\n",
    "# zfit.run.chunking.max_n_points = chunksize"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Build model and graphs\n",
    "## Create graphs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def formfactor(q2, subscript, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2): #returns real value\n",
    "    #check if subscript is viable\n",
    "\n",
    "    if subscript != \"0\" and subscript != \"+\" and subscript != \"T\":\n",
    "        raise ValueError('Wrong subscript entered, choose either 0, + or T')\n",
    "\n",
    "    #get constants\n",
    "\n",
    "    mK = ztf.constant(pdg['Ks_M'])\n",
    "    mbstar0 = ztf.constant(pdg[\"mbstar0\"])\n",
    "    mbstar = ztf.constant(pdg[\"mbstar\"])\n",
    "\n",
    "\n",
    "    mmu = ztf.constant(pdg['muon_M'])\n",
    "    mb = ztf.constant(pdg['bquark_M'])\n",
    "    ms = ztf.constant(pdg['squark_M'])\n",
    "    mB = ztf.constant(pdg['Bplus_M'])\n",
    "\n",
    "    #N comes from derivation in paper\n",
    "\n",
    "    N = 3\n",
    "\n",
    "    #some helperfunctions\n",
    "\n",
    "    tpos = (mB - mK)**2\n",
    "    tzero = (mB + mK)*(ztf.sqrt(mB)-ztf.sqrt(mK))**2\n",
    "\n",
    "    z_oben = ztf.sqrt(tpos - q2) - ztf.sqrt(tpos - tzero)\n",
    "    z_unten = ztf.sqrt(tpos - q2) + ztf.sqrt(tpos - tzero)\n",
    "    z = tf.divide(z_oben, z_unten)\n",
    "\n",
    "    #calculate f0\n",
    "\n",
    "    if subscript == \"0\":\n",
    "        prefactor = 1/(1 - q2/(mbstar0**2))\n",
    "        _sum = 0\n",
    "        b0 = [b0_0, b0_1, b0_2]\n",
    "\n",
    "        for i in range(N):\n",
    "            _sum += b0[i]*(tf.pow(z,i))\n",
    "\n",
    "        return ztf.to_complex(prefactor * _sum)\n",
    "\n",
    "    #calculate f+ or fT\n",
    "\n",
    "    else:\n",
    "        prefactor = 1/(1 - q2/(mbstar**2))\n",
    "        _sum = 0\n",
    "\n",
    "        if subscript == \"T\":\n",
    "            bT = [bT_0, bT_1, bT_2]\n",
    "            for i in range(N):\n",
    "                _sum += bT[i] * (tf.pow(z, i) - ((-1)**(i-N)) * (i/N) * tf.pow(z, N))\n",
    "        else:\n",
    "            bplus = [bplus_0, bplus_1, bplus_2]\n",
    "            for i in range(N):\n",
    "                _sum += bplus[i] * (tf.pow(z, i) - ((-1)**(i-N)) * (i/N) * tf.pow(z, N))\n",
    "\n",
    "        return ztf.to_complex(prefactor * _sum)\n",
    "\n",
    "def resonance(q, _mass, width, phase, scale):\n",
    "\n",
    "    q2 = tf.pow(q, 2)\n",
    "\n",
    "    mmu = ztf.constant(pdg['muon_M'])\n",
    "\n",
    "    p = 0.5 * ztf.sqrt(q2 - 4*(mmu**2))\n",
    "\n",
    "    p0 =  0.5 * ztf.sqrt(_mass**2 - 4*mmu**2)\n",
    "\n",
    "    gamma_j = tf.divide(p, q) * _mass * width / p0\n",
    "\n",
    "    #Calculate the resonance\n",
    "\n",
    "    _top = tf.complex(_mass * width, ztf.constant(0.0))\n",
    "\n",
    "    _bottom = tf.complex(_mass**2 - q2, -_mass*gamma_j)\n",
    "\n",
    "    com = _top/_bottom\n",
    "\n",
    "    #Rotate by the phase\n",
    "\n",
    "    r = ztf.to_complex(scale*tf.abs(com))\n",
    "\n",
    "    _phase = tf.angle(com)\n",
    "\n",
    "    _phase += phase\n",
    "\n",
    "    com = r * tf.exp(tf.complex(ztf.constant(0.0), _phase))\n",
    "\n",
    "    return com\n",
    "\n",
    "\n",
    "def axiv_nonres(q, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2):\n",
    "\n",
    "    GF = ztf.constant(pdg['GF'])\n",
    "    alpha_ew = ztf.constant(pdg['alpha_ew'])\n",
    "    Vtb = ztf.constant(pdg['Vtb'])\n",
    "    Vts = ztf.constant(pdg['Vts'])\n",
    "    C10eff = ztf.constant(pdg['C10eff'])\n",
    "\n",
    "    mmu = ztf.constant(pdg['muon_M'])\n",
    "    mb = ztf.constant(pdg['bquark_M'])\n",
    "    ms = ztf.constant(pdg['squark_M'])\n",
    "    mK = ztf.constant(pdg['Ks_M'])\n",
    "    mB = ztf.constant(pdg['Bplus_M'])\n",
    "\n",
    "    q2 = tf.pow(q, 2)\n",
    "\n",
    "    #Some helperfunctions\n",
    "\n",
    "    beta = 1. - 4. * mmu**2. / q2\n",
    "\n",
    "    kabs = ztf.sqrt(mB**2. + tf.pow(q2, 2)/mB**2. + mK**4./mB**2. - 2. * (mB**2. * mK**2. + mK**2. * q2 + mB**2. * q2) / mB**2.)\n",
    "\n",
    "    #prefactor in front of whole bracket\n",
    "\n",
    "    prefactor1 = GF**2. *alpha_ew**2. * (tf.abs(Vtb*Vts))**2. * kabs * beta / (128. * np.pi**5.)\n",
    "\n",
    "    #left term in bracket\n",
    "\n",
    "    bracket_left = 2./3. * tf.pow(kabs,2) * tf.pow(beta,2) * tf.pow(tf.abs(ztf.to_complex(C10eff)*formfactor(q2, \"+\", b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)),2)\n",
    "\n",
    "    #middle term in bracket\n",
    "\n",
    "    _top = 4. * mmu**2. * (mB**2. - mK**2.) * (mB**2. - mK**2.)\n",
    "\n",
    "    _under = q2 * mB**2.\n",
    "\n",
    "    bracket_middle = _top/_under *tf.pow(tf.abs(ztf.to_complex(C10eff) * formfactor(q2, \"0\", b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)), 2)\n",
    "    \n",
    "    #Note sqrt(q2) comes from derivation as we use q2 and plot q\n",
    "\n",
    "    return prefactor1 * (bracket_left + bracket_middle) * 2 * q\n",
    "\n",
    "def vec(q, funcs, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2):\n",
    "    \n",
    "    q2 = tf.pow(q, 2)\n",
    "\n",
    "    GF = ztf.constant(pdg['GF'])\n",
    "    alpha_ew = ztf.constant(pdg['alpha_ew'])\n",
    "    Vtb = ztf.constant(pdg['Vtb'])\n",
    "    Vts = ztf.constant(pdg['Vts'])\n",
    "    C7eff = ztf.constant(pdg['C7eff'])\n",
    "\n",
    "    mmu = ztf.constant(pdg['muon_M'])\n",
    "    mb = ztf.constant(pdg['bquark_M'])\n",
    "    ms = ztf.constant(pdg['squark_M'])\n",
    "    mK = ztf.constant(pdg['Ks_M'])\n",
    "    mB = ztf.constant(pdg['Bplus_M'])\n",
    "\n",
    "    #Some helperfunctions\n",
    "\n",
    "    beta = 1. - 4. * mmu**2. / q2\n",
    "\n",
    "    kabs = ztf.sqrt(mB**2. + tf.pow(q2, 2)/mB**2. + mK**4./mB**2. - 2 * (mB**2 * mK**2 + mK**2 * q2 + mB**2 * q2) / mB**2)\n",
    "    \n",
    "    #prefactor in front of whole bracket\n",
    "\n",
    "    prefactor1 = GF**2. *alpha_ew**2. * (tf.abs(Vtb*Vts))**2 * kabs * beta / (128. * np.pi**5.)\n",
    "\n",
    "    #right term in bracket\n",
    "\n",
    "    prefactor2 = tf.pow(kabs,2) * (1. - 1./3. * beta)\n",
    "\n",
    "    abs_bracket = tf.pow(tf.abs(c9eff(q, funcs) * formfactor(q2, \"+\", b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2) + ztf.to_complex(2.0 * C7eff * (mb + ms)/(mB + mK)) * formfactor(q2, \"T\", b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)),2)\n",
    "\n",
    "    bracket_right = prefactor2 * abs_bracket\n",
    "\n",
    "    #Note sqrt(q2) comes from derivation as we use q2 and plot q\n",
    "\n",
    "    return prefactor1 * bracket_right * 2 * q\n",
    "\n",
    "def c9eff(q, funcs):\n",
    "\n",
    "    C9eff_nr = ztf.to_complex(ztf.constant(pdg['C9eff']))\n",
    "\n",
    "    c9 = C9eff_nr + funcs\n",
    "\n",
    "    return c9"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def G(y):\n",
    "    \n",
    "    def inner_rect_bracket(q):\n",
    "        return tf.log(ztf.to_complex((1+tf.sqrt(q))/(1-tf.sqrt(q)))-tf.complex(ztf.constant(0), -1*ztf.constant(np.pi)))    \n",
    "    \n",
    "    def inner_right(q):\n",
    "        return ztf.to_complex(2 * tf.atan(1/tf.sqrt(tf.math.real(-q))))\n",
    "    \n",
    "    big_bracket = tf.where(tf.math.real(y) > ztf.constant(0.0), inner_rect_bracket(y), inner_right(y))\n",
    "    \n",
    "    return ztf.to_complex(tf.sqrt(tf.abs(y))) * big_bracket\n",
    "\n",
    "def h_S(m, q):\n",
    "    \n",
    "    return ztf.to_complex(2) - G(ztf.to_complex(1) - ztf.to_complex(4*tf.pow(m, 2)) / ztf.to_complex(tf.pow(q, 2)))\n",
    "\n",
    "def h_P(m, q):\n",
    "    \n",
    "    return ztf.to_complex(2/3) + (ztf.to_complex(1) - ztf.to_complex(4*tf.pow(m, 2)) / ztf.to_complex(tf.pow(q, 2))) * h_S(m,q)\n",
    "\n",
    "def two_p_ccbar(mD, m_D_bar, m_D_star, q):\n",
    "    \n",
    "    \n",
    "    #Load constants\n",
    "    nu_D_bar = ztf.to_complex(pdg[\"nu_D_bar\"])\n",
    "    nu_D = ztf.to_complex(pdg[\"nu_D\"])\n",
    "    nu_D_star = ztf.to_complex(pdg[\"nu_D_star\"])\n",
    "    \n",
    "    phase_D_bar = ztf.to_complex(pdg[\"phase_D_bar\"])\n",
    "    phase_D = ztf.to_complex(pdg[\"phase_D\"])\n",
    "    phase_D_star = ztf.to_complex(pdg[\"phase_D_star\"])\n",
    "    \n",
    "    #Calculation\n",
    "    left_part =  nu_D_bar * tf.exp(tf.complex(ztf.constant(0.0), phase_D_bar)) * h_S(m_D_bar, q) \n",
    "    \n",
    "    right_part_D = nu_D * tf.exp(tf.complex(ztf.constant(0.0), phase_D)) * h_P(m_D, q) \n",
    "    \n",
    "    right_part_D_star = nu_D_star * tf.exp(tf.complex(ztf.constant(0.0), phase_D_star)) * h_P(m_D_star, q) \n",
    "\n",
    "    return left_part + right_part_D + right_part_D_star"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Build pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "class total_pdf_cut(zfit.pdf.ZPDF):\n",
    "    _N_OBS = 1  # dimension, can be omitted\n",
    "    _PARAMS = ['b0_0', 'b0_1', 'b0_2', \n",
    "               'bplus_0', 'bplus_1', 'bplus_2', \n",
    "               'bT_0', 'bT_1', 'bT_2', \n",
    "               'rho_mass', 'rho_scale', 'rho_phase', 'rho_width',\n",
    "               'jpsi_mass', 'jpsi_scale', 'jpsi_phase', 'jpsi_width',\n",
    "               'psi2s_mass', 'psi2s_scale', 'psi2s_phase', 'psi2s_width',\n",
    "               'p3770_mass', 'p3770_scale', 'p3770_phase', 'p3770_width',\n",
    "               'p4040_mass', 'p4040_scale', 'p4040_phase', 'p4040_width',\n",
    "               'p4160_mass', 'p4160_scale', 'p4160_phase', 'p4160_width',\n",
    "               'p4415_mass', 'p4415_scale', 'p4415_phase', 'p4415_width',\n",
    "               'omega_mass', 'omega_scale', 'omega_phase', 'omega_width',\n",
    "               'phi_mass', 'phi_scale', 'phi_phase', 'phi_width',\n",
    "               'Dbar_mass', 'Dbar_scale', 'Dbar_phase',\n",
    "               'Dstar_mass', 'DDstar_scale', 'DDstar_phase', 'D_mass',\n",
    "               'tau_mass', 'C_tt']\n",
    "# the name of the parameters\n",
    "\n",
    "    def _unnormalized_pdf(self, x):\n",
    "        \n",
    "        x = x.unstack_x()\n",
    "        \n",
    "        b0 = [self.params['b0_0'], self.params['b0_1'], self.params['b0_2']]\n",
    "        bplus = [self.params['bplus_0'], self.params['bplus_1'], self.params['bplus_2']]\n",
    "        bT = [self.params['bT_0'], self.params['bT_1'], self.params['bT_2']]\n",
    "        \n",
    "        def rho_res(q):\n",
    "            return resonance(q, _mass = self.params['rho_mass'], scale = self.params['rho_scale'],\n",
    "                             phase = self.params['rho_phase'], width = self.params['rho_width'])\n",
    "    \n",
    "        def omega_res(q):\n",
    "            return resonance(q, _mass = self.params['omega_mass'], scale = self.params['omega_scale'],\n",
    "                             phase = self.params['omega_phase'], width = self.params['omega_width'])\n",
    "        \n",
    "        def phi_res(q):\n",
    "            return resonance(q, _mass = self.params['phi_mass'], scale = self.params['phi_scale'],\n",
    "                             phase = self.params['phi_phase'], width = self.params['phi_width'])\n",
    "\n",
    "        def jpsi_res(q):\n",
    "            return  ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['jpsi_mass'], 2)) * resonance(q, _mass = self.params['jpsi_mass'], \n",
    "                                                                                  scale = self.params['jpsi_scale'],\n",
    "                                                                                  phase = self.params['jpsi_phase'], \n",
    "                                                                                  width = self.params['jpsi_width'])\n",
    "        def psi2s_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['psi2s_mass'], 2)) * resonance(q, _mass = self.params['psi2s_mass'], \n",
    "                                                                                   scale = self.params['psi2s_scale'],\n",
    "                                                                                   phase = self.params['psi2s_phase'], \n",
    "                                                                                   width = self.params['psi2s_width'])\n",
    "        def p3770_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p3770_mass'], 2)) * resonance(q, _mass = self.params['p3770_mass'], \n",
    "                                                                                   scale = self.params['p3770_scale'],\n",
    "                                                                                   phase = self.params['p3770_phase'], \n",
    "                                                                                   width = self.params['p3770_width'])\n",
    "        \n",
    "        def p4040_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4040_mass'], 2)) * resonance(q, _mass = self.params['p4040_mass'], \n",
    "                                                                                   scale = self.params['p4040_scale'],\n",
    "                                                                                   phase = self.params['p4040_phase'], \n",
    "                                                                                   width = self.params['p4040_width'])\n",
    "        \n",
    "        def p4160_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4160_mass'], 2)) * resonance(q, _mass = self.params['p4160_mass'], \n",
    "                                                                                   scale = self.params['p4160_scale'],\n",
    "                                                                                   phase = self.params['p4160_phase'], \n",
    "                                                                                   width = self.params['p4160_width'])\n",
    "        \n",
    "        def p4415_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4415_mass'], 2)) * resonance(q, _mass = self.params['p4415_mass'], \n",
    "                                                                                   scale = self.params['p4415_scale'],\n",
    "                                                                                   phase = self.params['p4415_phase'], \n",
    "                                                                                   width = self.params['p4415_width'])\n",
    "        \n",
    "        def P2_D(q):\n",
    "            Dbar_contrib = ztf.to_complex(self.params['Dbar_scale'])*tf.exp(tf.complex(ztf.constant(0.0), self.params['Dbar_phase']))*ztf.to_complex(h_S(self.params['Dbar_mass'], q))\n",
    "            DDstar_contrib = ztf.to_complex(self.params['DDstar_scale'])*tf.exp(tf.complex(ztf.constant(0.0), self.params['DDstar_phase']))*(ztf.to_complex(h_P(self.params['Dstar_mass'], q)) + ztf.to_complex(h_P(self.params['D_mass'], q)))\n",
    "            return Dbar_contrib + DDstar_contrib\n",
    "        \n",
    "        def ttau_cusp(q):\n",
    "            return ztf.to_complex(self.params['C_tt'])*(ztf.to_complex((h_S(self.params['tau_mass'], q))) - ztf.to_complex(h_P(self.params['tau_mass'], q)))\n",
    "        \n",
    "\n",
    "        funcs = rho_res(x) + omega_res(x) + phi_res(x) + jpsi_res(x) + psi2s_res(x) + p3770_res(x) + p4040_res(x)+ p4160_res(x) + p4415_res(x) + P2_D(x) + ttau_cusp(x)\n",
    "\n",
    "        vec_f = vec(x, funcs, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)\n",
    "\n",
    "        axiv_nr = axiv_nonres(x, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)\n",
    "\n",
    "        tot = vec_f + axiv_nr\n",
    "        \n",
    "        #Cut out jpsi and psi2s\n",
    "        \n",
    "        tot = tf.where(tf.math.logical_or(x < ztf.constant(jpsi_mass-60.), x > ztf.constant(jpsi_mass+70.)), tot, 0.0*tot)\n",
    "        \n",
    "        tot = tf.where(tf.math.logical_or(x < ztf.constant(psi2s_mass-50.), x > ztf.constant(psi2s_mass+50.)), tot, 0.0*tot)\n",
    "        \n",
    "        return tot\n",
    "    \n",
    "class total_pdf_full(zfit.pdf.ZPDF):\n",
    "    _N_OBS = 1  # dimension, can be omitted\n",
    "    _PARAMS = ['b0_0', 'b0_1', 'b0_2', \n",
    "               'bplus_0', 'bplus_1', 'bplus_2', \n",
    "               'bT_0', 'bT_1', 'bT_2', \n",
    "               'rho_mass', 'rho_scale', 'rho_phase', 'rho_width',\n",
    "               'jpsi_mass', 'jpsi_scale', 'jpsi_phase', 'jpsi_width',\n",
    "               'psi2s_mass', 'psi2s_scale', 'psi2s_phase', 'psi2s_width',\n",
    "               'p3770_mass', 'p3770_scale', 'p3770_phase', 'p3770_width',\n",
    "               'p4040_mass', 'p4040_scale', 'p4040_phase', 'p4040_width',\n",
    "               'p4160_mass', 'p4160_scale', 'p4160_phase', 'p4160_width',\n",
    "               'p4415_mass', 'p4415_scale', 'p4415_phase', 'p4415_width',\n",
    "               'omega_mass', 'omega_scale', 'omega_phase', 'omega_width',\n",
    "               'phi_mass', 'phi_scale', 'phi_phase', 'phi_width',\n",
    "               'Dbar_mass', 'Dbar_scale', 'Dbar_phase',\n",
    "               'Dstar_mass', 'DDstar_scale', 'DDstar_phase', 'D_mass',\n",
    "               'tau_mass', 'C_tt']\n",
    "# the name of the parameters\n",
    "\n",
    "    def _unnormalized_pdf(self, x):\n",
    "        \n",
    "        x = x.unstack_x()\n",
    "        \n",
    "        b0 = [self.params['b0_0'], self.params['b0_1'], self.params['b0_2']]\n",
    "        bplus = [self.params['bplus_0'], self.params['bplus_1'], self.params['bplus_2']]\n",
    "        bT = [self.params['bT_0'], self.params['bT_1'], self.params['bT_2']]\n",
    "        \n",
    "        def rho_res(q):\n",
    "            return resonance(q, _mass = self.params['rho_mass'], scale = self.params['rho_scale'],\n",
    "                             phase = self.params['rho_phase'], width = self.params['rho_width'])\n",
    "    \n",
    "        def omega_res(q):\n",
    "            return resonance(q, _mass = self.params['omega_mass'], scale = self.params['omega_scale'],\n",
    "                             phase = self.params['omega_phase'], width = self.params['omega_width'])\n",
    "        \n",
    "        def phi_res(q):\n",
    "            return resonance(q, _mass = self.params['phi_mass'], scale = self.params['phi_scale'],\n",
    "                             phase = self.params['phi_phase'], width = self.params['phi_width'])\n",
    "\n",
    "        def jpsi_res(q):\n",
    "            return  ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['jpsi_mass'], 2)) * resonance(q, _mass = self.params['jpsi_mass'], \n",
    "                                                                                  scale = self.params['jpsi_scale'],\n",
    "                                                                                  phase = self.params['jpsi_phase'], \n",
    "                                                                                  width = self.params['jpsi_width'])\n",
    "        def psi2s_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['psi2s_mass'], 2)) * resonance(q, _mass = self.params['psi2s_mass'], \n",
    "                                                                                   scale = self.params['psi2s_scale'],\n",
    "                                                                                   phase = self.params['psi2s_phase'], \n",
    "                                                                                   width = self.params['psi2s_width'])\n",
    "        def p3770_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p3770_mass'], 2)) * resonance(q, _mass = self.params['p3770_mass'], \n",
    "                                                                                   scale = self.params['p3770_scale'],\n",
    "                                                                                   phase = self.params['p3770_phase'], \n",
    "                                                                                   width = self.params['p3770_width'])\n",
    "        \n",
    "        def p4040_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4040_mass'], 2)) * resonance(q, _mass = self.params['p4040_mass'], \n",
    "                                                                                   scale = self.params['p4040_scale'],\n",
    "                                                                                   phase = self.params['p4040_phase'], \n",
    "                                                                                   width = self.params['p4040_width'])\n",
    "        \n",
    "        def p4160_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4160_mass'], 2)) * resonance(q, _mass = self.params['p4160_mass'], \n",
    "                                                                                   scale = self.params['p4160_scale'],\n",
    "                                                                                   phase = self.params['p4160_phase'], \n",
    "                                                                                   width = self.params['p4160_width'])\n",
    "        \n",
    "        def p4415_res(q):\n",
    "            return ztf.to_complex(tf.pow(q, 2) / tf.pow(self.params['p4415_mass'], 2)) * resonance(q, _mass = self.params['p4415_mass'], \n",
    "                                                                                   scale = self.params['p4415_scale'],\n",
    "                                                                                   phase = self.params['p4415_phase'], \n",
    "                                                                                   width = self.params['p4415_width'])\n",
    "        \n",
    "        def P2_D(q):\n",
    "            Dbar_contrib = ztf.to_complex(self.params['Dbar_scale'])*tf.exp(tf.complex(ztf.constant(0.0), self.params['Dbar_phase']))*ztf.to_complex(h_S(self.params['Dbar_mass'], q))\n",
    "            DDstar_contrib = ztf.to_complex(self.params['DDstar_scale'])*tf.exp(tf.complex(ztf.constant(0.0), self.params['DDstar_phase']))*(ztf.to_complex(h_P(self.params['Dstar_mass'], q)) + ztf.to_complex(h_P(self.params['D_mass'], q)))\n",
    "            return Dbar_contrib + DDstar_contrib\n",
    "        \n",
    "        def ttau_cusp(q):\n",
    "            return ztf.to_complex(self.params['C_tt'])*(ztf.to_complex((h_S(self.params['tau_mass'], q))) - ztf.to_complex(h_P(self.params['tau_mass'], q)))\n",
    "        \n",
    "\n",
    "        funcs = rho_res(x) + omega_res(x) + phi_res(x) + jpsi_res(x) + psi2s_res(x) + p3770_res(x) + p4040_res(x)+ p4160_res(x) + p4415_res(x) + P2_D(x) + ttau_cusp(x)\n",
    "\n",
    "        vec_f = vec(x, funcs, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)\n",
    "\n",
    "        axiv_nr = axiv_nonres(x, b0_0, b0_1, b0_2, bplus_0, bplus_1, bplus_2, bT_0, bT_1, bT_2)\n",
    "\n",
    "        tot = vec_f + axiv_nr\n",
    "        \n",
    "        #Cut out jpsi and psi2s\n",
    "        \n",
    "#         tot = tf.where(tf.math.logical_or(x < ztf.constant(jpsi_mass-60.), x > ztf.constant(jpsi_mass+70.)), tot, 0.0*tot)\n",
    "        \n",
    "#         tot = tf.where(tf.math.logical_or(x < ztf.constant(psi2s_mass-50.), x > ztf.constant(psi2s_mass+50.)), tot, 0.0*tot)\n",
    "        \n",
    "        return tot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\resource_variable_ops.py:435: colocate_with (from tensorflow.python.framework.ops) is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "Colocations handled automatically by placer.\n"
     ]
    }
   ],
   "source": [
    "# formfactors\n",
    "\n",
    "b0_0 = zfit.Parameter(\"b0_0\", ztf.constant(0.292), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "b0_1 = zfit.Parameter(\"b0_1\", ztf.constant(0.281), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "b0_2 = zfit.Parameter(\"b0_2\", ztf.constant(0.150), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "\n",
    "bplus_0 = zfit.Parameter(\"bplus_0\", ztf.constant(0.466), lower_limit = -2.0, upper_limit= 2.0)\n",
    "bplus_1 = zfit.Parameter(\"bplus_1\", ztf.constant(-0.885), lower_limit = -2.0, upper_limit= 2.0)\n",
    "bplus_2 = zfit.Parameter(\"bplus_2\", ztf.constant(-0.213), lower_limit = -2.0, upper_limit= 2.0)\n",
    "\n",
    "bT_0 = zfit.Parameter(\"bT_0\", ztf.constant(0.460), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "bT_1 = zfit.Parameter(\"bT_1\", ztf.constant(-1.089), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "bT_2 = zfit.Parameter(\"bT_2\", ztf.constant(-1.114), floating = False) #, lower_limit = -2.0, upper_limit= 2.0)\n",
    "\n",
    "\n",
    "#rho\n",
    "\n",
    "rho_mass, rho_width, rho_phase, rho_scale = pdg[\"rho\"]\n",
    "\n",
    "rho_m = zfit.Parameter(\"rho_m\", ztf.constant(rho_mass), floating = False) #lower_limit = rho_mass - rho_width, upper_limit = rho_mass + rho_width)\n",
    "rho_w = zfit.Parameter(\"rho_w\", ztf.constant(rho_width), floating = False)\n",
    "rho_p = zfit.Parameter(\"rho_p\", ztf.constant(rho_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "rho_s = zfit.Parameter(\"rho_s\", ztf.constant(rho_scale), lower_limit=rho_scale-np.sqrt(rho_scale), upper_limit=rho_scale+np.sqrt(rho_scale))\n",
    "\n",
    "#omega\n",
    "\n",
    "omega_mass, omega_width, omega_phase, omega_scale = pdg[\"omega\"]\n",
    "\n",
    "omega_m = zfit.Parameter(\"omega_m\", ztf.constant(omega_mass), floating = False)\n",
    "omega_w = zfit.Parameter(\"omega_w\", ztf.constant(omega_width), floating = False)\n",
    "omega_p = zfit.Parameter(\"omega_p\", ztf.constant(omega_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "omega_s = zfit.Parameter(\"omega_s\", ztf.constant(omega_scale), lower_limit=omega_scale-np.sqrt(omega_scale), upper_limit=omega_scale+np.sqrt(omega_scale))\n",
    "\n",
    "\n",
    "#phi\n",
    "\n",
    "phi_mass, phi_width, phi_phase, phi_scale = pdg[\"phi\"]\n",
    "\n",
    "phi_m = zfit.Parameter(\"phi_m\", ztf.constant(phi_mass), floating = False)\n",
    "phi_w = zfit.Parameter(\"phi_w\", ztf.constant(phi_width), floating = False)\n",
    "phi_p = zfit.Parameter(\"phi_p\", ztf.constant(phi_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "phi_s = zfit.Parameter(\"phi_s\", ztf.constant(phi_scale), lower_limit=phi_scale-np.sqrt(phi_scale), upper_limit=phi_scale+np.sqrt(phi_scale))\n",
    "\n",
    "#jpsi\n",
    "\n",
    "jpsi_mass, jpsi_width, jpsi_phase, jpsi_scale = pdg[\"jpsi\"]\n",
    "\n",
    "jpsi_m = zfit.Parameter(\"jpsi_m\", ztf.constant(jpsi_mass), floating = False)\n",
    "jpsi_w = zfit.Parameter(\"jpsi_w\", ztf.constant(jpsi_width), floating = False)\n",
    "jpsi_p = zfit.Parameter(\"jpsi_p\", ztf.constant(jpsi_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "jpsi_s = zfit.Parameter(\"jpsi_s\", ztf.constant(jpsi_scale), floating = False) #, lower_limit=jpsi_scale-np.sqrt(jpsi_scale), upper_limit=jpsi_scale+np.sqrt(jpsi_scale))\n",
    "\n",
    "#psi2s\n",
    "\n",
    "psi2s_mass, psi2s_width, psi2s_phase, psi2s_scale = pdg[\"psi2s\"]\n",
    "\n",
    "psi2s_m = zfit.Parameter(\"psi2s_m\", ztf.constant(psi2s_mass), floating = False)\n",
    "psi2s_w = zfit.Parameter(\"psi2s_w\", ztf.constant(psi2s_width), floating = False)\n",
    "psi2s_p = zfit.Parameter(\"psi2s_p\", ztf.constant(psi2s_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "psi2s_s = zfit.Parameter(\"psi2s_s\", ztf.constant(psi2s_scale), floating = False) #, lower_limit=psi2s_scale-np.sqrt(psi2s_scale), upper_limit=psi2s_scale+np.sqrt(psi2s_scale))\n",
    "\n",
    "#psi(3770)\n",
    "\n",
    "p3770_mass, p3770_width, p3770_phase, p3770_scale = pdg[\"p3770\"]\n",
    "\n",
    "p3770_m = zfit.Parameter(\"p3770_m\", ztf.constant(p3770_mass), floating = False)\n",
    "p3770_w = zfit.Parameter(\"p3770_w\", ztf.constant(p3770_width), floating = False)\n",
    "p3770_p = zfit.Parameter(\"p3770_p\", ztf.constant(p3770_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "p3770_s = zfit.Parameter(\"p3770_s\", ztf.constant(p3770_scale), lower_limit=p3770_scale-np.sqrt(p3770_scale), upper_limit=p3770_scale+np.sqrt(p3770_scale))\n",
    "\n",
    "#psi(4040)\n",
    "\n",
    "p4040_mass, p4040_width, p4040_phase, p4040_scale = pdg[\"p4040\"]\n",
    "\n",
    "p4040_m = zfit.Parameter(\"p4040_m\", ztf.constant(p4040_mass), floating = False)\n",
    "p4040_w = zfit.Parameter(\"p4040_w\", ztf.constant(p4040_width), floating = False)\n",
    "p4040_p = zfit.Parameter(\"p4040_p\", ztf.constant(p4040_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "p4040_s = zfit.Parameter(\"p4040_s\", ztf.constant(p4040_scale), lower_limit=p4040_scale-np.sqrt(p4040_scale), upper_limit=p4040_scale+np.sqrt(p4040_scale))\n",
    "\n",
    "#psi(4160)\n",
    "\n",
    "p4160_mass, p4160_width, p4160_phase, p4160_scale = pdg[\"p4160\"]\n",
    "\n",
    "p4160_m = zfit.Parameter(\"p4160_m\", ztf.constant(p4160_mass), floating = False)\n",
    "p4160_w = zfit.Parameter(\"p4160_w\", ztf.constant(p4160_width), floating = False)\n",
    "p4160_p = zfit.Parameter(\"p4160_p\", ztf.constant(p4160_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "p4160_s = zfit.Parameter(\"p4160_s\", ztf.constant(p4160_scale), lower_limit=p4160_scale-np.sqrt(p4160_scale), upper_limit=p4160_scale+np.sqrt(p4160_scale))\n",
    "\n",
    "#psi(4415)\n",
    "\n",
    "p4415_mass, p4415_width, p4415_phase, p4415_scale = pdg[\"p4415\"]\n",
    "\n",
    "p4415_m = zfit.Parameter(\"p4415_m\", ztf.constant(p4415_mass), floating = False)\n",
    "p4415_w = zfit.Parameter(\"p4415_w\", ztf.constant(p4415_width), floating = False)\n",
    "p4415_p = zfit.Parameter(\"p4415_p\", ztf.constant(p4415_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)\n",
    "p4415_s = zfit.Parameter(\"p4415_s\", ztf.constant(p4415_scale), lower_limit=p4415_scale-np.sqrt(p4415_scale), upper_limit=p4415_scale+np.sqrt(p4415_scale))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dynamic generation of 2 particle contribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "m_c = 1300\n",
    "\n",
    "Dbar_phase = 0.0\n",
    "DDstar_phase = 0.0\n",
    "Dstar_mass = pdg['Dst_M']\n",
    "Dbar_mass = pdg['D0_M']\n",
    "D_mass = pdg['D0_M']\n",
    "\n",
    "Dbar_s = zfit.Parameter(\"Dbar_s\", ztf.constant(0.0), lower_limit=-0.3, upper_limit=0.3)\n",
    "Dbar_m = zfit.Parameter(\"Dbar_m\", ztf.constant(Dbar_mass), floating = False)\n",
    "Dbar_p = zfit.Parameter(\"Dbar_p\", ztf.constant(Dbar_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)#, floating = False)\n",
    "DDstar_s = zfit.Parameter(\"DDstar_s\", ztf.constant(0.0), lower_limit=-0.3, upper_limit=0.3)#, floating = False)\n",
    "Dstar_m = zfit.Parameter(\"Dstar_m\", ztf.constant(Dstar_mass), floating = False)\n",
    "D_m = zfit.Parameter(\"D_m\", ztf.constant(D_mass), floating = False)\n",
    "DDstar_p = zfit.Parameter(\"DDstar_p\", ztf.constant(DDstar_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)#, floating = False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tau parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "tau_m = zfit.Parameter(\"tau_m\", ztf.constant(pdg['tau_M']), floating = False)\n",
    "Ctt = zfit.Parameter(\"Ctt\", ztf.constant(0.0), lower_limit=-1.5, upper_limit=1.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "x_min = 2*pdg['muon_M']\n",
    "x_max = (pdg[\"Bplus_M\"]-pdg[\"Ks_M\"]-0.1)\n",
    "\n",
    "# # Full spectrum\n",
    "\n",
    "obs_toy = zfit.Space('q', limits = (x_min, x_max))\n",
    "\n",
    "# Jpsi and Psi2s cut out\n",
    "\n",
    "obs1 = zfit.Space('q', limits = (x_min, jpsi_mass - 60.))\n",
    "obs2 = zfit.Space('q', limits = (jpsi_mass + 70., psi2s_mass - 50.))\n",
    "obs3 = zfit.Space('q', limits = (psi2s_mass + 50., x_max))\n",
    "\n",
    "obs_fit = obs1 + obs2 + obs3\n",
    "\n",
    "# with open(r\"./data/slim_points/slim_points_toy_0_range({0}-{1}).pkl\".format(int(x_min), int(x_max)), \"rb\") as input_file:\n",
    "#     part_set = pkl.load(input_file)\n",
    "\n",
    "# x_part = part_set['x_part']\n",
    "\n",
    "# x_part = x_part.astype('float64')\n",
    "\n",
    "# data = zfit.data.Data.from_numpy(array=x_part, obs=obs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "total_f = total_pdf_cut(obs=obs_toy, jpsi_mass = jpsi_m, jpsi_scale = jpsi_s, jpsi_phase = jpsi_p, jpsi_width = jpsi_w,\n",
    "                    psi2s_mass = psi2s_m, psi2s_scale = psi2s_s, psi2s_phase = psi2s_p, psi2s_width = psi2s_w,\n",
    "                    p3770_mass = p3770_m, p3770_scale = p3770_s, p3770_phase = p3770_p, p3770_width = p3770_w,\n",
    "                    p4040_mass = p4040_m, p4040_scale = p4040_s, p4040_phase = p4040_p, p4040_width = p4040_w,\n",
    "                    p4160_mass = p4160_m, p4160_scale = p4160_s, p4160_phase = p4160_p, p4160_width = p4160_w,\n",
    "                    p4415_mass = p4415_m, p4415_scale = p4415_s, p4415_phase = p4415_p, p4415_width = p4415_w,\n",
    "                    rho_mass = rho_m, rho_scale = rho_s, rho_phase = rho_p, rho_width = rho_w,\n",
    "                    omega_mass = omega_m, omega_scale = omega_s, omega_phase = omega_p, omega_width = omega_w,\n",
    "                    phi_mass = phi_m, phi_scale = phi_s, phi_phase = phi_p, phi_width = phi_w,\n",
    "                    Dstar_mass = Dstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p, D_mass = D_m,\n",
    "                    Dbar_mass = Dbar_m, Dbar_scale = Dbar_s, Dbar_phase = Dbar_p,\n",
    "                    tau_mass = tau_m, C_tt = Ctt, b0_0 = b0_0, b0_1 = b0_1, b0_2 = b0_2,\n",
    "                    bplus_0 = bplus_0, bplus_1 = bplus_1, bplus_2 = bplus_2,\n",
    "                    bT_0 = bT_0, bT_1 = bT_1, bT_2 = bT_2)\n",
    "\n",
    "total_f_fit = total_pdf_full(obs=obs_fit, jpsi_mass = jpsi_m, jpsi_scale = jpsi_s, jpsi_phase = jpsi_p, jpsi_width = jpsi_w,\n",
    "                    psi2s_mass = psi2s_m, psi2s_scale = psi2s_s, psi2s_phase = psi2s_p, psi2s_width = psi2s_w,\n",
    "                    p3770_mass = p3770_m, p3770_scale = p3770_s, p3770_phase = p3770_p, p3770_width = p3770_w,\n",
    "                    p4040_mass = p4040_m, p4040_scale = p4040_s, p4040_phase = p4040_p, p4040_width = p4040_w,\n",
    "                    p4160_mass = p4160_m, p4160_scale = p4160_s, p4160_phase = p4160_p, p4160_width = p4160_w,\n",
    "                    p4415_mass = p4415_m, p4415_scale = p4415_s, p4415_phase = p4415_p, p4415_width = p4415_w,\n",
    "                    rho_mass = rho_m, rho_scale = rho_s, rho_phase = rho_p, rho_width = rho_w,\n",
    "                    omega_mass = omega_m, omega_scale = omega_s, omega_phase = omega_p, omega_width = omega_w,\n",
    "                    phi_mass = phi_m, phi_scale = phi_s, phi_phase = phi_p, phi_width = phi_w,\n",
    "                    Dstar_mass = Dstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p, D_mass = D_m,\n",
    "                    Dbar_mass = Dbar_m, Dbar_scale = Dbar_s, Dbar_phase = Dbar_p,\n",
    "                    tau_mass = tau_m, C_tt = Ctt, b0_0 = b0_0, b0_1 = b0_1, b0_2 = b0_2,\n",
    "                    bplus_0 = bplus_0, bplus_1 = bplus_1, bplus_2 = bplus_2,\n",
    "                    bT_0 = bT_0, bT_1 = bT_1, bT_2 = bT_2)\n",
    "                   \n",
    "# print(total_pdf.obs)\n",
    "\n",
    "# print(calcs_test)\n",
    "\n",
    "# for param in total_f.get_dependents():\n",
    "#     print(zfit.run(param))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor 'normalization/hook_integrate/hook_numeric_integrate/mul_1:0' shape=() dtype=float64>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "total_f_fit.normalization(obs_toy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Test if graphs actually work and compute values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# def total_test_tf(xq):\n",
    "\n",
    "#     def jpsi_res(q):\n",
    "#         return resonance(q, jpsi_m, jpsi_s, jpsi_p, jpsi_w)\n",
    "\n",
    "#     def psi2s_res(q):\n",
    "#         return resonance(q, psi2s_m, psi2s_s, psi2s_p, psi2s_w)\n",
    "\n",
    "#     def cusp(q):\n",
    "#         return bifur_gauss(q, cusp_m, sig_L, sig_R, cusp_s)\n",
    "\n",
    "#     funcs = jpsi_res(xq) + psi2s_res(xq) + cusp(xq)\n",
    "\n",
    "#     vec_f = vec(xq, funcs)\n",
    "\n",
    "#     axiv_nr = axiv_nonres(xq)\n",
    "\n",
    "#     tot = vec_f + axiv_nr\n",
    "    \n",
    "#     return tot\n",
    "\n",
    "# def jpsi_res(q):\n",
    "#     return resonance(q, jpsi_m, jpsi_s, jpsi_p, jpsi_w)\n",
    "\n",
    "# calcs = zfit.run(total_test_tf(x_part))\n",
    "\n",
    "test_q = np.linspace(x_min, x_max, int(2e6))\n",
    "\n",
    "probs = total_f_fit.pdf(test_q, norm_range=False)\n",
    "\n",
    "calcs_test = zfit.run(probs)\n",
    "# res_y = zfit.run(jpsi_res(test_q))\n",
    "# b0 = [b0_0, b0_1, b0_2]\n",
    "# bplus = [bplus_0, bplus_1, bplus_2]\n",
    "# bT = [bT_0, bT_1, bT_2]\n",
    "# f0_y = zfit.run(tf.math.real(formfactor(test_q,\"0\", b0, bplus, bT)))\n",
    "# fplus_y = zfit.run(tf.math.real(formfactor(test_q,\"+\", b0, bplus, bT)))\n",
    "# fT_y = zfit.run(tf.math.real(formfactor(test_q,\"T\", b0, bplus, bT)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel_launcher.py:12: UserWarning: Creating legend with loc=\"best\" can be slow with large amounts of data.\n",
      "  if sys.path[0] == '':\n",
      "C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\pylabtools.py:128: UserWarning: Creating legend with loc=\"best\" can be slow with large amounts of data.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaUAAAD4CAYAAABMtfkzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO29eXzc1XX3/z4zo9Eu2dq825JtGSMbsxmbHWM2k5CYJqRxElJSaMmTkrZJnjxNaJukTeH3lLZPaNKsvGISsmEoIcEhBIdgUwIYY4MNeMWyvMmbVkvWrpm5vz++d0bj0ayypBlJ5/16+eWZO/d77p2vre9H59xzzxVjDIqiKIqSCbjSPQFFURRFCaKipCiKomQMKkqKoihKxqCipCiKomQMKkqKoihKxuBJ9wQyjbKyMlNZWZnuaShKxrP35BkKsj3MnJx7Vvvprj6OtnazYEoh2Z70/t6781gbZYXZTC3KoaWzj2Onuzl/WhEel6R1XuORN998s8kYU36udlSUIqisrGTbtm3pnoaiZDzLHvwDKxdW8K8fXnJW+zM7jvG363bwqy9cx/yKgjTNzmHe3z/H/7puLv/nloX8YssR/v5X77Lh729gSlFOWuc1HhGRw8NhR8N3iqIMiYABVxSPw23bAhmwB9IYg0uc+QSn6g+kf15KbFSUFEUZEgFjiBYFc1sRyISHf8BAcIquDBJLJTYqSoqiDImAMSEBCifoKaVblILVaiTkKYltT9uUlCTQNSVFUYaEP2BCD/xwMkeUnL+DU0xX+K6/v5/6+np6enpGddyRIicnh5kzZ5KVlTUi9lWUFEUZEsYMCFA4wTCZP80uSXB0sQG8dK111dfXU1hYSGVlZVQRH0sYY2hubqa+vp6qqqoRGUPDd4qiDAl/IP6aUiDtnpIzfnCOQUEY7Wn19PRQWlo65gUJnHtYWlo6ol6fipKiKEMiYEzU7DtPhoTvAjHCd+lIdBgPghRkpL+LipKiKEMiEJZuHY4rQ0TJcHaiQ8iD00yHjEZFSVGUIeEPJMi+S/eaUoSnJBmUqp6JvPTSS9x2220A9Pb2cuONN3LRRRfxxBNPjOo8NNFBUZQhEWvzrCtDHv4hUYpIdFBHKTHbt2+nv7+fHTt2jPrYSXlKIrJKRPaJSK2IfDnK59ki8oT9fIuIVIZ9dr9t3ycitySyKSJV1sZ+a9MbbwwRKRWRTSLSISLfjjH/9SKyM7lboihKIoJJDNESHTwZskk1GL5zZcCaUro5dOgQCxcu5K677mLJkiXccccddHV18fzzz7Nw4UKuvvpqnn76aQAaGhq488472bFjBxdddBEHDhwY1bkm9JRExA18B7gJqAe2ish6Y8zusG73AK3GmPkisgZ4CPioiNQAa4BFwHTgDyKywF4Ty+ZDwMPGmHUi8n1r+3uxxgB6gK8Ai+2fyPl/COhI6a4oihKX4IM9XvjO58+0RIf0e3D//Jtd7D7ePqw2a6YX8bUPLErYb9++faxdu5arrrqKu+++m2984xv84Ac/YOPGjcyfP5+PfvSjAFRUVPDDH/6Q//iP/+DZZ58d1rkmQzKe0jKg1hhTZ4zpA9YBqyP6rAYes6+fAm4QJ4C7GlhnjOk1xhwEaq29qDbtNSutDazN2+ONYYzpNMa8giNOZyEiBcAXgAeS+J6KoiRJcL0oXvgu3R7JQEq4regQ8uDSNqW0MmvWLK666ioA7rzzTrZt20ZVVRXV1dWICHfeeWeaZ+iQzJrSDOBo2Pt6YHmsPsYYn4i0AaW2/fWIa2fY19FslgKnjTG+KP1jjdEUZ+7/Avw/oCveFxSRe4F7AWbPnh2vq6IoDKzLRMu+G6joMJozGkyk+AT106RRLJPxaEaKyFTutra2jExVT8ZTijbryH/VWH2Gqz3ZeQxMSOQiYL4x5lex+oSMGPOIMWapMWZpefk5HweiKOMef5w1pUzJviMUvju79l26EzDSxZEjR9i8eTMAjz/+ODfeeCMHDx4MrRk9/vjj6ZxeiGREqR6YFfZ+JnA8Vh8R8QDFQEuca2O1NwGTrI3IsWKNEYsrgEtF5BDwCrBARF6K+00VRUkKn32wRyszNOAppddVGpzoMLHDd+effz6PPfYYS5YsoaWlhc9//vM88sgjvP/97+fqq69mzpw56Z4ikFz4bitQLSJVwDGcxIWPR/RZD9wFbAbuADYaY4yIrAd+ISLfwEl0qAbewPF6Btm012yyNtZZm8/EGyPWpI0x38NJkMBm6j1rjFmRxPdVFCUBwey7aCe4DhxdMapTGkQo0cG+n8jZdwAul4vvf//7Z7WtWrWKvXv3Duq7YsUKVqxYMUozO5uEomTXbz4LbADcwKPGmF0i8nVgmzFmPbAW+KmI1OJ4L2vstbtE5ElgN+AD7jPG+AGi2bRDfglYJyIPANutbWKNYW0dAooAr4jcDtwckR2oKMowEs9Tctn4S8bUvnOltyCrkhpJbZ41xjwHPBfR9tWw1z3AR2Jc+yDwYDI2bXsdTnZeZHu8MSoTzP8QUdLFFUUZGqGUcNfgFQCPbUv3mlKkp5SugqyZQGVlJTt3jo2tmlpmSFGUlBnwlAZ/FtQpX7o9pYFMByAsfJeGeaUz42+4GenvoqKkKErKBAKxPaVMOboiqElBMUpX+C4nJ4fm5uZxIUzB85RycnJGbAytfacoSsrE85Qy5eTZgfDd2Snhoz2tmTNnUl9fT2Nj4+gOPEIET54dKVSUFEVJGX88TylDEgoiU8KD+0RHWyyzsrJG7JTW8YiG7xRFSZmQKMWrfZcpnlJE+G48hNHGMypKiqKkjD9eSniGVE4Iik9k+C7dWYFKfFSUFEVJmXiiFArfpV2UnL8jq4SnO/9CiY+KkqIoKRP0NuJWdEj3mtKg2nfBdlWlTEZFSVGUlAnWtYt6dIVLEMmA8F2M2nfpnpcSHxUlRVFSJljXLlqiQ7A93Q//WIkOqkmZjYqSoigp47OeUrQ1JXC8pfSH785OdAiKU7pT1ZX4qCgpipIywVMpYomSxyXpT3Swfw9KdFBXKaNRUVIUJWUSeUpukbTvUwp5ShJZJTxtU1KSQEVJUZSUGagSHjt8l26PZODIdufvUEUHDd9lNCpKiqKkTDDRIVpKODhile6Hf2Ttu2BShqaEZzYqSoqipEwoJTxW9p0r/dl3MY9D1/hdRqOipChKyoQ8JXd0UcpyCf3+NHtKdo6RiQ5pnpaSABUlRVFSxpfAU/K4XfiCypUmzED+HTBw+KCG7zIbFSVFUVImEKfMEECWW+hPd/guItFhoPadilImo6KkKErK+Pzxs++y3C76fWn2lAbVvguWGUrXjJRkSEqURGSViOwTkVoR+XKUz7NF5An7+RYRqQz77H7bvk9EbklkU0SqrI391qY33hgiUioim0SkQ0S+HWYnT0R+KyJ7RWSXiPxr6rdHUZRoJEoJ97gzYJ9SZKKDfdqpp5TZJBQlEXED3wFuBWqAj4lITUS3e4BWY8x84GHgIXttDbAGWASsAr4rIu4ENh8CHjbGVAOt1nbMMYAe4CvAF6NM/z+MMQuBi4GrROTWRN9XUZTE+OIcXQHgcbnoT7NLEln7zmNVyaeZDhlNMp7SMqDWGFNnjOkD1gGrI/qsBh6zr58CbhDHZ14NrDPG9BpjDgK11l5Um/aaldYG1ubt8cYwxnQaY17BEacQxpguY8wm+7oPeAsYuYPlFWUCEUggSlluSfvDP7L2nTtUvVzjd5lMMqI0Azga9r7etkXtY4zxAW1AaZxrY7WXAqetjcixYo2REBGZBHwAeDHG5/eKyDYR2dbY2JiMSUWZ0PjiHIcOjlfiS/PDP7L2HTiJGekOKyrxSUaUov2vi/xXjdVnuNqTnccgRMQDPA58yxhTF62PMeYRY8xSY8zS8vLyRCYVZcITOnk2xj4ljzv9+5Qia99BZmzqVeKTjCjVA7PC3s8EjsfqY0WgGGiJc22s9iZgkrUROVasMRLxCLDfGPOfSfRVFCUJ/Ak8pSx3BnhKoTJDAzgenIpSJpOMKG0Fqm1WnBcncWF9RJ/1wF329R3ARuP8mrIeWGMz56qAauCNWDbtNZusDazNZxKMERMReQBHvD6XxPdUFCVJ/Imy71zpX1PyR1n3Uk8p8/Ek6mCM8YnIZ4ENgBt41BizS0S+DmwzxqwH1gI/FZFaHO9ljb12l4g8CewGfMB9xhg/QDSbdsgvAeusoGy3tok1hrV1CCgCvCJyO3Az0A78A7AXeMu68N82xvww9dukKEo4/mT2KWVI9l141QlnTUkTHTKZhKIEYIx5Dnguou2rYa97gI/EuPZB4MFkbNr2OpzsvMj2eGNUxph69J8YRVHOiZCnFLPMUPoTCoL7kcJ1Uz2lzEcrOiiKkjL+gEHEOTcpGh6XK+3hu2gbfDMhrKjER0VJUZSU8QdMTC8JbO27NIfvgh7RWdl3GeDBKfFRUVIUJWV8ARPz2ArIjPBdMA3qbE9Js+8yHRUlRVFSpt8fIMsd+/GRCWWGgp6SK2LzrFZ0yGxUlBRFSZlEouT1ZM6akiti82y656XER0VJUZSU8flNzLOUIDNSr6OJkset2XeZjoqSEpOTbT387PXD6Z6GkoH0+0388J3bRb/fpPWU10CUNSW3rillPEntU1ImJnf/eCu7T7RzU80UphTlpHs6SgbhCwTiJjpkuYIH6sVPiBhJYq8pqShlMuopKTFp7eoD0B9iZRAJw3fWi0qnVxIK30WUGUp3WFGJj4qSoigpkyjRIct6R+nMwIu6pqSeUsajoqQoSsokTgl3hCCdmW5Bh8gdmX1nRckYw7HT3emYmhIHFSVFUVIm8eZZ59HSn8ZQmT90ntJAW3iZoWffOcFV/7qRNw+3pmN6SgxUlJSEaLBDiaTfHyDLlTh8l05PyUSpfReefbdpXwMABxo6Rn9ySkxUlBRFSRmfP4GnZAUrnaIUXM4KX1PKcg9UdMj2uIH0enPKYFSUlITo+R9KJP0BEwrRRSMoWOl84Mc6uiLoKXnC0taVzEFFSUmI/sgqkfj8gdBepGh4g2tKacy+M1FSwsOz74JhvT6fekqZhIqSoigpkyh85/VYUfKlM3wXrfbdQE2+oKfU0+8f/ckpMVFRUhKi4TslkmQKsgL0+tL3wA+VGYqxTyn4eVefilImoaKkJETDd0ok/YH4ohRMIuhNY2gsuKYkYdMMP+Svz++IUbd6ShmFipKiKCmTqMxQtvWU0rleEzoOfZCn5MwpGFrU8F1mkZQoicgqEdknIrUi8uUon2eLyBP28y0iUhn22f22fZ+I3JLIpohUWRv7rU1vvDFEpFRENolIh4h8O2Jel4rIu/aab4nEOb9ZiYneNCWSfn/87LtMCN9FSwkPz77rsx16+jXRIZNIKEoi4ga+A9wK1AAfE5GaiG73AK3GmPnAw8BD9toaYA2wCFgFfFdE3AlsPgQ8bIypBlqt7ZhjAD3AV4AvRpn+94B7gWr7Z1Wi76sMRsN3SiS+QCC0QTYa2SFRSr+nFL7HN7yiQ1CUunVNKaNIxlNaBtQaY+qMMX3AOmB1RJ/VwGP29VPADdYrWQ2sM8b0GmMOArXWXlSb9pqV1gbW5u3xxjDGdBpjXsERpxAiMg0oMsZsNk5u6E/CbClJoB6SEgsnfJeMp5RGUYqSfedxu0JVwoOhRV1TyiySEaUZwNGw9/W2LWofY4wPaANK41wbq70UOG1tRI4Va4x4865PMG8AROReEdkmItsaGxvjmJxYqIekxMLJvovnKWVCooPzd/iaktcePhgImJAo6ZpSZpGMKEX7nxf5vIrVZ7jak51HMnMa3GjMI8aYpcaYpeXl5XFMTkzUY1IiSZQSnp1lPaU0PvCjFWQNenB9/kBoY6+KUmaRjCjVA7PC3s8EjsfqIyIeoBhoiXNtrPYmYJK1ETlWrDHizXtmgnkrSaAekxJOIGAIGOJvnnUPPPzThTEGl0B4flN2mCgNeEqa6JBJJCNKW4FqmxXnxUlcWB/RZz1wl319B7DRruOsB9bYzLkqnGSDN2LZtNdssjawNp9JMEZUjDEngDMicrldq/qzMFuKogyRYD27+PuUgp5SGo+uCJiz1pNgYM79vgFPSdeUMgtPog7GGJ+IfBbYALiBR40xu0Tk68A2Y8x6YC3wUxGpxfFe1thrd4nIk8BuwAfcZ4zxA0SzaYf8ErBORB4AtlvbxBrD2joEFAFeEbkduNkYsxv4DPBjIBf4nf2jpIiG75Rw+iPK9ERDRPC6XWlfU3JFzPHs8J3zPVSUMouEogRgjHkOeC6i7athr3uAj8S49kHgwWRs2vY6nOy8yPZ4Y1TGaN8GLI72maIoQyMY9gp6Q7HI9rjSvnk2UjdDYUVfIFRuSNeUMgut6KAkRNeUlHBCopTljtsvO8uV3tp3AXNW5h2EeUq+QCgMqaKUWagoKQmJs3SnTECCQuONs6YU/Dyd4Tu/ib2m1OcPhDbR9vsNvlFKyGjr6h+VccYyKkqKoqREb8hTShC+y3KnNXxnzNnp4HB2Tb5wIeoZhXn+7t0TXPj137OlrnnExxrLqCgpCVFHSQknKDTJeUrpC435olQyDw/f+cJOnB2NEN4f9jQA8OoBFaV4qCgpMVExUqIRFJrk1pTSmxLudsVICfcbfAFDjvX2RqP+3fHT3QAcbOoc8bHGMipKiqKkRG+SnlK6s++iHa8xkBLup98foCA7CxgdT6m9x1lPOtnWPeJjjWVUlJSYBOPx6jEp4SS7puT1ZICnFFF1Ijwl3Oc3FOU4u2JGo6pDW7cjSsdP9yToObFRUVJiomKkRCPZNaVsjzut6da+wOBK5uHVy/0BQ4EVpdHYQNtuRelUe49mtMZBRUlJiNGdSkoYQe8nJ4GnlOt1p7Vagi8QGLSm5A1bU+oPBCjIDnpKIzvPQMBwptdHvteNL2Bo7/YlvmiCoqKkKEpKDHhK8RMd8rLcaT1AL96aUk+/H2OgcJQ8pTO9PoyBqvJ8AJo7e0d0vLGMipKSEI00KOEMZN/Ff3zked109qbPI/AHzKBK5kFR6upz5jVaiQ7B0F1VWQEALZ19IzreWEZFSUmIapISTrK173K9njSH7wzuGGtKXdaDK8wZnfBdMMmhqjQPgGYVpZioKCmKkhKhlPAEopTvdTtrN2k6U8kfGBy+C56WGylKIx1mDHlKNnynnlJsVJSUhGimkBJOstl3uV5nzakrTetKsRIdROBMTzB8Zz2lEU5db4sI3zV36JpSLFSUFEVJiV6fH7dL8CQQpTyv88APrt+MNtESHUSE3Cw3HXatKz97lDwlu3G2rMBLQbZHw3dxUFFSEqJ+khJOny+QcD0JnEQHSKenZKIKZ57XzRkrEl63i2yPi54RrtEX9JSKc7Moyfdq+C4OKkpKQjR6p4TT6wskXE+CAVFKV1p4tDUlcMKKwTUej1vI9brpGfE1JR8ugXyvR0UpASpKiqKkRE+/P0lPKRi+S5+nFLmmBJCX5aHdril53C7ystx0jvAc27r7KcrNwuUSygq8NHWoKMVCRUlJAnWVlAG6+vwhwYlHMNGhM01rSv5AILGn5BKKcrNC4byRor2nn6IcZ09UaX62JjrEQUVJiYmWF1Ki0dPvJzfBsRWQ/vCdzx/DU/K6Q4kHHpdQlJMVWvMZKdq6+ynOdUSpvDCb5s4+AgH9+YpGUqIkIqtEZJ+I1IrIl6N8ni0iT9jPt4hIZdhn99v2fSJySyKbIlJlbey3Nr3nMMbnRWSXiOwUkcdFJCe126OArikpZ+N4SolFKT8DwnfRPKU8rztUFTw7y01RbhZtI1yLzgnfOfejrMCLP2Bo7dIQXjQSipKIuIHvALcCNcDHRKQmots9QKsxZj7wMPCQvbYGWAMsAlYB3xURdwKbDwEPG2OqgVZreyhjzAD+BlhqjFkMuG0/JUkE5wdaNUkJp6vPHwrNxWNgn1K6wnfRs+9yw0KP2R4XRbmeUDhvpGgP85TKCrMBdF0pBsl4SsuAWmNMnTGmD1gHrI7osxp4zL5+CrhBRMS2rzPG9BpjDgK11l5Um/aaldYG1ubtQxwDwAPkiogHyAOOJ/F9FYuG75RoJBu+y8+2a0q96ds8G9VTCpt7tsdFUU7WiItSW7cvtKZUVhAUJV1XikYyojQDOBr2vt62Re1jjPEBbUBpnGtjtZcCp62NyLFSGsMYcwz4D+AIcAJoM8b8PtoXFJF7RWSbiGxrbGyMeSMmKhq+U8JJNnyXm+XG7RI6ekf2gR+LaMehA2d5edkeN8W5WZzp9eEfwTWe9p4wT0lFKS7JiNLgf9XBEZ1YfYarPeUxRGQyjhdVBUwH8kXkzih9McY8YoxZaoxZWl5eHq2LoiiWZMN3IkJRjidtZwfFWlM6S5SyXBRZsRipDLyefj99vkBonHIrSo1nVJSikYwo1QOzwt7PZHAYLNTHhsqKgZY418ZqbwImWRuRY6U6xo3AQWNMozGmH3gauDKJ76tEoGE8JRwnfJc4JRygKDcrlOk22vT7AoOqhINTKDZItscV8mBGSjyDmX1BUSrK9eB1u3RNKQbJiNJWoNpmxXlxkgXWR/RZD9xlX98BbDROFc/1wBqbOVcFVANvxLJpr9lkbWBtPjPEMY4Al4tInl17ugHYk9xtUcLR8J0SxBhDV58vqfAdMCrrNbHo95uolSeK87yh19keN0W2Uvjp7pERiaAoTbKiJCKUFng1fBeDhL/uGGN8IvJZYANOBtujxphdIvJ1YJsxZj2wFvipiNTieC9r7LW7RORJYDfgA+4zxvgBotm0Q34JWCciDwDbrW2GMMYWEXkKeMu2bwceGeqNUhTFKTEUMCQVvgPHKxjpPUDRMMbQ549eDqkkXJSyXGHZcCMjEuF174KUFWSrKMUgKR/cGPMc8FxE21fDXvcAH4lx7YPAg8nYtO11DGTPhbcPZYyvAV+Ldo2SPOopKUGCG2GTyb4Dx1M61T76D98+f+yDCCfnD4hDvtdDReHIrvG0dUUTJS8NuqYUFa3ooMRExUiJJHiSbLLhu+Lc9ITv4p35VJI/4Cm5XRLKhmsYIfE8HQzf5Q2IUkVhTlrEeiygoqQkRBMdlCDB6gzJh+/Sk+jQF+d03PDwHUBOlpMWPlKeS7Tw3fRJuTR19I74MexjERUlJSHqMSlBgtUZkinIClCU46GnP0DvCJ9XFEkwfBdNlIKeUfAodHDq0Y1Y+M6KUmHOgCjNmJwLwIm2nhEZcyyT3P8sRVEUBo4RD3+gx2NgD5CP7ILkvKvhIF74zuUSfvGXy88K41UUZtNwZmQEoq2rj6Icz1kbeWdMckTp+OluqsryR2TcsYp6SoqiJE1wg2nSomS9g9HOwIsXvgO4cl4ZC6cWhd5PLcrh5Ah5LW3d/RSHrScBzLSe0rHW7hEZcyyjoqQoStIED8cryslK0NNhsvVGRvuk1d4EohTJ7NI8TrT3jMgaT2tXP5Mj1rGmFucgAvWnVZQiUVFSEqJrSkqQjhTDd2UFzsN4tA+164+zphSNOaV5GAP1rV3DPpfGM72h0kJBstwuphTmqKcUBRUlRVGSJrimVJCdnCiF6ryNckmdYPguO8qaUjTmlDrrOoebh1+UGs70UlGUPah95uTcERHBsY6KkpIQTQlXgpzp6SfP6456TlE0guG70faU4mXfRaPSitKhYRYlf8DQ0jnYUwKYW57PgcbOYR1vPKCipCREw3dKkDM9vqS9JHDCVJPzska9pE6iRIdIJudlUZTj4UBjx7DOo7mzl4BxUs4jmVdeQFNHb6jig+KgoqQoStKc6e1Pej0pSGlBNk1n0hO+S1aURISa6UXsOt6e8ljGGL74329z5w+30NF7dqXxYJWI8sKcQdfNrygAoHaYhXCso6KkJEQdJSXImR7fWZtAk6GswEtz5+h6Sj12s26OJ/m9UYunF7P3RDs+G/pLljcPt/LUm/W8UtvE2j8ePOuzoy1OOHBWSe6g6+aVO6J0oEFFKRwVJSUmKkZKJO09vpQ9Jaci9uh6SsEj2JOt0QewaEYRvb5Ayus82w63AnDpnMn89PVDZ1WvCK5RBRMpwplVkofX41JPKQIVJSUhRheVFEt7d/9ZNdySoaLQ2Zg6mv+PulOs0QdwwYxJAGw/0prSWLUNHZQXZvM3N1TT1NHH7949GfrscHMnZQXeqOtwbpcwr7yAPSdSDxmOZ1SUlJgEi6KoJClBmjt6Kc33Ju4YxozJuXT3+2kdxQX9YOHYZGv0Acwrz6e8MJvXDjSnNNb+hg6qKwq4Zn4ZVWX5/GTzodBne06eCYXporFkRjHvHmvTX/zCUFFSYqI/Jko4/f4A7T0+SvIHZ5LFI1jnbTQ3inb1+/B6XGfVm0uEiHDlvFJeO9CctEgYYzjQ0MH8igJcLuHOy+fw1pHT7DzWRq/Pz57j7Vw0e1LM65fMKuZ0Vz9HW3QTbRAVJSUh+kucAtDa5awLleSnFr4L1Xk7PXobRbv7/CmtJwW5an4ZTR297E4ypHaqvZeOXh/VNpPujktnkpvl5iebD7H5QDN9/gCXzSmJef0SGzJ859jplOc6XlFRUhQlKVo7nfBbqp5SUJTqR9NT6vOTl+TpuOHcsLACt0v47Tsnkuq/v+EMAPOsKBXnZnHHpTP55VvH+OozuyjOzeLq6rKY1583tRCv28U79W0pz3W8oqKkJIG6SgqhtO7JKXpKxblZ5HvdHBvF4qNdfb6UkhyClBZkc9X8Mn7zzvGkQni1Np27uqIw1PZ3q87jghnFnGjr5msfqCEnjjh6PS6WzCxmy8GWhGP9YfcpNu49lcS3GNuoKCkJ0fCdAgOeUmmKnpKIMHNyXmjPzmjQ1edPKckhnA9eOJ2jLd1srkuc8LC/oYPi3KxQ4VlwDvP71V9dyc5/voUPXTIzoY0r55Xybv3puCf0tnb28Rc/2cbdP94WOmhxvJKUKInIKhHZJyK1IvLlKJ9ni8gT9vMtIlIZ9tn9tn2fiNySyKaIVFkb+61N7zmMMUlEnhKRvSKyR0SuSO32KIoSpGWInhI41Qv2j+Im0a4+/5A8JYDblkyjJN/Lo68cSti31mbeiZydUCEiZCe5cfeKeWUEDLxRF9tbert+YM1p66HUUtbHGglFSUTcwHeAWx/S3TUAACAASURBVIEa4GMiUhPR7R6g1RgzH3gYeMheWwOsARYBq4Dviog7gc2HgIeNMdVAq7Wd8hj2mm8CzxtjFgIXAnuSvTHKAOooKQBNHX2IMOhsoGSonlLAkZauETmvKBrt3f1Jn/kUSU6Wm08sn82Le09Ra9eMYlFrM+/OhYtnTyLb4+LVA00x+4R7me/Wj++kiGQ8pWVArTGmzhjTB6wDVkf0WQ08Zl8/Bdwgzq8Oq4F1xpheY8xBoNbai2rTXrPS2sDavH0oY4hIEXAtsBbAGNNnjBnf/5rDjIbtlHBOtfdQVpBNVpIVwsOprijEmIE1mJGmvbufSXlDEyWAT11ZSb7Xw0PP74vZp6mjl5bOPhZMKYzZJxlystwsn1vKS/saY65jHW7uItvjorI0b9wnRSTzv2sGcDTsfb1ti9rHGOMD2oDSONfGai8FTlsbkWOlOsZcoBH4kYhsF5EfisjgWh+AiNwrIttEZFtjY2PsOzFBUXFSAE609TCteHBh0WRYMMXxJvYn8DyGi9NDqDwRTmlBNp9ZMY8Xdp9ic4zNtO+ddL7LeVPPTZQAbqqZwsGmzpiifaSli9kleSyaXsyek+O7AkQyohRt91nkYypWn+FqH8oYHuAS4HvGmIuBTmDQehiAMeYRY8xSY8zS8vLyaF0mJMEwue42VwBOtvUwpWhoolRZlo/X7WLPiZEXpX5/gK4+/zmJEsDdV1UxuySPLz/9TtTkgn2nnO9yrp4SwM01UwDYsOtk1M+PtHQxpzSP86cVcrSlO25SxFgnGVGqB2aFvZ8JHI/VR0Q8QDHQEufaWO1NwCRrI3KsoYxRb4zZYtufwhEpJUlUi5RwTrYP3VPKcrtYPKOItw6P/CJ9W7fzwD6X8B04dfP+/Y4lHG7u4oHfDl6O3nH0NBWF2Wdl3g2VKUU5XDRrEht2DU75NsZwpKWLWSV51EwvAmDvKIh7ukhGlLYC1TYrzouTVLA+os964C77+g5go3F+vV4PrLGZc1VANfBGLJv2mk3WBtbmM0MZwxhzEjgqIufZa24AdifxfZUIVJuUrj4fbd39TB2iKAFcMnsy7x5rC511NFIERelcPSWA5XNL+fS1c/nFliM8sfVIqN0Yw5a6FpZVlQzKvBsqty6eyrvH2qiLqBre3NlHV5+fOSV5nD/NEaXxXMQ1oSjZ9ZvPAhtwsteeNMbsEpGvi8gHbbe1QKmI1AJfwIbJjDG7gCdxxOB54D5jjD+WTWvrS8AXrK1SazvlMew1fw38XETeAS4C/r9Ub5CiKE7oDmDqEMN3AJfMmUyvLzDiD9TWTqcc0qQhZAlG4//cch7XVJfxj7/eyYt7HE9mx9HTnGzv4drq4Qv3337xDFwCT7917Kz2I6EzmfKYWpTDpLyscS1KSe0uM8Y8BzwX0fbVsNc9wEdiXPsg8GAyNm17HU52XmT7UMbYASyNdo2SPBrGU4IPxpmT84ZsY2nlZABePdDEhbNiFyk9V07Z016nFKW2yTcWHreLb3/sEj756BY+/dM3ue/6+Ww+0ExBtodVF0wdljHACeFdU13O02/V84WbFuCyxWSPhM5kynNOyJ1WlHRtvrGIVnRQEmI0gDfhOdTkHHxXWTZ0UaoozOGCGcVs3NMwXNOKysn2c/fqIinOy+Lnf7Gcm2qm8M0X97PtcAv/9MFFQ94LFYsPXzqT4209Zx2fcbi5C5GBXwjOn1bEvpNnUj4hd6wwtDociqJMKA41d5HvdVNecG7ex/ULK/j2xv20dPZRkuK5TMnS0N6D1+MaljWlcApzsvjenZdyuLkTr8fFtOLBR5yfKzfXTGFSXhY/2XwoVMj1SEsXU4tyQjX0zp/mnJB7qLmT+RXnnvmXaainpCRGHaUJz8GmTirL8s95Uf/mmikEDPz2ncgE3uHjZHsPU4tyhi0BIZI5pfkjIkjgbKS9c/kcXthzKuSd1jZ2UFU2sMWyxiY77Do+PkN4KkqKoiTkULMjSufK4hnFLJpexLqtRxN3HiKHm7uYVTIyojEa/NkVc8hyuVj7ykH6/U5iyOIZxaHP51cUkOWWUdnzlQ5UlJSEqKM0senu83O0pSvusd6psOayWew63s7WQ4mPa0iV4EmwwzXXdFBRlMOHL53J428cYd3Wo/T5Alwclhji9biYV14wbjPwVJSUhGj23cRm78l2AgYW2Y2b58odl86irCCbb/z+vWGxF07jmV7O9PrGtCgBfPHmBRTkePjKr3dSmO3huvPOTj2vmV6koqQoysQkmH4cXMs4V3K9bu67fh6b65p57t3kTnhNlnePOcVKFw5DPbp0UlqQzU/uXsb7L5jGtz528aCzoWqmFdFwppemjt40zXDkUFFSEqIp4ROb3cfbKcrxhI41Hw4+efkcLphRzD/+eicn2obvRNpth1vJcsuI7oMaLZbMnMR3PnEJ1y+sGPTZeK7soKKkxEHFSHEOmFs0vXhYs9k8bhff+NML6fMF+PMfbeV0V9+w2H35vUaWzJwU9wjy8cDi6cWIwPYj4+80HhUlJSG6pjRxae/pZ/fxdpZVlQy77eophXz3E5dQ19jJh7732jmftVTbcIZdx9u5dfHwVVnIVIrzslg4tYgtBxMf2T7WUFFS4uD8ZqyaNHF581ArAQPLR0CUAK5dUM5P71lGS2cf7/vWH/nGC+8N2Wv6r4215GS5WH1R5HFv45PlVSW8ebh1xAvcjjYqSkocVI4mOq/XNZPlFi6ePXnExlg+t5Tff/5abqqZwrde3M+V/7qRzz+xg+d3ngwVV42HMYYfv3qQZ3Yc5y+vmUt54fDUvMt0Lp9bQk9/gHePja8QnpYZUhKih/xNXF7Yc4rlVaXkekd2jaaiMIfvfPwS/nplOz965RDP7zrJr7Y71bIrS/OYX1HI7JI8pk/KoSDbQ67XTU+/nxNtPWza28Db9W3csLCCv7mhekTnmUksqyoF4PW6Fi6dMzKebDpQUVIUJSq1DR3UNXbyqSsrR23MhVOLeOiOJTzwJ4t583Arbx1p5e2jpznc3MVrB5ro6vMPuqZmWhEP3L6Yjy2bjds1MqWFMpGSfC8LphTwel0z910/P93TGTZUlJSYBKyDpH7SxOR3dg/RTfao7tEky+3i8rmlXD63NNRmjKG9x0dXn4/OXj85WS5K8r2D9vBMJK6eX87Ptxymp98/bjIOdU1JiUkobKeqNOEIBAxPbDvKFXNLR6z4aKqICMW5WUwrzmV+RQEzJ+dNaEECuO68cnp9AV6vGz9ZeCpKSkxUiyYuf6xtor61m48tn53uqShxWF5VQk6Wi/95rzHdUxk2VJSUmAw4SipPE43vv3SAsoJsblk0+qE7JXlystxcPreU/9mnoqRMADTrbmKypa6ZzXXNfGbFPLI942OdYjxz3YJy6po6Q8emj3VUlJSYhDwl1aYJg88f4F9+u5uKwmw+oaG7McF1C5wK4v+zf3x4S0mJkoisEpF9IlIrIl+O8nm2iDxhP98iIpVhn91v2/eJyC2JbIpIlbWx39r0DnUM+5lbRLaLyLPJ3xYFBtaUVJQmDj9+7RA7j7XztQ8sGjfZXOOdqrJ8ZpXk8j/7GtI9lWEhoSiJiBv4DnArUAN8TERqIrrdA7QaY+YDDwMP2WtrgDXAImAV8F0rEvFsPgQ8bIypBlqt7ZTHCJvb3wJ7krsdSjgavptYvFN/mn97fh83LKzgfReM//px4wUR4frzKni1tple3+B9XGONZDylZUCtMabOGNMHrANWR/RZDTxmXz8F3CBOSeHVwDpjTK8x5iBQa+1FtWmvWWltYG3ePsQxEJGZwPuBHyZ3O5RwdJ/SxKGhvYfP/Owtyguz+fePXDisFcGVkWfFeeV09/t54+Dwn+Y72iQjSjOAo2Hv621b1D7GGB/QBpTGuTZWeylw2tqIHCvVMQD+E/g7IG7FQhG5V0S2ici2xsbxEZcdDjTrbmLQ0tnHnWu30NrVx/fuvISSfG+6p6SkyBVzy/B6XGzaO/afX8mIUrRfmSKfVrH6DFd7ymOIyG1AgzHmzSifn93ZmEeMMUuNMUvLy8sTdZ8wDCQ6qDiNV+pbu1jzyGYON3fxw7uWsmTm2D8cbyKS63VzxdxSXnpv7K8rJSNK9cCssPczgeOx+oiIBygGWuJcG6u9CZhkbUSOleoYVwEfFJFDOOHBlSLysyS+r2IxEX8r44s3D7fyJ999jRNtPfzoU5dx5byydE9JOQeuP6+cusZODjd3pnsq50QyorQVqLZZcV6cpIL1EX3WA3fZ13cAG43z6/V6YI3NnKsCqoE3Ytm012yyNrA2nxnKGMaY+40xM40xldb+RmPMnUneFwVUjcYp/oDhWy/u509/sJmcLBdPf+ZKrpyvgjTWWXGec2z6S2N8I23CwlHGGJ+IfBbYALiBR40xu0Tk68A2Y8x6YC3wUxGpxfFe1thrd4nIk8BuwAfcZ4zxA0SzaYf8ErBORB4AtlvbDGUM5dwI2LCdRu/GD+/Wt/GPz+zk7aOnWX3RdP7l9sUU5WSle1rKMFBZlk9VWT6b9jVw1yhWdh9ukqpmaIx5Dnguou2rYa97gI/EuPZB4MFkbNr2Omz2XER7ymOEff4S8FKsz5XoqBaNH1o6+3j4hff42ZbDlOZn8801F02YE1onEivOK+cXW46M6arhE7vErhKXgQQHlaexSlt3P2v/WMfaVw7S3e/nU1dW8vmbFqh3NE65/rwKfvTqITbXNXO9DeeNNVSUlJhoRYexy+muPn6y+TA//GMd7T0+3n/BND53YzXVUwrTPTVlBFlWVUJulpuX9jaoKCnjDxWjscfh5k7WvnKQ/95WT3e/nxsWVvCFmxewaHpxuqemjAI5WW6unFfKpn2N/JMxY3ITtIqSkhDVpszGGMO2w62s/eNBNuw+icclrL5oBvdcXcX504rSPT1llFmxsIIX9zZQ19TJvPKCdE8nZVSUFGWM0tbVz9Pb6/nFliPsb+igODeLv1oxj7uuqKSiKCfd01PSxApbNfylfY0qSsr4RMN4mYMxhreOtPLzLUf47Tsn6PUFuHDWJP7tw0u47cJpE/54cAVmleRRXVHAS/sauOfqqnRPJ2X0f7CSEK2Bl34az/Sy/u3jPLn1KPtOnaEg28Mdl87k48tn63qRMojrF1bw41cP0dnrIz97bD3mx9ZsFWUC0dPv58U9DfzyrXr+571G/AHDkpnF/N8PXcAHL5w+5h42yuixYkE5j7xcx2sHmrmpZmwdaa//q5WEaPhu9AiG55568xjPvnOcMz0+phbl8JfXzOVDl8xggaZ0K0mwtLKEfK+bl/Y1qCgpipI6R5q7+NX2Yzy9vZ7DzV3kZrlZtXgqH7pkBlfOK8PtGnupvUr68HpcXF1dxkv7GjFjLDVcRUlJiDpKI8PJth6efec4v3nnBG8fPY0IXDG3lL9eWc2qxVMp0PCccg5cf14FG3adYn9Dx5jysPV/vZIQPU9p+Gju6OW5nSf5zdvH2XqoBWNg8Ywi7r91IbddOJ0Zk3LTPUVlnHDdeU5q+Ka9DSpKiqIM0Nbdz4ZdjhC9dqAZf8Awv6KAz9+4gNuWTGPuGNxLomQ+04pzWTi1kE37Gvj0dfPSPZ2kUVFSlBGgo9fHi3tO8Zu3T/Dye430+QPMLsnjf103lw9cOJ3zphSOqTi/Mja5bkE5j756kO4+P7nesVE1XEVJUYaJtq5+Xthziud3nuDl/U30+QJMLcrhz66YwwcunM6SmcUqRMqocuX8Mn7wch3bDrdwTXV5uqeTFCpKSlTC15F0SSk2TR29/H7XKX638wSbDzTjCximF+fwieWzuXXxNJbOmYxLM+eUNHFZ5WQ8LuG1A80qSsrYxh8IEyXNvzuLk209bNh1kt/tPMEbB1sIGJhTmsc911Rx6+JpXKgekZIh5Hk9XDx7Eq8daE73VJJGRUmJii+gQhTO0ZYunt/pCNFbR04DUF1RwGevn8+qxdM4f5quESmZyRXzyvj2xv209/SPicMdVZSUqPT7A6HXEzF85w8Y3q4/zaa9Dby4p4HdJ9oBWDS9iC/evIBVi6cxv0Kz5pTM58p5pXzrxf1sqWsZE9UdVJTGOZ9bt51rqsv58KUzU7qu3z/xlOhMTz9/3N/Ei3saeGlfA82dfbgEls4p4f5bF3Lr4mnMLs1L9zQVJSUunj2JbI+L1w40jR9REpFVwDcBN/BDY8y/RnyeDfwEuBRoBj5qjDlkP7sfuAfwA39jjNkQz6aIVAHrgBLgLeCTxpi+VMcQkVm2/1QgADxijPlmqjdorPPrHcf59Y7jKYtSd78/9Ho8e0p1jR1s3NvAxr0NvHGwBV/AUJybxYrzylm5sILrFpQzKc+b7mkqypDJ9ri5rLKEzWNkXSmhKImIG/gOcBNQD2wVkfXGmN1h3e4BWo0x80VkDfAQ8FERqQHWAIuA6cAfRGSBvSaWzYeAh40x60Tk+9b294Ywhg/438aYt0SkEHhTRF6ImPe4xhcWgkuV7j5/4k5jkD5fgG2HWnjRCtHBpk4AFkwp4C+umcvKhRVcMnsSHrcrzTNVlOHjinml/PuGfTR19FJWkJ3u6cQlGU9pGVBrjKkDEJF1wGog/OG+Gvgn+/op4NvirPquBtYZY3qBgyJSa+0RzaaI7AFWAh+3fR6zdr+X6hjGmM3ACQBjzBlre0bEvMc1necgLD3hntJwTCaNNHX08tK+RjbuPcXL7zXR0evD63ZxxbxSPnVlJSsXVjCrRMNyyvjlynmlALxe18xtS6aneTbxSUaUZgBHw97XA8tj9THG+ESkDSi17a9HXDvDvo5msxQ4bYzxRek/lDEAEJFK4GJgS7QvKCL3AvcCzJ49O1qXMUlnry9xpxh09YWH78aWLBlj2H2inY17Gti4r4EdR09jDFQUZnPbkmmsXFjBVfPL9DwiZcJwwYxiCrI9vHZgfIhStDzXyKdUrD6x2qPFRuL1H8oYzkUiBcAvgc8ZY9qj9MUY8wjwCMDSpUvH1hM4Dl19Qxel8DWlsUB3n59Xa5vYuK+BjXsaONneA8CFM4v53A0LuOH8CmqmFelGVmVC4nG7uHTOZLYebEn3VBKSjCjVA7PC3s8EjsfoUy8iHqAYaElwbbT2JmCSiHistxTeP+UxRCQLR5B+box5OonvOq7o6B26sHSFeVmZqtLHTnc7SQp7TvHagWZ6fQHyvW6uqS5n5fkVrDivnIrCnHRPU1EygmVVJfz7hn20dPZRkp+5yTvJiNJWoNpmxR3DSSr4eESf9cBdwGbgDmCjMcaIyHrgFyLyDZwkhGrgDRzvZpBNe80ma2OdtfnMUMaw601rgT3GmG+kemPGA13nEL5r6eobxpkMD/6AYcfRVl7c4yQp7D15BoDZJXl8fPlsVi6sYFlVCdmesVF4UlFGk2VVJQBsPdTCLYumpnk2sUkoSnb95rPABpz07UeNMbtE5OvANmPMepyH/09tkkELjshg+z2Jk1zgA+4zxvgBotm0Q34JWCciDwDbrW1SHUNErgY+CbwrIjusjb83xjw3tFs19uiwouT1pJ5J1toZJkppdJXauvt5+b1GNu519g61dvXjdgmXVU7mH953PtcvrGBeeb5WU1CUBCyZWYzX42LrwTEuSgD2Qf5cRNtXw173AB+Jce2DwIPJ2LTtdQxk6IW3pzSGMeYVoq83TRiCyQrZQxClls7+0OvRrH1njOFAYycb957ixT0NbDvcij9gmJyXxfXnVXD9wgquXVBOcW7ml0tRlEwi2+PmolmTeONQZq8rafrROKa9xxGWnKzUw1kn2rrJyXLR0z/0vU7J0uvz88bBltAm1sPNXQAsnFrIp6+dyw3nV3DRrMm4NUlBUc6J5VUlfPelA3T0+ijI0OzTzJyVMiw0nekFoHAI//kONnUyt6yA3SfaR6SiQ+OZXjbZTLk/7m+ks89PtsfFlfNKQ5tY9WhwRRleLqsswR+o5a3DrVy7IDOPslBRGsc0djiiFEhRVfp8AQ42dXJTzZRQIdJzxRjDruPtNknhFG/XtwEwtSiH1RfP4IaFFVw5r2zMnI6pKGORS+ZMxiVOsoOKkjLq1Ld2A6kXV32n/jS9vgCXVZbw7Dsnhryi1NPvZ3NdM3/Y7awPnWzvQQQumjWJL968gOsXOnuHNElBUUaHgmwPi2cU80YG71dSURrH7LMp075AautCz75zAq/bxeVzndIkqThaLZ19bNzbwB92n+Ll/Y109fnJ87q5prqML5y/gJULKzK+9paijGcuqyzhp68fptfnz8jtEypK45T9p87QYNeUUvGUTrR18+S2o9x6wVQm5SWX4dbW3c/vd53kN++c4NXaJvwBw5SibP7k4hncWDOFK+aWDinZQlGU4WdZVQlrXznIu/VtLK0sSfd0BqGiNA7p6ffz0PN78Xpc3FQzhZffa0zqus5eH3/187cIGMMXbz4vlE8fbU0qEDC8dqCZn285zIt7GujzB5hVksu9187l1sVTuWCGHgmuKJnIZVaIthxsUVFSRpauPh+/2HKEH7xcR+OZXr56Ww0n2rrxJeEp7T7ezuee2E5tQwff/cSlzCrJ47St6tDnGwj/+QOG9W8f47821lLX2MnkvCw+ecUcPnjhdJbMVCFSlEynJN9LdUUBWzN0v5KK0jigrbufx147xI9ePUhrVz9Xzivlvz52MZfPLeXfnt8bd03pQGMHD7/wHs++c4KSfC8/uXs5V1eXAYTizb1WlPafOsPnn9zBzmPtnD+tiIc/eiG3Lp6moTlFGWNcVlXCb3Ycxx8wGbf/T0VpDNN4ppe1rxzkZ68fpqPXxw0LK7hv5XwumT051MfjdtHvNxhjzvJijrZ08a0X9/PLt+rJyXJz3/XzuPeaeRSHrSMFyxP1+QJsqWvmz3+8lZwsN99ccxEfWDJdK24ryhhleVUJv9hyhN3H27lgZnG6p3MWKkpjkENNnax95SBPbDuKzx/gfRdM469WzKdmetGgvllWOHwBQ5ZbaDzTy39t3M/jbxxBRPjzq6r4zIp5UTPi3C7B4xIONXey9pU6pk/K5ed/sZwpRVp5W1HGMsHM2tfrmlWUlKHR6/OzYdcp1r1xhNcONON1u/jwpTP49LXzqCzLj3ldlvV2+v0Bfv76Yf5twz56fQE+etks/mZlNVOL4wtMtsfFr7YfI8/rZu1dS1WQFGUcMKUoh6qyfLYcbOYvr52b7umchYpSBtPr8/NabTPP7zzJht0nOd3Vz8zJuXzx5gX86dJZVCQhEB7rKX3l17v45Vv1rDivnK/eVsPc8oKk5hA8Uv0Ty2czpzS2+CmKMra4fK6zOT7T1pVUlDKMzl4fL+1r5PldJ9m0t4GOXh+F2R5Wnl/BHZfO5Kp5ZSmt5WS5HU/pl2/Vc/dVVXzltvOHlCH351dVpXyNoiiZy+VzS3n8jaPsOdHO4hmZE8JTUcoATnf18Yc9DTy/8yQv72+kzxegNN/LbUumccviqVw5r3TIO6+Dx1Z43S7uf9/ClAXp0U8tpbmjj+laHFVRxhXLqwbWlVSUFBrae9iw+xQbdp5kc10z/oBhenEOn1g+m1WLprK0smRYXOpL50wmyy383w9dEPKaUmHlwinnPAdFUTKPqcU5VJbm8XpdC39xTeasK6kojSLHTnfz/M6T/O7dE7x5pBVjYG55Pp++di6rRqgKQvWUQnZ/fdWQBElRlPHN5XNLee7dzFpXUlEaYbr7/Pxu5wme2HqULbYy78KphXz+xgXcungq1VMKR3wOKkiKokRj+dwS1m3NrHUlFaURorWzjx+9epAfv3aI9h4flaV5fPHmBdy2ZHrcFG5FUZTRIrhf6dXaJhWl8Yo/YPjZ64f59w376Oj1sWrRVD51VSXLq0q0LpyiKBnFtOJcFk0v4ve7T/Hp6+alezoAJBXXEZFVIrJPRGpF5MtRPs8WkSfs51tEpDLss/tt+z4RuSWRTRGpsjb2W5ve4R5jpDjd1cddj77B19bv4uLZk/j956/l+5+8lMvnlqogKYqSkdxcM5W3jrTScKYn3VMBkhAlEXED3wFuBWqAj4lITUS3e4BWY8x84GHgIXttDbAGWASsAr4rIu4ENh8CHjbGVAOt1vZwjzHstHb28ZHvb+aNgy089OEL+Mndy1gwCutFiqIo58LNi6ZgDPzm7RPpngqQnKe0DKg1xtQZY/qAdcDqiD6rgcfs66eAG8RxDVYD64wxvcaYg0CttRfVpr1mpbWBtXn7cI6R3G1JjT5fgE/9eCuHW7r48d2X8dHLZqtnpCjKmGDh1EKWVZbw7Y37OX66O93TSWpNaQZwNOx9PbA8Vh9jjE9E2oBS2/56xLUz7OtoNkuB08YYX5T+wzXGIETkXuBe+7ZDRJqBpmh9E3HVg0O5KmMpY4j3YRyi92IAvRcO4+4+zPjakC8tA+YMxxySEaVov/JHnhoXq0+s9mgeWrz+wznG4EZjHgEeCb4XkW3GmKXR+k4k9D4MoPdiAL0XDnofBrD3onI4bCUTvqsHZoW9nwkcj9VHRDxAMdAS59pY7U3AJGsjcqzhGkNRFEXJUJIRpa1Atc2K8+IkFayP6LMeuMu+vgPYaIwxtn2NzZyrAqqBN2LZtNdssjawNp8ZzjGSuy2KoihKOkgYvrPrN58FNgBu4FFjzC4R+TqwzRizHlgL/FREanG8lzX22l0i8iSwG/AB9xlj/ADRbNohvwSsE5EHgO3WNsM8RiIeSdxlQqD3YQC9FwPovXDQ+zDAsN0LcZwNRVEURUk/WhRNURRFyRhUlBRFUZSMQUUpjNEuS5QORORREWkQkZ1hbSUi8oIt7fSCiEy27SIi37L34x0RuSTsmrts//0icle0sTIZEZklIptEZI+I7BKRv7XtE/Fe5IjIGyLytr0X/2zbh63k11jCVoTZLiLP2vcT9T4cEpF3RWSHiGyzbSP/82GM0T/OupobOADMBbzA20BNuuc1At/zJNTnoQAAAw1JREFUWuASYGdY278BX7avvww8ZF+/D/gdzl6wy4Ettr0EqLN/T7avJ6f7u6V4H6YBl9jXhcB7OOWoJuK9EKDAvs4Cttjv+CSwxrZ/H/iMff1XwPft6zXAE/Z1jf25yQaq7M+TO93fbwj34wvAL4Bn7fuJeh8OAWURbSP+86Ge0gCjVpYonRhjXsbJXgwnvIRTZGmnnxiH13H2kE0DbgFeMMa0GGNagRdw6g6OGYwxJ4wxb9nXZ4A9OJVAJuK9MMaYDvs2y/4xDF/JrzGDiMwE3g/80L4fztJn44ER//lQURogWjmlGTH6jjemGGNOgPOwBipse6x7Mq7ulQ27XIzjIUzIe2FDVjuABpwHxwGSLPkFhJf8Guv34j+BvwMC9n3Spc8YX/cBnF9Mfi8ib4pTig1G4edDz1MaIJlyShONVEs7jTlEpAD4JfA5Y0y7xC6kO67vhXH29l0kIpOAXwHnR+tm/x6X90JEbgMajDFvisiKYHOUruP6PoRxlTHmuIhUAC+IyN44fYftXqinNMBELkt0yrra2L8bbPu4LuEkIlk4gvRzY8zTtnlC3osgxpjTwEs46wLDVfJrrHAV8EEROYQTvl+J4zlNtPsAgDHmuP27AecXlWWMws+HitIAE7ksUXgJp8jSTn9mM2suB9qsy74BuFlEJtvsm5tt25jBxv7XAnuMMd8I+2gi3oty6yEhIrnAjThrbMNV8mtMYIy53xgz0ziFRdfgfK9PMMHuA4CI5ItIYfA1zv/rnYzGz0e6Mzwy6Q9OBsl7OPH0f0j3fEboOz4OnAD6cX6LuQcnDv4isN/+XWL7Cs5BiQeAd4GlYXbuxlnArQX+PN3fawj34WqcMMI7wA77530T9F4swSnp9Y598HzVts/FeZjWAv8NZNv2HPu+1n4+N8zWP9h7tA+4Nd3f7RzuyQoGsu8m3H2w3/lt+2dX8Hk4Gj8fWmZIURRFyRg0fKcoiqJkDCpKiqIoSsagoqQoiqJkDCpKiqIoSsagoqQoiqJkDCpKiqIoSsagoqQoiqJkDP8/XqW/tC349s4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.clf()\n",
    "# plt.plot(x_part, calcs, '.')\n",
    "plt.plot(test_q, calcs_test, label = 'pdf')\n",
    "# plt.plot(test_q, f0_y, label = '0')\n",
    "# plt.plot(test_q, fT_y, label = 'T')\n",
    "# plt.plot(test_q, fplus_y, label = '+')\n",
    "# plt.plot(test_q, res_y, label = 'res')\n",
    "plt.legend()\n",
    "plt.ylim(0.0, 1.5e-6)\n",
    "# plt.yscale('log')\n",
    "# plt.xlim(770, 785)\n",
    "plt.savefig('test.png')\n",
    "# print(jpsi_width)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "# probs = mixture.prob(test_q)\n",
    "# probs_np = zfit.run(probs)\n",
    "# probs_np *= np.max(calcs_test) / np.max(probs_np)\n",
    "# plt.figure()\n",
    "# plt.semilogy(test_q, probs_np,label=\"importance sampling\")\n",
    "# plt.semilogy(test_q, calcs_test, label = 'pdf')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 0.213/(0.00133+0.213+0.015)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Adjust scaling of different parts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "total_f.update_integration_options(draws_per_dim=2000000, mc_sampler=None)\n",
    "# inte = total_f.integrate(limits = (950., 1050.), norm_range=False)\n",
    "# inte_fl = zfit.run(inte)\n",
    "# print(inte_fl/4500)\n",
    "# print(pdg[\"jpsi_BR\"]/pdg[\"NR_BR\"], inte_fl*pdg[\"psi2s_auc\"]/pdg[\"NR_auc\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # print(\"jpsi:\", inte_fl)\n",
    "# # print(\"Increase am by factor:\", np.sqrt(pdg[\"jpsi_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "# # print(\"New amp:\", pdg[\"jpsi\"][3]*np.sqrt(pdg[\"jpsi_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "\n",
    "# # print(\"psi2s:\", inte_fl)\n",
    "# # print(\"Increase am by factor:\", np.sqrt(pdg[\"psi2s_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "# # print(\"New amp:\", pdg[\"psi2s\"][3]*np.sqrt(pdg[\"psi2s_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "\n",
    "# name = \"phi\"\n",
    "\n",
    "# print(name+\":\", inte_fl)\n",
    "# print(\"Increase am by factor:\", np.sqrt(pdg[name+\"_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "# print(\"New amp:\", pdg[name][0]*np.sqrt(pdg[name+\"_BR\"]/pdg[\"NR_BR\"]*pdg[\"NR_auc\"]/inte_fl))\n",
    "\n",
    "\n",
    "# print(x_min)\n",
    "# print(x_max)\n",
    "# # total_f.update_integration_options(draws_per_dim=2000000, mc_sampler=None)\n",
    "# total_f.update_integration_options(mc_sampler=lambda dim, num_results,\n",
    "#                                     dtype: tf.random_uniform(maxval=1., shape=(num_results, dim), dtype=dtype),\n",
    "#                                    draws_per_dim=1000000)\n",
    "# # _ = []\n",
    "\n",
    "# # for i in range(10):\n",
    "\n",
    "# #     inte = total_f.integrate(limits = (x_min, x_max))\n",
    "# #     inte_fl = zfit.run(inte)\n",
    "# #     print(inte_fl)\n",
    "# #     _.append(inte_fl)\n",
    "\n",
    "# # print(\"mean:\", np.mean(_))\n",
    "\n",
    "# _ = time.time()\n",
    "\n",
    "# inte = total_f.integrate(limits = (x_min, x_max))\n",
    "# inte_fl = zfit.run(inte)\n",
    "# print(inte_fl)\n",
    "# print(\"Time taken: {}\".format(display_time(int(time.time() - _))))\n",
    "\n",
    "# print(pdg['NR_BR']/pdg['NR_auc']*inte_fl)\n",
    "# print(0.25**2*4.2/1000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sampling\n",
    "## Mixture distribution for sampling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "    \n",
    "# print(list_of_borders[:9])\n",
    "# print(list_of_borders[-9:])\n",
    "\n",
    "\n",
    "class UniformSampleAndWeights(zfit.util.execution.SessionHolderMixin):\n",
    "    def __call__(self, limits, dtype, n_to_produce):\n",
    "        # n_to_produce = tf.cast(n_to_produce, dtype=tf.int32)\n",
    "        low, high = limits.limit1d\n",
    "        low = tf.cast(low, dtype=dtype)\n",
    "        high = tf.cast(high, dtype=dtype)\n",
    "#         uniform = tfd.Uniform(low=low, high=high)\n",
    "#         uniformjpsi = tfd.Uniform(low=tf.constant(3080, dtype=dtype), high=tf.constant(3112, dtype=dtype))\n",
    "#         uniformpsi2s = tfd.Uniform(low=tf.constant(3670, dtype=dtype), high=tf.constant(3702, dtype=dtype))\n",
    "\n",
    "#         list_of_borders = []\n",
    "#         _p = []\n",
    "#         splits = 10\n",
    "\n",
    "#         _ = np.linspace(x_min, x_max, splits)\n",
    "\n",
    "#         for i in range(splits):\n",
    "#             list_of_borders.append(tf.constant(_[i], dtype=dtype))\n",
    "#             _p.append(tf.constant(1/splits, dtype=dtype))\n",
    "    \n",
    "#         mixture = tfd.MixtureSameFamily(mixture_distribution=tfd.Categorical(probs=_p[:(splits-1)]),\n",
    "#                                         components_distribution=tfd.Uniform(low=list_of_borders[:(splits-1)], \n",
    "#                                                                             high=list_of_borders[-(splits-1):]))\n",
    "        mixture = tfd.MixtureSameFamily(mixture_distribution=tfd.Categorical(probs=[tf.constant(0.05, dtype=dtype),\n",
    "                                                                                    tf.constant(0.93, dtype=dtype),\n",
    "                                                                                    tf.constant(0.05, dtype=dtype),\n",
    "                                                                                    tf.constant(0.065, dtype=dtype),\n",
    "                                                                                    tf.constant(0.04, dtype=dtype),\n",
    "                                                                                    tf.constant(0.05, dtype=dtype)]),\n",
    "                                        components_distribution=tfd.Uniform(low=[tf.constant(x_min, dtype=dtype), \n",
    "                                                                                 tf.constant(3090, dtype=dtype),\n",
    "                                                                                 tf.constant(3681, dtype=dtype), \n",
    "                                                                                 tf.constant(3070, dtype=dtype),\n",
    "                                                                                 tf.constant(1000, dtype=dtype),\n",
    "                                                                                 tf.constant(3660, dtype=dtype)], \n",
    "                                                                            high=[tf.constant(x_max, dtype=dtype),\n",
    "                                                                                  tf.constant(3102, dtype=dtype), \n",
    "                                                                                  tf.constant(3691, dtype=dtype),\n",
    "                                                                                  tf.constant(3110, dtype=dtype),\n",
    "                                                                                  tf.constant(1040, dtype=dtype),\n",
    "                                                                                  tf.constant(3710, dtype=dtype)]))\n",
    "#         dtype = tf.float64\n",
    "#         mixture = tfd.MixtureSameFamily(mixture_distribution=tfd.Categorical(probs=[tf.constant(0.04, dtype=dtype),\n",
    "#                                                                                     tf.constant(0.90, dtype=dtype),\n",
    "#                                                                                     tf.constant(0.02, dtype=dtype),\n",
    "#                                                                                     tf.constant(0.07, dtype=dtype),\n",
    "#                                                                                     tf.constant(0.02, dtype=dtype)]),\n",
    "#                                         components_distribution=tfd.Uniform(low=[tf.constant(x_min, dtype=dtype), \n",
    "#                                                                                  tf.constant(3089, dtype=dtype),\n",
    "#                                                                                  tf.constant(3103, dtype=dtype), \n",
    "#                                                                                  tf.constant(3681, dtype=dtype),\n",
    "#                                                                                  tf.constant(3691, dtype=dtype)], \n",
    "#                                                                             high=[tf.constant(3089, dtype=dtype),\n",
    "#                                                                                   tf.constant(3103, dtype=dtype), \n",
    "#                                                                                   tf.constant(3681, dtype=dtype),\n",
    "#                                                                                   tf.constant(3691, dtype=dtype), \n",
    "#                                                                                   tf.constant(x_max, dtype=dtype)]))\n",
    "#         mixture = tfd.Uniform(tf.constant(x_min, dtype=dtype), tf.constant(x_max, dtype=dtype))\n",
    "#         sample = tf.random.uniform((n_to_produce, 1), dtype=dtype)\n",
    "        sample = mixture.sample((n_to_produce, 1))\n",
    "#         sample = tf.random.uniform((n_to_produce, 1), dtype=dtype)\n",
    "        weights = mixture.prob(sample)[:,0]\n",
    "#         weights = tf.broadcast_to(tf.constant(1., dtype=dtype), shape=(n_to_produce,))\n",
    "        # sample = tf.expand_dims(sample, axis=-1)\n",
    "#         print(sample, weights)\n",
    "        \n",
    "#         weights = tf.ones(shape=(n_to_produce,), dtype=dtype)\n",
    "        weights_max = None\n",
    "        thresholds = tf.random_uniform(shape=(n_to_produce,), dtype=dtype)\n",
    "        return sample, thresholds, weights, weights_max, n_to_produce"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# total_f._sample_and_weights = UniformSampleAndWeights"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 0.00133/(0.00133+0.213+0.015)*(x_max-3750)/(x_max-x_min)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "# zfit.settings.set_verbosity(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# # zfit.run.numeric_checks = False   \n",
    "\n",
    "# nr_of_toys = 1\n",
    "# nevents = int(pdg[\"number_of_decays\"])\n",
    "# nevents = pdg[\"number_of_decays\"]\n",
    "# event_stack = 1000000\n",
    "# # zfit.settings.set_verbosity(10)\n",
    "# calls = int(nevents/event_stack + 1)\n",
    "\n",
    "# total_samp = []\n",
    "\n",
    "# start = time.time()\n",
    "\n",
    "# sampler = total_f.create_sampler(n=event_stack)\n",
    "\n",
    "# for toy in range(nr_of_toys):\n",
    "    \n",
    "#     dirName = 'data/zfit_toys/toy_{0}'.format(toy)\n",
    "    \n",
    "#     if not os.path.exists(dirName):\n",
    "#         os.mkdir(dirName)\n",
    "#         print(\"Directory \" , dirName ,  \" Created \")\n",
    "\n",
    "#     for call in range(calls):\n",
    "\n",
    "#         sampler.resample(n=event_stack)\n",
    "#         s = sampler.unstack_x()\n",
    "#         sam = zfit.run(s)\n",
    "# #         clear_output(wait=True)\n",
    "\n",
    "#         c = call + 1\n",
    "        \n",
    "#         print(\"{0}/{1} of Toy {2}/{3}\".format(c, calls, toy+1, nr_of_toys))\n",
    "#         print(\"Time taken: {}\".format(display_time(int(time.time() - start))))\n",
    "#         print(\"Projected time left: {}\".format(display_time(int((time.time() - start)/(c+calls*(toy))*((nr_of_toys-toy)*calls-c)))))\n",
    "\n",
    "#         with open(\"data/zfit_toys/toy_{0}/{1}.pkl\".format(toy, call), \"wb\") as f:\n",
    "#             pkl.dump(sam, f, pkl.HIGHEST_PROTOCOL)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with open(r\"data/zfit_toys/toy_0/0.pkl\", \"rb\") as input_file:\n",
    "#     sam = pkl.load(input_file)\n",
    "# print(sam[:10])\n",
    "\n",
    "# with open(r\"data/zfit_toys/toy_0/1.pkl\", \"rb\") as input_file:\n",
    "#     sam2 = pkl.load(input_file)\n",
    "# print(sam2[:10])\n",
    "\n",
    "# print(np.sum(sam-sam2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(\"Time to generate full toy: {} s\".format(int(time.time()-start)))\n",
    "\n",
    "# total_samp = []\n",
    "\n",
    "# for call in range(calls):\n",
    "#     with open(r\"data/zfit_toys/toy_0/{}.pkl\".format(call), \"rb\") as input_file:\n",
    "#         sam = pkl.load(input_file)\n",
    "#         total_samp = np.append(total_samp, sam)\n",
    "\n",
    "# total_samp = total_samp.astype('float64')\n",
    "\n",
    "# data2 = zfit.data.Data.from_numpy(array=total_samp[:int(nevents)], obs=obs)\n",
    "\n",
    "# data3 = zfit.data.Data.from_numpy(array=total_samp, obs=obs)\n",
    "\n",
    "# print(total_samp[:nevents].shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plt.clf()\n",
    "\n",
    "# bins = int((x_max-x_min)/7)\n",
    "\n",
    "# # calcs = zfit.run(total_test_tf(samp))\n",
    "# print(total_samp[:nevents].shape)\n",
    "\n",
    "# plt.hist(total_samp[:nevents], bins = bins, range = (x_min,x_max), label = 'data')\n",
    "# # plt.plot(test_q, calcs_test*nevents , label = 'pdf')\n",
    "\n",
    "# # plt.plot(sam, calcs, '.')\n",
    "# # plt.plot(test_q, calcs_test)\n",
    "# # plt.yscale('log')\n",
    "# plt.ylim(0, 200)\n",
    "# # plt.xlim(3080, 3110)\n",
    "\n",
    "# plt.legend()\n",
    "\n",
    "# plt.savefig('test2.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sampler = total_f.create_sampler(n=nevents)\n",
    "# nll = zfit.loss.UnbinnedNLL(model=total_f, data=sampler, fit_range = (x_min, x_max))\n",
    "\n",
    "# # for param in pdf.get_dependents():\n",
    "# #     param.set_value(initial_value)\n",
    "\n",
    "# sampler.resample(n=nevents)\n",
    "\n",
    "# # Randomise initial values\n",
    "# # for param in pdf.get_dependents():\n",
    "# #     param.set_value(random value here)\n",
    "\n",
    "# # Minimise the NLL\n",
    "# minimizer = zfit.minimize.MinuitMinimizer(verbosity = 10)\n",
    "# minimum = minimizer.minimize(nll)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "# jpsi_width"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plt.hist(sample, weights=1 / prob(sample))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fitting"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "# start = time.time()\n",
    "\n",
    "# for param in total_f.get_dependents():\n",
    "#     param.randomize()\n",
    "    \n",
    "# # for param in total_f.get_dependents():\n",
    "# #     print(zfit.run(param))\n",
    "    \n",
    "# nll = zfit.loss.UnbinnedNLL(model=total_f, data=data2, fit_range = (x_min, x_max))\n",
    "\n",
    "# minimizer = zfit.minimize.MinuitMinimizer(verbosity = 5)\n",
    "# # minimizer._use_tfgrad = False\n",
    "# result = minimizer.minimize(nll)\n",
    "\n",
    "# # param_errors = result.error()\n",
    "\n",
    "# # for var, errors in param_errors.items():\n",
    "# #     print('{}: ^{{+{}}}_{{{}}}'.format(var.name, errors['upper'], errors['lower']))\n",
    "\n",
    "# print(\"Function minimum:\", result.fmin)\n",
    "# # print(\"Results:\", result.params)\n",
    "# print(\"Hesse errors:\", result.hesse())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(\"Time taken for fitting: {}\".format(display_time(int(time.time()-start))))\n",
    "\n",
    "# # probs = total_f.pdf(test_q)\n",
    "\n",
    "# calcs_test = zfit.run(probs)\n",
    "# res_y = zfit.run(jpsi_res(test_q))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plt.clf()\n",
    "# # plt.plot(x_part, calcs, '.')\n",
    "# plt.plot(test_q, calcs_test, label = 'pdf')\n",
    "# # plt.plot(test_q, res_y, label = 'res')\n",
    "# plt.legend()\n",
    "# plt.ylim(0.0, 10e-6)\n",
    "# # plt.yscale('log')\n",
    "# # plt.xlim(3080, 3110)\n",
    "# plt.savefig('test3.png')\n",
    "# # print(jpsi_width)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "# _tot = 4.37e-7+6.02e-5+4.97e-6\n",
    "# _probs = []\n",
    "# _probs.append(6.02e-5/_tot)\n",
    "# _probs.append(4.97e-6/_tot)\n",
    "# _probs.append(4.37e-7/_tot)\n",
    "# print(_probs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [],
   "source": [
    "# dtype = 'float64'\n",
    "# # mixture = tfd.Uniform(tf.constant(x_min, dtype=dtype), tf.constant(x_max, dtype=dtype))\n",
    "# mixture = tfd.MixtureSameFamily(mixture_distribution=tfd.Categorical(probs=[tf.constant(0.007, dtype=dtype),\n",
    "#                                                                             tf.constant(0.917, dtype=dtype),\n",
    "#                                                                             tf.constant(0.076, dtype=dtype)]),\n",
    "#                                 components_distribution=tfd.Uniform(low=[tf.constant(x_min, dtype=dtype), \n",
    "#                                                                          tf.constant(3080, dtype=dtype),\n",
    "#                                                                          tf.constant(3670, dtype=dtype)], \n",
    "#                                                                     high=[tf.constant(x_max, dtype=dtype),\n",
    "#                                                                           tf.constant(3112, dtype=dtype), \n",
    "#                                                                           tf.constant(3702, dtype=dtype)]))\n",
    "# # for i in range(10):\n",
    "# #     print(zfit.run(mixture.prob(mixture.sample((10, 1)))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print((zfit.run(jpsi_p)%(2*np.pi))/np.pi)\n",
    "# print((zfit.run(psi2s_p)%(2*np.pi))/np.pi)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "#         def jpsi_res(q):\n",
    "#             return resonance(q, _mass = jpsi_mass, scale = jpsi_scale,\n",
    "#                              phase = jpsi_phase, width = jpsi_width)\n",
    "\n",
    "#         def psi2s_res(q):\n",
    "#             return resonance(q, _mass = psi2s_mass, scale = psi2s_scale,\n",
    "#                              phase = psi2s_phase, width = psi2s_width)\n",
    "        \n",
    "#         def p3770_res(q):\n",
    "#             return resonance(q, _mass = p3770_mass, scale = p3770_scale,\n",
    "#                              phase = p3770_phase, width = p3770_width)\n",
    "        \n",
    "#         def p4040_res(q):\n",
    "#             return resonance(q, _mass = p4040_mass, scale = p4040_scale,\n",
    "#                              phase = p4040_phase, width = p4040_width)\n",
    "        \n",
    "#         def p4160_res(q):\n",
    "#             return resonance(q, _mass = p4160_mass, scale = p4160_scale,\n",
    "#                              phase = p4160_phase, width = p4160_width)\n",
    "        \n",
    "#         def p4415_res(q):\n",
    "#             return resonance(q, _mass = p4415_mass, scale = p4415_scale,\n",
    "#                              phase = p4415_phase, width = p4415_width)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 0.15**2*4.2/1000\n",
    "# result.hesse()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Constraints"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 1. Constraint - Real part of sum of Psi contrib and D contribs\n",
    "\n",
    "sum_list = []\n",
    "\n",
    "sum_list.append(ztf.to_complex(jpsi_s) * tf.exp(tf.complex(ztf.constant(0.0), jpsi_p)) * ztf.to_complex(jpsi_w / (tf.pow(jpsi_m,3))))\n",
    "sum_list.append(ztf.to_complex(psi2s_s) * tf.exp(tf.complex(ztf.constant(0.0), psi2s_p)) * ztf.to_complex(psi2s_w / (tf.pow(psi2s_m,3))))\n",
    "sum_list.append(ztf.to_complex(p3770_s) * tf.exp(tf.complex(ztf.constant(0.0), p3770_p)) * ztf.to_complex(p3770_w / (tf.pow(p3770_m,3))))\n",
    "sum_list.append(ztf.to_complex(p4040_s) * tf.exp(tf.complex(ztf.constant(0.0), p4040_p)) * ztf.to_complex(p4040_w / (tf.pow(p4040_m,3))))\n",
    "sum_list.append(ztf.to_complex(p4160_s) * tf.exp(tf.complex(ztf.constant(0.0), p4160_p)) * ztf.to_complex(p4160_w / (tf.pow(p4160_m,3))))\n",
    "sum_list.append(ztf.to_complex(p4415_s) * tf.exp(tf.complex(ztf.constant(0.0), p4415_p)) * ztf.to_complex(p4415_w / (tf.pow(p4415_m,3))))\n",
    "sum_list.append(ztf.to_complex(DDstar_s) * tf.exp(tf.complex(ztf.constant(0.0), DDstar_p)) * (ztf.to_complex(1.0 / (10.0*tf.pow(Dstar_m,2)) + 1.0 / (10.0*tf.pow(D_m,2)))))\n",
    "sum_list.append(ztf.to_complex(Dbar_s) * tf.exp(tf.complex(ztf.constant(0.0), Dbar_p)) * ztf.to_complex(1.0 / (6.0*tf.pow(Dbar_m,2))))\n",
    "\n",
    "sum_ru_1 = ztf.to_complex(ztf.constant(0.0))\n",
    "\n",
    "for part in sum_list:\n",
    "    sum_ru_1 += part\n",
    "\n",
    "sum_1 = tf.math.real(sum_ru_1)\n",
    "# constraint1 = zfit.constraint.GaussianConstraint(params = sum_1, mu = ztf.constant(1.7*10**-8), \n",
    "#                                                  sigma = ztf.constant(2.2*10**-8))\n",
    "\n",
    "constraint1 = tf.pow((sum_1-ztf.constant(1.7*10**-8))/ztf.constant(2.2*10**-8),2)/ztf.constant(2.)\n",
    "\n",
    "# 2. Constraint - Abs. of sum of Psi contribs and D contribs\n",
    "\n",
    "sum_2 = tf.abs(sum_ru_1)\n",
    "constraint2 = tf.cond(tf.greater_equal(sum_2, 5.0e-8), lambda: 100000., lambda: 0.)\n",
    "\n",
    "# 3. Constraint - Maximum eta of D contribs\n",
    "\n",
    "constraint3_0 = tf.cond(tf.greater_equal(tf.abs(Dbar_s), 0.2), lambda: 100000., lambda: 0.)\n",
    "\n",
    "constraint3_1 = tf.cond(tf.greater_equal(tf.abs(DDstar_s), 0.2), lambda: 100000., lambda: 0.)\n",
    "\n",
    "# 4. Constraint - Formfactor multivariant gaussian covariance fplus\n",
    "\n",
    "Cov_matrix = [[ztf.constant(   1.), ztf.constant( 0.45), ztf.constant( 0.19), ztf.constant(0.857), ztf.constant(0.598), ztf.constant(0.531), ztf.constant(0.752), ztf.constant(0.229), ztf.constant(0,117)],\n",
    "              [ztf.constant( 0.45), ztf.constant(   1.), ztf.constant(0.677), ztf.constant(0.708), ztf.constant(0.958), ztf.constant(0.927), ztf.constant(0.227), ztf.constant(0.443), ztf.constant(0.287)],\n",
    "              [ztf.constant( 0.19), ztf.constant(0.677), ztf.constant(   1.), ztf.constant(0.595), ztf.constant(0.770), ztf.constant(0.819),ztf.constant(-0.023), ztf.constant( 0.07), ztf.constant(0.196)],\n",
    "              [ztf.constant(0.857), ztf.constant(0.708), ztf.constant(0.595), ztf.constant(   1.), ztf.constant( 0.83), ztf.constant(0.766), ztf.constant(0.582), ztf.constant(0.237), ztf.constant(0.192)],\n",
    "              [ztf.constant(0.598), ztf.constant(0.958), ztf.constant(0.770), ztf.constant( 0.83), ztf.constant(   1.), ztf.constant(0.973), ztf.constant(0.324), ztf.constant(0.372), ztf.constant(0.272)],\n",
    "              [ztf.constant(0.531), ztf.constant(0.927), ztf.constant(0.819), ztf.constant(0.766), ztf.constant(0.973), ztf.constant(   1.), ztf.constant(0.268), ztf.constant(0.332), ztf.constant(0.269)],\n",
    "              [ztf.constant(0.752), ztf.constant(0.227),ztf.constant(-0.023), ztf.constant(0.582), ztf.constant(0.324), ztf.constant(0.268), ztf.constant(   1.), ztf.constant( 0.59), ztf.constant(0.515)],\n",
    "              [ztf.constant(0.229), ztf.constant(0.443), ztf.constant( 0.07), ztf.constant(0.237), ztf.constant(0.372), ztf.constant(0.332), ztf.constant( 0.59), ztf.constant(   1.), ztf.constant(0.897)],\n",
    "              [ztf.constant(0.117), ztf.constant(0.287), ztf.constant(0.196), ztf.constant(0.192), ztf.constant(0.272), ztf.constant(0.269), ztf.constant(0.515), ztf.constant(0.897), ztf.constant(   1.)]]\n",
    "\n",
    "def triGauss(val1,val2,val3,m = Cov_matrix):\n",
    "\n",
    "    mean1 = ztf.constant(0.466)\n",
    "    mean2 = ztf.constant(-0.885)\n",
    "    mean3 = ztf.constant(-0.213)\n",
    "    sigma1 = ztf.constant(0.014/3.)\n",
    "    sigma2 = ztf.constant(0.128/3.)\n",
    "    sigma3 = ztf.constant(0.548/3.)\n",
    "    x1 = (val1-mean1)/sigma1\n",
    "    x2 = (val2-mean2)/sigma2\n",
    "    x3 = (val3-mean3)/sigma3\n",
    "    rho12 = m[0][1]\n",
    "    rho13 = m[0][2]\n",
    "    rho23 = m[1][2]\n",
    "    w = x1*x1*(rho23*rho23-1) + x2*x2*(rho13*rho13-1)+x3*x3*(rho12*rho12-1)+2*(x1*x2*(rho12-rho13*rho23)+x1*x3*(rho13-rho12*rho23)+x2*x3*(rho23-rho12*rho13))\n",
    "    d = 2*(rho12*rho12+rho13*rho13+rho23*rho23-2*rho12*rho13*rho23-1)\n",
    "    \n",
    "    fcn = -w/d\n",
    "    chisq = -2*fcn\n",
    "    return chisq\n",
    "\n",
    "constraint4 = triGauss(bplus_0, bplus_1, bplus_2)\n",
    "\n",
    "# mean1 = ztf.constant(0.466)\n",
    "# mean2 = ztf.constant(-0.885)\n",
    "# mean3 = ztf.constant(-0.213)\n",
    "# sigma1 = ztf.constant(0.014)\n",
    "# sigma2 = ztf.constant(0.128)\n",
    "# sigma3 = ztf.constant(0.548)\n",
    "# constraint4_0 = tf.pow((bplus_0-mean1)/sigma1,2)/ztf.constant(2.)\n",
    "# constraint4_1 = tf.pow((bplus_1-mean2)/sigma2,2)/ztf.constant(2.)\n",
    "# constraint4_2 = tf.pow((bplus_2-mean3)/sigma3,2)/ztf.constant(2.)\n",
    "\n",
    "# 5. Constraint - Abs. of sum of light contribs\n",
    "\n",
    "sum_list_5 = []\n",
    "\n",
    "sum_list_5.append(rho_s*rho_w/rho_m)\n",
    "sum_list_5.append(omega_s*omega_w/omega_m)\n",
    "sum_list_5.append(phi_s*phi_w/phi_m)\n",
    "\n",
    "\n",
    "sum_ru_5 = ztf.constant(0.0)\n",
    "\n",
    "for part in sum_list_5:\n",
    "    sum_ru_5 += part\n",
    "\n",
    "constraint5 = tf.cond(tf.greater_equal(tf.abs(sum_ru_5), ztf.constant(0.02)), lambda: 100000., lambda: 0.)\n",
    "\n",
    "# 6. Constraint on phases of Jpsi and Psi2s for cut out fit\n",
    "\n",
    "\n",
    "# constraint6_0 = zfit.constraint.GaussianConstraint(params = jpsi_p, mu = ztf.constant(pdg[\"jpsi_phase_unc\"]),\n",
    "#                                                    sigma = ztf.constant(jpsi_phase))\n",
    "# constraint6_1 = zfit.constraint.GaussianConstraint(params = psi2s_p, mu = ztf.constant(pdg[\"psi2s_phase_unc\"]),\n",
    "#                                                    sigma = ztf.constant(psi2s_phase))\n",
    "\n",
    "constraint6_0  =  tf.pow((jpsi_p-ztf.constant(jpsi_phase))/ztf.constant(pdg[\"jpsi_phase_unc\"]),2)/ztf.constant(2.)\n",
    "constraint6_1  =  tf.pow((psi2s_p-ztf.constant(psi2s_phase))/ztf.constant(pdg[\"psi2s_phase_unc\"]),2)/ztf.constant(2.)\n",
    "\n",
    "# 7. Constraint on Ctt with higher limits\n",
    "\n",
    "constraint7 = tf.cond(tf.greater_equal(Ctt*Ctt, 0.25), lambda: 100000., lambda: 0.)\n",
    "\n",
    "constraint7dtype = tf.float64\n",
    "\n",
    "# zfit.run(constraint6_0)\n",
    "\n",
    "# ztf.convert_to_tensor(constraint6_0)\n",
    "\n",
    "#List of all constraints\n",
    "\n",
    "constraints = [constraint1, constraint2, constraint3_0, constraint3_1,# constraint4, #constraint4_0, constraint4_1, constraint4_2,\n",
    "               constraint6_0, constraint6_1]#, constraint7]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reset params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "def reset_param_values():   \n",
    "    jpsi_m.set_value(jpsi_mass)\n",
    "    jpsi_s.set_value(jpsi_scale)\n",
    "    jpsi_p.set_value(jpsi_phase)\n",
    "    jpsi_w.set_value(jpsi_width)\n",
    "    psi2s_m.set_value(psi2s_mass)\n",
    "    psi2s_s.set_value(psi2s_scale)\n",
    "    psi2s_p.set_value(psi2s_phase)\n",
    "    psi2s_w.set_value(psi2s_width)\n",
    "    p3770_m.set_value(p3770_mass)\n",
    "    p3770_s.set_value(p3770_scale)\n",
    "    p3770_p.set_value(p3770_phase)\n",
    "    p3770_w.set_value(p3770_width)\n",
    "    p4040_m.set_value(p4040_mass)\n",
    "    p4040_s.set_value(p4040_scale)\n",
    "    p4040_p.set_value(p4040_phase)\n",
    "    p4040_w.set_value(p4040_width)\n",
    "    p4160_m.set_value(p4160_mass)\n",
    "    p4160_s.set_value(p4160_scale)\n",
    "    p4160_p.set_value(p4160_phase)\n",
    "    p4160_w.set_value(p4160_width)\n",
    "    p4415_m.set_value(p4415_mass)\n",
    "    p4415_s.set_value(p4415_scale)\n",
    "    p4415_p.set_value(p4415_phase)\n",
    "    p4415_w.set_value(p4415_width)\n",
    "    rho_m.set_value(rho_mass)\n",
    "    rho_s.set_value(rho_scale)\n",
    "    rho_p.set_value(rho_phase)\n",
    "    rho_w.set_value(rho_width)\n",
    "    omega_m.set_value(omega_mass)\n",
    "    omega_s.set_value(omega_scale)\n",
    "    omega_p.set_value(omega_phase)\n",
    "    omega_w.set_value(omega_width)\n",
    "    phi_m.set_value(phi_mass)\n",
    "    phi_s.set_value(phi_scale)\n",
    "    phi_p.set_value(phi_phase)\n",
    "    phi_w.set_value(phi_width)\n",
    "    Dstar_m.set_value(Dstar_mass)\n",
    "    DDstar_s.set_value(0.0)\n",
    "    DDstar_p.set_value(0.0)\n",
    "    D_m.set_value(D_mass)\n",
    "    Dbar_m.set_value(Dbar_mass)\n",
    "    Dbar_s.set_value(0.0)\n",
    "    Dbar_p.set_value(0.0)\n",
    "    tau_m.set_value(pdg['tau_M'])\n",
    "    Ctt.set_value(0.0)\n",
    "    b0_0.set_value(0.292)\n",
    "    b0_1.set_value(0.281)\n",
    "    b0_2.set_value(0.150)\n",
    "    bplus_0.set_value(0.466)\n",
    "    bplus_1.set_value(-0.885)\n",
    "    bplus_2.set_value(-0.213)\n",
    "    bT_0.set_value(0.460)\n",
    "    bT_1.set_value(-1.089)\n",
    "    bT_2.set_value(-1.114)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
<<<<<<< HEAD
   "outputs": [],
=======
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py:163: to_int32 (from tensorflow.python.ops.math_ops) is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "Use tf.cast instead.\n",
      "Toy 0: Generating data...\n",
      "Toy 0: Data generation finished\n",
      "Toy 0: Loading data...\n",
      "Toy 0: Loading data finished\n",
      "Toy 0: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.979E+05               |     Ncalls=616 (627 total)     |\n",
      "| EDM = 1.12E+06 (Goal: 5E-06)  |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   True    |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   False   |   False   |  True  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297878.03533313115\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |  -0.017   |   0.133   |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |     4     |    11     |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |   19.27   |   0.29    |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   1.780   |   0.021   |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   5.93    |   0.24    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |   -4.8    |    0.4    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   2.54    |   0.05    |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |   4.45    |   0.17    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |  -2.268   |   0.017   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |   0.258   |   0.008   |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |  -0.233   |   0.016   |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   | -2.996E-1 |  0.005E-1 |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |   -0.41   |    0.20   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |   -3.57   |    0.17   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |   8.64    |   0.17    |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |  -0.647   |   0.003   |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   2.94    |   0.09    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   1.04    |   0.03    |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |   3.72    |   0.05    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |   1.365   |   0.028   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |   0.321   |   0.014   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |    0.5    |    0.3    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   0.515   |   0.026   |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |   3.95    |   0.07    |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000    0.965    0.862    0.207    0.933    0.950    0.074    0.882    0.631    0.797   -0.593   -0.035    0.707    0.703   -0.904   -0.323    0.781    0.492    0.179    0.878    0.961   -0.962   -0.116    0.384 |\n",
      "|  omega_p |    0.965    1.000    0.894    0.210    0.967    0.988    0.083    0.914    0.654    0.826   -0.638   -0.036    0.723    0.723   -0.937   -0.337    0.815    0.506    0.200    0.913    0.996   -0.997   -0.123    0.408 |\n",
      "|    phi_s |    0.862    0.894    1.000    0.187    0.864    0.883    0.074    0.817    0.585    0.738   -0.576   -0.032    0.646    0.646   -0.837   -0.305    0.728    0.452    0.179    0.816    0.890   -0.891   -0.110    0.364 |\n",
      "|  p4040_s |    0.207    0.210    0.187    1.000    0.203    0.207    0.024    0.192    0.135    0.171   -0.137   -0.008    0.151    0.147   -0.196   -0.056    0.155    0.095    0.050    0.192    0.210   -0.209   -0.023    0.070 |\n",
      "|    rho_p |    0.933    0.967    0.864    0.203    1.000    0.955    0.080    0.884    0.632    0.799   -0.615   -0.035    0.699    0.700   -0.906   -0.324    0.788    0.489    0.193    0.883    0.964   -0.964   -0.119    0.394 |\n",
      "|   Dbar_p |    0.950    0.988    0.883    0.207    0.955    1.000    0.077    0.903    0.657    0.813   -0.617   -0.036    0.706    0.713   -0.925   -0.313    0.798    0.498    0.195    0.898    0.986   -0.985   -0.120    0.397 |\n",
      "|  p3770_s |    0.074    0.083    0.074    0.024    0.080    0.077    1.000    0.076    0.077    0.068   -0.039   -0.004    0.048    0.064   -0.078   -0.021    0.072    0.047    0.013    0.070    0.083   -0.083   -0.010    0.030 |\n",
      "|    phi_p |    0.882    0.914    0.817    0.192    0.884    0.903    0.076    1.000    0.598    0.755   -0.582   -0.033    0.662    0.662   -0.857   -0.306    0.745    0.463    0.183    0.835    0.912   -0.912   -0.112    0.373 |\n",
      "|  psi2s_p |    0.631    0.654    0.585    0.135    0.632    0.657    0.077    0.598    1.000    0.544   -0.418   -0.020    0.542    0.472   -0.613   -0.238    0.547    0.331    0.172    0.597    0.650   -0.652   -0.085    0.287 |\n",
      "| DDstar_s |    0.797    0.826    0.738    0.171    0.799    0.813    0.068    0.755    0.544    1.000   -0.521   -0.030    0.593    0.593   -0.774   -0.258    0.669    0.414    0.168    0.751    0.825   -0.823   -0.101    0.336 |\n",
      "|  bplus_2 |   -0.593   -0.638   -0.576   -0.137   -0.615   -0.617   -0.039   -0.582   -0.418   -0.521    1.000    0.022   -0.482   -0.467    0.596    0.134   -0.494   -0.318   -0.111   -0.567   -0.643    0.637    0.066   -0.231 |\n",
      "|   Dbar_s |   -0.035   -0.036   -0.032   -0.008   -0.035   -0.036   -0.004   -0.033   -0.020   -0.030    0.022    1.000   -0.032   -0.026    0.034    0.013   -0.029   -0.019   -0.006   -0.032   -0.036    0.036    0.004   -0.014 |\n",
      "| DDstar_p |    0.707    0.723    0.646    0.151    0.699    0.706    0.048    0.662    0.542    0.593   -0.482   -0.032    1.000    0.512   -0.677   -0.200    0.610    0.356    0.170    0.700    0.723   -0.721   -0.094    0.340 |\n",
      "|  p4415_p |    0.703    0.723    0.646    0.147    0.700    0.713    0.064    0.662    0.472    0.593   -0.467   -0.026    0.512    1.000   -0.678   -0.210    0.598    0.366    0.152    0.657    0.724   -0.721   -0.053    0.318 |\n",
      "|  omega_s |   -0.904   -0.937   -0.837   -0.196   -0.906   -0.925   -0.078   -0.857   -0.613   -0.774    0.596    0.034   -0.677   -0.678    1.000    0.315   -0.763   -0.474   -0.187   -0.856   -0.934    0.934    0.115   -0.382 |\n",
      "|  bplus_1 |   -0.323   -0.337   -0.305   -0.056   -0.324   -0.313   -0.021   -0.306   -0.238   -0.258    0.134    0.013   -0.200   -0.210    0.315    1.000   -0.250   -0.146   -0.081   -0.280   -0.352    0.336    0.042   -0.139 |\n",
      "|  p4040_p |    0.781    0.815    0.728    0.155    0.788    0.798    0.072    0.745    0.547    0.669   -0.494   -0.029    0.610    0.598   -0.763   -0.250    1.000    0.424    0.193    0.738    0.814   -0.812   -0.096    0.292 |\n",
      "|  p4160_s |    0.492    0.506    0.452    0.095    0.489    0.498    0.047    0.463    0.331    0.414   -0.318   -0.019    0.356    0.366   -0.474   -0.146    0.424    1.000    0.112    0.459    0.506   -0.504   -0.054    0.175 |\n",
      "|  p3770_p |    0.179    0.200    0.179    0.050    0.193    0.195    0.013    0.183    0.172    0.168   -0.111   -0.006    0.170    0.152   -0.187   -0.081    0.193    0.112    1.000    0.175    0.198   -0.199   -0.028    0.096 |\n",
      "|   jpsi_p |    0.878    0.913    0.816    0.192    0.883    0.898    0.070    0.835    0.597    0.751   -0.567   -0.032    0.700    0.657   -0.856   -0.280    0.738    0.459    0.175    1.000    0.912   -0.910   -0.112    0.372 |\n",
      "|  bplus_0 |    0.961    0.996    0.890    0.210    0.964    0.986    0.083    0.912    0.650    0.825   -0.643   -0.036    0.723    0.724   -0.934   -0.352    0.814    0.506    0.198    0.912    1.000   -0.993   -0.122    0.406 |\n",
      "|    rho_s |   -0.962   -0.997   -0.891   -0.209   -0.964   -0.985   -0.083   -0.912   -0.652   -0.823    0.637    0.036   -0.721   -0.721    0.934    0.336   -0.812   -0.504   -0.199   -0.910   -0.993    1.000    0.122   -0.407 |\n",
      "|  p4415_s |   -0.116   -0.123   -0.110   -0.023   -0.119   -0.120   -0.010   -0.112   -0.085   -0.101    0.066    0.004   -0.094   -0.053    0.115    0.042   -0.096   -0.054   -0.028   -0.112   -0.122    0.122    1.000   -0.050 |\n",
      "|  p4160_p |    0.384    0.408    0.364    0.070    0.394    0.397    0.030    0.373    0.287    0.336   -0.231   -0.014    0.340    0.318   -0.382   -0.139    0.292    0.175    0.096    0.372    0.406   -0.407   -0.050    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.13307582874104307}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 11.312544525700241}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 0.29435979835990267}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.020773255462037232}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.2436242306495524}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.42032440023747597}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.04542769991233664}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 0.16705000022500105}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.01701667878468638}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.008402722883520175}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.015655568513703977}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.00046055330791114346}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.20138193109194136}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.17099716522674058}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 0.16848322966260731}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.0025357515987360912}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.08785001542580062}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.03182752693261037}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.04623833317689208}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.02804294337792035}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.01433694235612526}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.34127837859046173}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.02591065046357091}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.0690782454106591})])\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel_launcher.py:196: UserWarning: Creating legend with loc=\"best\" can be slow with large amounts of data.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 1/10\n",
      "Time taken: 1 min, 40 s\n",
      "Projected time left: 15 min\n",
      "Toy 1: Generating data...\n",
      "Toy 1: Data generation finished\n",
      "Toy 1: Loading data...\n",
      "Toy 1: Loading data finished\n",
      "Toy 1: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.979E+05               |    Ncalls=1287 (1287 total)    |\n",
      "| EDM = 6.42E-05 (Goal: 5E-06)  |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     True      |     True      |   False   |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   False   |   False   |  True  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297909.6307224464\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |   0.62    |   0.19    |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |   -5.6    |    0.4    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |   22.0    |    1.0    |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   1.09    |   0.18    |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   -0.30   |    0.71   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |    0.9    |    0.3    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   1.89    |   0.24    |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |   -5.53   |    0.15   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |   -2.08   |    0.03   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |   0.300   |   0.024   |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |   0.12    |   0.07    |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   |  -0.300   |   0.015   |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |    4.0    |    0.4    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |   -2.54   |    0.13   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |     8     |     4     |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |   0.73    |   0.06    |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   -2.70   |    0.17   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   2.60    |   0.17    |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |   4.24    |   0.14    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |   4.580   |   0.025   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |  -0.376   |   0.028   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |    0.6    |    0.4    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   1.73    |   0.19    |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |   -2.13   |    0.09   |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000   -0.014   -0.053    0.288    0.090   -0.352   -0.037   -0.015   -0.436    0.012   -0.639    0.002   -0.141    0.135   -0.012    0.102   -0.020    0.228   -0.274    0.467    0.026   -0.005    0.075   -0.120 |\n",
      "|  omega_p |   -0.014    1.000   -0.012   -0.000   -0.141    0.005    0.004   -0.005   -0.007    0.002    0.022    0.001    0.018   -0.001    0.827    0.177    0.004    0.001   -0.006   -0.041   -0.198    0.419    0.001   -0.006 |\n",
      "|    phi_s |   -0.053   -0.012    1.000   -0.006   -0.108   -0.014   -0.032    0.457    0.032   -0.001    0.158   -0.003   -0.007   -0.041    0.012   -0.146   -0.036   -0.039    0.016    0.056    0.161   -0.124    0.020   -0.017 |\n",
      "|  p4040_s |    0.288   -0.000   -0.006    1.000    0.005    0.025    0.156   -0.001   -0.200   -0.000   -0.001    0.011   -0.134   -0.143   -0.002   -0.021   -0.166   -0.070   -0.159    0.066    0.004    0.004    0.139   -0.417 |\n",
      "|    rho_p |    0.090   -0.141   -0.108    0.005    1.000    0.007    0.037   -0.294   -0.044    0.000   -0.253    0.003   -0.020    0.046    0.019   -0.216    0.031    0.047   -0.012   -0.049    0.235   -0.045   -0.027    0.024 |\n",
      "|   Dbar_p |   -0.352    0.005   -0.014    0.025    0.007    1.000   -0.057   -0.005    0.023   -0.006    0.026    0.012   -0.164   -0.152   -0.002   -0.074   -0.174   -0.002   -0.061   -0.232    0.002    0.023    0.058   -0.220 |\n",
      "|  p3770_s |   -0.037    0.004   -0.032    0.156    0.037   -0.057    1.000   -0.011    0.265   -0.004   -0.000   -0.016    0.076    0.005   -0.002   -0.020    0.033    0.090   -0.131   -0.082   -0.022    0.020    0.070   -0.048 |\n",
      "|    phi_p |   -0.015   -0.005    0.457   -0.001   -0.294   -0.005   -0.011    1.000    0.009    0.001    0.063   -0.000   -0.001   -0.010   -0.081    0.053   -0.009   -0.012    0.002    0.015   -0.052    0.036    0.006   -0.006 |\n",
      "|  psi2s_p |   -0.436   -0.007    0.032   -0.200   -0.044    0.023    0.265    0.009    1.000    0.019    0.156    0.026   -0.044    0.091    0.001    0.019    0.183   -0.094    0.238   -0.094    0.006   -0.026   -0.136    0.264 |\n",
      "| DDstar_s |    0.012    0.002   -0.001   -0.000    0.000   -0.006   -0.004    0.001    0.019    1.000    0.016   -0.001    0.032   -0.000    0.000   -0.032   -0.003    0.000    0.013    0.037    0.006    0.006    0.010    0.012 |\n",
      "|  bplus_2 |   -0.639    0.022    0.158   -0.001   -0.253    0.026   -0.000    0.063    0.156    0.016    1.000   -0.006    0.262   -0.234    0.023   -0.231   -0.096   -0.080    0.082   -0.264    0.032   -0.026    0.231   -0.163 |\n",
      "|   Dbar_s |    0.002    0.001   -0.003    0.011    0.003    0.012   -0.016   -0.000    0.026   -0.001   -0.006    1.000   -0.002   -0.001   -0.000   -0.011   -0.008    0.004    0.004    0.030    0.002    0.003    0.011   -0.003 |\n",
      "| DDstar_p |   -0.141    0.018   -0.007   -0.134   -0.020   -0.164    0.076   -0.001   -0.044    0.032    0.262   -0.002    1.000   -0.347   -0.001   -0.271   -0.165   -0.078    0.079   -0.296    0.030    0.056    0.049   -0.240 |\n",
      "|  p4415_p |    0.135   -0.001   -0.041   -0.143    0.046   -0.152    0.005   -0.010    0.091   -0.000   -0.234   -0.001   -0.347    1.000   -0.006    0.169    0.188   -0.050    0.066    0.149   -0.069    0.012   -0.183    0.383 |\n",
      "|  omega_s |   -0.012    0.827    0.012   -0.002    0.019   -0.002   -0.002   -0.081    0.001    0.000    0.023   -0.000   -0.001   -0.006    1.000   -0.000   -0.004   -0.005    0.001   -0.017    0.000    0.070    0.002   -0.003 |\n",
      "|  bplus_1 |    0.102    0.177   -0.146   -0.021   -0.216   -0.074   -0.020    0.053    0.019   -0.032   -0.231   -0.011   -0.271    0.169   -0.000    1.000    0.114   -0.027   -0.016    0.047   -0.942    0.437   -0.147    0.172 |\n",
      "|  p4040_p |   -0.020    0.004   -0.036   -0.166    0.031   -0.174    0.033   -0.009    0.183   -0.003   -0.096   -0.008   -0.165    0.188   -0.004    0.114    1.000    0.372    0.193    0.027   -0.072    0.022    0.074    0.090 |\n",
      "|  p4160_s |    0.228    0.001   -0.039   -0.070    0.047   -0.002    0.090   -0.012   -0.094    0.000   -0.080    0.004   -0.078   -0.050   -0.005   -0.027    0.372    1.000   -0.027    0.019   -0.011    0.020    0.330   -0.147 |\n",
      "|  p3770_p |   -0.274   -0.006    0.016   -0.159   -0.012   -0.061   -0.131    0.002    0.238    0.013    0.082    0.004    0.079    0.066    0.001   -0.016    0.193   -0.027    1.000   -0.069    0.020   -0.020   -0.089    0.192 |\n",
      "|   jpsi_p |    0.467   -0.041    0.056    0.066   -0.049   -0.232   -0.082    0.015   -0.094    0.037   -0.264    0.030   -0.296    0.149   -0.017    0.047    0.027    0.019   -0.069    1.000    0.124   -0.089   -0.060    0.092 |\n",
      "|  bplus_0 |    0.026   -0.198    0.161    0.004    0.235    0.002   -0.022   -0.052    0.006    0.006    0.032    0.002    0.030   -0.069    0.000   -0.942   -0.072   -0.011    0.020    0.124    1.000   -0.496    0.044   -0.050 |\n",
      "|    rho_s |   -0.005    0.419   -0.124    0.004   -0.045    0.023    0.020    0.036   -0.026    0.006   -0.026    0.003    0.056    0.012    0.070    0.437    0.022    0.020   -0.020   -0.089   -0.496    1.000   -0.004   -0.010 |\n",
      "|  p4415_s |    0.075    0.001    0.020    0.139   -0.027    0.058    0.070    0.006   -0.136    0.010    0.231    0.011    0.049   -0.183    0.002   -0.147    0.074    0.330   -0.089   -0.060    0.044   -0.004    1.000   -0.107 |\n",
      "|  p4160_p |   -0.120   -0.006   -0.017   -0.417    0.024   -0.220   -0.048   -0.006    0.264    0.012   -0.163   -0.003   -0.240    0.383   -0.003    0.172    0.090   -0.147    0.192    0.092   -0.050   -0.010   -0.107    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.18748768755829082}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 0.40323634554077215}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 0.985024000483012}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.17586720794456895}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.7105697644993554}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.32118284401931696}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.2372102719369873}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 0.14503339498061774}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.03185662332362904}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.02421642931361845}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.07483773848310715}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.01520721320902993}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.39795186227992385}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.13390670241756997}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 3.8706357513686984}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.057866643089661274}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.16823464853939285}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.171806968358724}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.138889598592427}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.02515049586844853}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.02826194125557635}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.39800274915520395}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.19429670305431568}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.08939828920188164})])\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 2/10\n",
      "Time taken: 3 min, 54 s\n",
      "Projected time left: 15 min, 36 s\n",
      "Toy 2: Generating data...\n",
      "Toy 2: Data generation finished\n",
      "Toy 2: Loading data...\n",
      "Toy 2: Loading data finished\n",
      "Toy 2: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.977E+05               |    Ncalls=1481 (1481 total)    |\n",
      "| EDM = 5.87E-05 (Goal: 5E-06)  |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     True      |     True      |   False   |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   True    |   True    | False  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297742.96509176795\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |  -0.017   |   0.193   |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |   6.28    |   0.09    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |   23.0    |    0.9    |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   1.56    |   0.17    |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   -0.55   |    0.27   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |    1.5    |    0.6    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   0.919   |   0.024   |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |   -5.27   |    0.12   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |   -2.04   |    0.03   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |   0.30    |   0.05    |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |   0.23    |   0.10    |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   |   -0.30   |    0.08   |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |  -0.011   |   0.472   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |    4.4    |    0.5    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |    6.0    |    0.9    |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |   1.02    |   0.08    |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   2.86    |   0.11    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   0.717   |   0.012   |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |  -6.283   |   0.019   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |   4.613   |   0.030   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |   -0.55   |    0.04   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |   1.37    |   0.31    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   0.40    |   0.19    |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |  -6.283   |   0.011   |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000   -0.001   -0.001    0.031    0.025   -0.498   -0.008   -0.001   -0.402   -0.011   -0.589    0.022    0.398    0.240    0.026    0.041    0.334    0.004    0.001    0.479    0.157   -0.079    0.022   -0.006 |\n",
      "|  omega_p |   -0.001    1.000   -0.002    0.000    0.021    0.001   -0.000    0.005    0.000   -0.000    0.001    0.000   -0.001   -0.000   -0.081    0.003   -0.000   -0.000    0.000    0.000   -0.004    0.014    0.001    0.000 |\n",
      "|    phi_s |   -0.001   -0.002    1.000   -0.024   -0.143   -0.082    0.000    0.331   -0.000    0.003    0.120   -0.020    0.003   -0.020   -0.004   -0.093   -0.005    0.001   -0.001    0.029    0.117   -0.034   -0.012    0.000 |\n",
      "|  p4040_s |    0.031    0.000   -0.024    1.000    0.037    0.373   -0.006   -0.007   -0.081    0.000   -0.015    0.075   -0.081   -0.213   -0.013   -0.067   -0.269   -0.007    0.005    0.017   -0.057    0.028    0.104    0.010 |\n",
      "|    rho_p |    0.025    0.021   -0.143    0.037    1.000    0.123   -0.000   -0.248   -0.013   -0.001   -0.253    0.034    0.007    0.041    0.159    0.008    0.008   -0.001    0.001   -0.051   -0.048    0.161    0.008   -0.000 |\n",
      "|   Dbar_p |   -0.498    0.001   -0.082    0.373    0.123    1.000    0.002   -0.029    0.275    0.032   -0.054    0.078   -0.300   -0.262   -0.030   -0.082   -0.417   -0.014    0.013   -0.132   -0.192    0.070    0.184    0.007 |\n",
      "|  p3770_s |   -0.008   -0.000    0.000   -0.006   -0.000    0.002    1.000    0.000    0.009   -0.001    0.001    0.002    0.004    0.002    0.000   -0.001    0.005    0.000   -0.002    0.002    0.001   -0.000   -0.002   -0.000 |\n",
      "|    phi_p |   -0.001    0.005    0.331   -0.007   -0.248   -0.029    0.000    1.000    0.009   -0.001    0.073   -0.006    0.014   -0.011   -0.117    0.021    0.003    0.000   -0.000    0.023   -0.010    0.131    0.002    0.000 |\n",
      "|  psi2s_p |   -0.402    0.000   -0.000   -0.081   -0.013    0.275    0.009    0.009    1.000    0.022    0.188    0.133    0.236   -0.040   -0.017   -0.020    0.042    0.006   -0.006    0.122   -0.073    0.041   -0.030    0.001 |\n",
      "| DDstar_s |   -0.011   -0.000    0.003    0.000   -0.001    0.032   -0.001   -0.001    0.022    1.000   -0.049    0.006    0.027    0.029    0.007    0.014    0.038    0.000    0.000    0.012    0.030   -0.017    0.003   -0.001 |\n",
      "|  bplus_2 |   -0.589    0.001    0.120   -0.015   -0.253   -0.054    0.001    0.073    0.188   -0.049    1.000   -0.047   -0.071   -0.245   -0.049    0.058   -0.042   -0.003   -0.003   -0.221   -0.180    0.123    0.161    0.001 |\n",
      "|   Dbar_s |    0.022    0.000   -0.020    0.075    0.034    0.078    0.002   -0.006    0.133    0.006   -0.047    1.000   -0.092   -0.038   -0.003   -0.010   -0.000   -0.001    0.003    0.104   -0.030    0.006    0.055    0.001 |\n",
      "| DDstar_p |    0.398   -0.001    0.003   -0.081    0.007   -0.300    0.004    0.014    0.236    0.027   -0.071   -0.092    1.000    0.139    0.006    0.015    0.483    0.015    0.000    0.635    0.033   -0.017    0.036   -0.005 |\n",
      "|  p4415_p |    0.240   -0.000   -0.020   -0.213    0.041   -0.262    0.002   -0.011   -0.040    0.029   -0.245   -0.038    0.139    1.000    0.013    0.068    0.076    0.011   -0.003    0.106    0.034   -0.035   -0.149    0.001 |\n",
      "|  omega_s |    0.026   -0.081   -0.004   -0.013    0.159   -0.030    0.000   -0.117   -0.017    0.007   -0.049   -0.003    0.006    0.013    1.000   -0.122    0.006    0.001   -0.000   -0.002    0.132   -0.460   -0.017   -0.000 |\n",
      "|  bplus_1 |    0.041    0.003   -0.093   -0.067    0.008   -0.082   -0.001    0.021   -0.020    0.014    0.058   -0.010    0.015    0.068   -0.122    1.000    0.065    0.002   -0.002   -0.041   -0.941    0.268   -0.050   -0.002 |\n",
      "|  p4040_p |    0.334   -0.000   -0.005   -0.269    0.008   -0.417    0.005    0.003    0.042    0.038   -0.042   -0.000    0.483    0.076    0.006    0.065    1.000    0.001   -0.009    0.265    0.011   -0.017    0.092   -0.004 |\n",
      "|  p4160_s |    0.004   -0.000    0.001   -0.007   -0.001   -0.014    0.000    0.000    0.006    0.000   -0.003   -0.001    0.015    0.011    0.001    0.002    0.001    1.000   -0.000    0.009    0.003   -0.001   -0.001   -0.002 |\n",
      "|  p3770_p |    0.001    0.000   -0.001    0.005    0.001    0.013   -0.002   -0.000   -0.006    0.000   -0.003    0.003    0.000   -0.003   -0.000   -0.002   -0.009   -0.000    1.000    0.004   -0.001    0.001    0.002    0.000 |\n",
      "|   jpsi_p |    0.479    0.000    0.029    0.017   -0.051   -0.132    0.002    0.023    0.122    0.012   -0.221    0.104    0.635    0.106   -0.002   -0.041    0.265    0.009    0.004    1.000    0.103   -0.017    0.019   -0.002 |\n",
      "|  bplus_0 |    0.157   -0.004    0.117   -0.057   -0.048   -0.192    0.001   -0.010   -0.073    0.030   -0.180   -0.030    0.033    0.034    0.132   -0.941    0.011    0.003   -0.001    0.103    1.000   -0.299   -0.080    0.000 |\n",
      "|    rho_s |   -0.079    0.014   -0.034    0.028    0.161    0.070   -0.000    0.131    0.041   -0.017    0.123    0.006   -0.017   -0.035   -0.460    0.268   -0.017   -0.001    0.001   -0.017   -0.299    1.000    0.041    0.000 |\n",
      "|  p4415_s |    0.022    0.001   -0.012    0.104    0.008    0.184   -0.002    0.002   -0.030    0.003    0.161    0.055    0.036   -0.149   -0.017   -0.050    0.092   -0.001    0.002    0.019   -0.080    0.041    1.000   -0.005 |\n",
      "|  p4160_p |   -0.006    0.000    0.000    0.010   -0.000    0.007   -0.000    0.000    0.001   -0.001    0.001    0.001   -0.005    0.001   -0.000   -0.002   -0.004   -0.002    0.000   -0.002    0.000    0.000   -0.005    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.19338376498028875}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 0.09355253342183101}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 0.8591335437952221}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.1745741995718021}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.26975749374548696}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.5675491117289342}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.023645611184699322}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 0.11694285135085192}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.034894129483006964}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.049562852297974136}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.09992913474301979}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.0804709270108857}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.4721007315609622}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.4971009105610804}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 0.9373573080054802}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.07869165370384468}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.11078527371686864}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.011981558188349617}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.01896914517553716}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.029718397699545207}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.0430274851462954}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.3078557703077741}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.1904672232500821}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.010943271603894189})])\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 3/10\n",
      "Time taken: 6 min, 17 s\n",
      "Projected time left: 14 min, 35 s\n",
      "Toy 3: Generating data...\n",
      "Toy 3: Data generation finished\n",
      "Toy 3: Loading data...\n",
      "Toy 3: Loading data finished\n",
      "Toy 3: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.978E+05               |    Ncalls=1302 (1302 total)    |\n",
      "| EDM = 3.73E-05 (Goal: 5E-06)  |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     True      |     True      |   False   |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   False   |   False   |  True  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297798.8924556829\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |   0.53    |   0.24    |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |   -5.63   |    0.23   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |   17.9    |    1.3    |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   0.98    |   0.17    |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   -6.28   |    0.13   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |   -2.4    |    0.8    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   0.919   |   0.023   |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |   0.51    |   0.29    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |   -2.18   |    0.04   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |   -0.30   |    0.41   |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |   0.04    |   0.10    |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   |   -0.30   |    0.07   |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |   -1.4    |    0.6    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |   4.10    |   0.21    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |    8.4    |    1.0    |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |   -0.95   |    0.06   |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   4.07    |   0.20    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   2.34    |   0.17    |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |  -6.283   |   0.018   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |   -4.72   |    0.04   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |   0.502   |   0.023   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |    0.7    |    0.4    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   1.11    |   0.20    |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |   4.41    |   0.12    |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000    0.013   -0.024    0.274    0.010    0.598   -0.001    0.010   -0.033   -0.039    0.439    0.141   -0.097    0.182   -0.018   -0.114    0.038    0.178    0.006    0.311    0.027    0.012    0.106    0.112 |\n",
      "|  omega_p |    0.013    1.000    0.036    0.003   -0.057   -0.013    0.001    0.031    0.015   -0.031    0.009    0.001   -0.077   -0.016    0.565   -0.011   -0.010    0.012    0.000    0.010    0.104    0.448    0.019   -0.009 |\n",
      "|    phi_s |   -0.024    0.036    1.000   -0.012   -0.012   -0.009    0.000    0.803    0.009   -0.001   -0.076    0.004   -0.013   -0.025    0.017    0.015   -0.016   -0.020    0.000   -0.033   -0.011   -0.046    0.006   -0.012 |\n",
      "|  p4040_s |    0.274    0.003   -0.012    1.000    0.003    0.006   -0.005   -0.003   -0.222   -0.024    0.095   -0.040    0.105   -0.150   -0.008   -0.025   -0.188    0.021   -0.003   -0.036    0.021    0.009    0.176   -0.415 |\n",
      "|    rho_p |    0.010   -0.057   -0.012    0.003    1.000    0.006   -0.000   -0.022   -0.008    0.004    0.025   -0.002    0.014    0.011    0.011    0.008    0.006    0.004   -0.000    0.001   -0.029   -0.062   -0.006    0.005 |\n",
      "|   Dbar_p |    0.598   -0.013   -0.009    0.006    0.006    1.000    0.009    0.005    0.357    0.066   -0.079    0.181   -0.038    0.320    0.003   -0.135    0.384   -0.036   -0.002    0.702   -0.032   -0.055   -0.076    0.559 |\n",
      "|  p3770_s |   -0.001    0.001    0.000   -0.005   -0.000    0.009    1.000    0.000    0.015   -0.006   -0.007    0.002   -0.014    0.003    0.000    0.007    0.010   -0.001   -0.002    0.011    0.001    0.002   -0.002    0.012 |\n",
      "|    phi_p |    0.010    0.031    0.803   -0.003   -0.022    0.005    0.000    1.000    0.003   -0.009   -0.011    0.001   -0.015   -0.004   -0.017    0.006   -0.003   -0.002    0.000   -0.012   -0.001   -0.028    0.001   -0.000 |\n",
      "|  psi2s_p |   -0.033    0.015    0.009   -0.222   -0.008    0.357    0.015    0.003    1.000   -0.296   -0.322    0.002   -0.425    0.114    0.004    0.214    0.303   -0.122   -0.008    0.425    0.024    0.035   -0.103    0.457 |\n",
      "| DDstar_s |   -0.039   -0.031   -0.001   -0.024    0.004    0.066   -0.006   -0.009   -0.296    1.000    0.136   -0.086    0.324   -0.015    0.005   -0.226   -0.016   -0.014   -0.003   -0.280   -0.055   -0.068   -0.074   -0.142 |\n",
      "|  bplus_2 |    0.439    0.009   -0.076    0.095    0.025   -0.079   -0.007   -0.011   -0.322    0.136    1.000   -0.070    0.388    0.230   -0.046   -0.497    0.006    0.143    0.002   -0.263    0.100    0.013   -0.174   -0.046 |\n",
      "|   Dbar_s |    0.141    0.001    0.004   -0.040   -0.002    0.181    0.002    0.001    0.002   -0.086   -0.070    1.000   -0.176    0.055    0.002    0.034    0.099   -0.016   -0.002    0.119   -0.003    0.002   -0.035    0.135 |\n",
      "| DDstar_p |   -0.097   -0.077   -0.013    0.105    0.014   -0.038   -0.014   -0.015   -0.425    0.324    0.388   -0.176    1.000    0.241    0.010   -0.637    0.072    0.055   -0.004   -0.310   -0.134   -0.183   -0.144   -0.047 |\n",
      "|  p4415_p |    0.182   -0.016   -0.025   -0.150    0.011    0.320    0.003   -0.004    0.114   -0.015    0.230    0.055    0.241    1.000   -0.009   -0.246    0.260   -0.108   -0.003    0.163   -0.022   -0.049   -0.211    0.403 |\n",
      "|  omega_s |   -0.018    0.565    0.017   -0.008    0.011    0.003    0.000   -0.017    0.004    0.005   -0.046    0.002    0.010   -0.009    1.000    0.005   -0.003   -0.014   -0.000   -0.016   -0.019   -0.109   -0.001   -0.001 |\n",
      "|  bplus_1 |   -0.114   -0.011    0.015   -0.025    0.008   -0.135    0.007    0.006    0.214   -0.226   -0.497    0.034   -0.637   -0.246    0.005    1.000   -0.122   -0.012    0.002    0.110   -0.610   -0.032    0.197   -0.098 |\n",
      "|  p4040_p |    0.038   -0.010   -0.016   -0.188    0.006    0.384    0.010   -0.003    0.303   -0.016    0.006    0.099    0.072    0.260   -0.003   -0.122    1.000    0.300   -0.005    0.249   -0.014   -0.034   -0.014    0.348 |\n",
      "|  p4160_s |    0.178    0.012   -0.020    0.021    0.004   -0.036   -0.001   -0.002   -0.122   -0.014    0.143   -0.016    0.055   -0.108   -0.014   -0.012    0.300    1.000   -0.001   -0.077    0.051    0.028    0.336   -0.136 |\n",
      "|  p3770_p |    0.006    0.000    0.000   -0.003   -0.000   -0.002   -0.002    0.000   -0.008   -0.003    0.002   -0.002   -0.004   -0.003   -0.000    0.002   -0.005   -0.001    1.000    0.001    0.000    0.001    0.001   -0.003 |\n",
      "|   jpsi_p |    0.311    0.010   -0.033   -0.036    0.001    0.702    0.011   -0.012    0.425   -0.280   -0.263    0.119   -0.310    0.163   -0.016    0.110    0.249   -0.077    0.001    1.000    0.057    0.025   -0.018    0.459 |\n",
      "|  bplus_0 |    0.027    0.104   -0.011    0.021   -0.029   -0.032    0.001   -0.001    0.024   -0.055    0.100   -0.003   -0.134   -0.022   -0.019   -0.610   -0.014    0.051    0.000    0.057    1.000    0.268    0.044   -0.020 |\n",
      "|    rho_s |    0.012    0.448   -0.046    0.009   -0.062   -0.055    0.002   -0.028    0.035   -0.068    0.013    0.002   -0.183   -0.049   -0.109   -0.032   -0.034    0.028    0.001    0.025    0.268    1.000    0.050   -0.035 |\n",
      "|  p4415_s |    0.106    0.019    0.006    0.176   -0.006   -0.076   -0.002    0.001   -0.103   -0.074   -0.174   -0.035   -0.144   -0.211   -0.001    0.197   -0.014    0.336    0.001   -0.018    0.044    0.050    1.000   -0.115 |\n",
      "|  p4160_p |    0.112   -0.009   -0.012   -0.415    0.005    0.559    0.012   -0.000    0.457   -0.142   -0.046    0.135   -0.047    0.403   -0.001   -0.098    0.348   -0.136   -0.003    0.459   -0.020   -0.035   -0.115    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.23746033613438966}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 0.22827854391421543}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 1.3242924020583473}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.17019599757914744}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.1332345481328261}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.7962709805312216}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.02308522313498257}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 0.2876950629670567}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.03769146402222612}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.40560761040792503}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.10018804096358758}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.06870386529805039}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.6338643163749156}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.2091282373871879}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 1.0387903247081258}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.05705684069783579}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.20352338541949067}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.1706073621912625}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.018358184856082094}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.03880632726508182}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.02288158646680416}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.3649741145457665}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.19896385950009027}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.11669689762358448})])\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 4/10\n",
      "Time taken: 8 min, 31 s\n",
      "Projected time left: 12 min, 42 s\n",
      "Toy 4: Generating data...\n",
      "Toy 4: Data generation finished\n",
      "Toy 4: Loading data...\n",
      "Toy 4: Loading data finished\n",
      "Toy 4: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.98E+05                |    Ncalls=1319 (1319 total)    |\n",
      "| EDM = 1.7E-05 (Goal: 5E-06)   |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     True      |     True      |   False   |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   False   |   False   |  True  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297957.0655247679\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |  -0.017   |   0.204   |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |    0.8    |    0.3    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |    18     |     6     |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   1.39    |   0.16    |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   -0.53   |    0.30   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |    1.1    |    0.3    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   1.62    |   0.22    |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |    -6     |     7     |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |   -2.03   |    0.03   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |  -0.300   |   0.029   |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |   -0.20   |    0.08   |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   |  -0.300   |   0.016   |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |    2.6    |    0.3    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |    4.6    |    0.5    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |    8.9    |    1.3    |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |   -0.93   |    0.03   |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   3.02    |   0.11    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   0.717   |   0.011   |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |   -1.85   |    0.15   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |  -1.752   |   0.032   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |   0.480   |   0.016   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |   1.32    |   0.31    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   0.40    |   0.18    |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |  -6.283   |   0.012   |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000   -0.019   -0.070    0.243    0.092   -0.423   -0.091   -0.067   -0.344    0.026    0.652    0.010    0.348    0.098   -0.006   -0.129    0.105   -0.002   -0.209    0.574   -0.023   -0.018    0.145   -0.003 |\n",
      "|  omega_p |   -0.019    1.000    0.306    0.001    0.031    0.007   -0.003    0.305    0.000   -0.001   -0.058    0.001    0.021   -0.017    0.720   -0.043   -0.015   -0.000    0.004   -0.040    0.058    0.489    0.002    0.000 |\n",
      "|    phi_s |   -0.070    0.306    1.000   -0.010   -0.245   -0.016   -0.014    0.995   -0.000    0.000   -0.203   -0.001    0.001   -0.050   -0.036   -0.019   -0.053   -0.001    0.015   -0.112   -0.001    0.145   -0.011    0.001 |\n",
      "|  p4040_s |    0.243    0.001   -0.010    1.000    0.011    0.079    0.084   -0.010   -0.204    0.006    0.013    0.014    0.052   -0.165   -0.001    0.072   -0.247   -0.006   -0.191    0.070    0.024    0.016    0.038    0.010 |\n",
      "|    rho_p |    0.092    0.031   -0.245    0.011    1.000    0.037    0.039   -0.250   -0.038    0.000    0.288    0.007    0.059    0.072    0.071    0.004    0.057    0.001   -0.027   -0.020    0.045    0.172    0.009   -0.001 |\n",
      "|   Dbar_p |   -0.423    0.007   -0.016    0.079    0.037    1.000   -0.073   -0.015    0.096   -0.000   -0.063    0.012   -0.099   -0.123    0.000    0.192   -0.316   -0.007   -0.017   -0.239    0.073    0.037   -0.006    0.004 |\n",
      "|  p3770_s |   -0.091   -0.003   -0.014    0.084    0.039   -0.073    1.000   -0.012    0.245   -0.012    0.018   -0.020   -0.109    0.012    0.000   -0.033    0.144   -0.001   -0.153   -0.118    0.021   -0.012    0.057   -0.002 |\n",
      "|    phi_p |   -0.067    0.305    0.995   -0.010   -0.250   -0.015   -0.012    1.000   -0.001    0.000   -0.194   -0.001    0.003   -0.048   -0.042   -0.020   -0.051   -0.001    0.014   -0.113    0.001    0.149   -0.011    0.001 |\n",
      "|  psi2s_p |   -0.344    0.000   -0.000   -0.204   -0.038    0.096    0.245   -0.001    1.000    0.017   -0.201    0.027    0.248    0.036   -0.000    0.091    0.149    0.007    0.238    0.025    0.021    0.017   -0.098   -0.000 |\n",
      "| DDstar_s |    0.026   -0.001    0.000    0.006    0.000   -0.000   -0.012    0.000    0.017    1.000    0.027   -0.001    0.039    0.015   -0.000   -0.030    0.041    0.001   -0.001    0.049   -0.009   -0.006    0.012   -0.001 |\n",
      "|  bplus_2 |    0.652   -0.058   -0.203    0.013    0.288   -0.063    0.018   -0.194   -0.201    0.027    1.000    0.016    0.059    0.282   -0.008   -0.237    0.067    0.002   -0.100    0.254   -0.031   -0.069   -0.100   -0.002 |\n",
      "|   Dbar_s |    0.010    0.001   -0.001    0.014    0.007    0.012   -0.020   -0.001    0.027   -0.001    0.016    1.000   -0.017   -0.007   -0.000    0.006   -0.003   -0.000    0.004    0.028    0.005    0.002    0.010    0.000 |\n",
      "| DDstar_p |    0.348    0.021    0.001    0.052    0.059   -0.099   -0.109    0.003    0.248    0.039    0.059   -0.017    1.000   -0.042   -0.000    0.234    0.160    0.006    0.110    0.588    0.082    0.056    0.069   -0.000 |\n",
      "|  p4415_p |    0.098   -0.017   -0.050   -0.165    0.072   -0.123    0.012   -0.048    0.036    0.015    0.282   -0.007   -0.042    1.000   -0.002   -0.126   -0.008    0.006    0.062    0.045   -0.010   -0.028   -0.095    0.003 |\n",
      "|  omega_s |   -0.006    0.720   -0.036   -0.001    0.071    0.000    0.000   -0.042   -0.000   -0.000   -0.008   -0.000   -0.000   -0.002    1.000    0.000   -0.002   -0.000    0.001   -0.009    0.001    0.104   -0.000    0.000 |\n",
      "|  bplus_1 |   -0.129   -0.043   -0.019    0.072    0.004    0.192   -0.033   -0.020    0.091   -0.030   -0.237    0.006    0.234   -0.126    0.000    1.000   -0.007   -0.002    0.014    0.022   -0.822   -0.127    0.125    0.000 |\n",
      "|  p4040_p |    0.105   -0.015   -0.053   -0.247    0.057   -0.316    0.144   -0.051    0.149    0.041    0.067   -0.003    0.160   -0.008   -0.002   -0.007    1.000   -0.005    0.178    0.095    0.020   -0.011    0.138   -0.003 |\n",
      "|  p4160_s |   -0.002   -0.000   -0.001   -0.006    0.001   -0.007   -0.001   -0.001    0.007    0.001    0.002   -0.000    0.006    0.006   -0.000   -0.002   -0.005    1.000    0.006    0.004   -0.000   -0.000    0.002   -0.001 |\n",
      "|  p3770_p |   -0.209    0.004    0.015   -0.191   -0.027   -0.017   -0.153    0.014    0.238   -0.001   -0.100    0.004    0.110    0.062    0.001    0.014    0.178    0.006    1.000    0.009    0.001    0.004   -0.065   -0.002 |\n",
      "|   jpsi_p |    0.574   -0.040   -0.112    0.070   -0.020   -0.239   -0.118   -0.113    0.025    0.049    0.254    0.028    0.588    0.045   -0.009    0.022    0.095    0.004    0.009    1.000   -0.014    0.008    0.031   -0.001 |\n",
      "|  bplus_0 |   -0.023    0.058   -0.001    0.024    0.045    0.073    0.021    0.001    0.021   -0.009   -0.031    0.005    0.082   -0.010    0.001   -0.822    0.020   -0.000    0.001   -0.014    1.000    0.196    0.041   -0.000 |\n",
      "|    rho_s |   -0.018    0.489    0.145    0.016    0.172    0.037   -0.012    0.149    0.017   -0.006   -0.069    0.002    0.056   -0.028    0.104   -0.127   -0.011   -0.000    0.004    0.008    0.196    1.000    0.021    0.000 |\n",
      "|  p4415_s |    0.145    0.002   -0.011    0.038    0.009   -0.006    0.057   -0.011   -0.098    0.012   -0.100    0.010    0.069   -0.095   -0.000    0.125    0.138    0.002   -0.065    0.031    0.041    0.021    1.000   -0.006 |\n",
      "|  p4160_p |   -0.003    0.000    0.001    0.010   -0.001    0.004   -0.002    0.001   -0.000   -0.001   -0.002    0.000   -0.000    0.003    0.000    0.000   -0.003   -0.001   -0.002   -0.001   -0.000    0.000   -0.006    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.20356743696537594}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 0.3225986671314396}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 5.601509675112658}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.1629207983331552}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.3048294644485092}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.3470474381158444}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.22325462827629838}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 7.06425982044021}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.03353790591960104}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.028565280512765218}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.08327237332492732}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.01630754606864146}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.34158714572292137}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.47996263337494316}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 1.2607557762653512}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.033473910114156125}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.11311693903164866}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.011131439035562851}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.1531198411151733}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.031531811376411145}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.016322193507993354}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.31173277656278864}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.18417991434425623}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.011534763945483029})])\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 5/10\n",
      "Time taken: 11 min, 1 \n",
      "Projected time left: 11 min\n",
      "Toy 5: Generating data...\n",
      "Toy 5: Data generation finished\n",
      "Toy 5: Loading data...\n",
      "Toy 5: Loading data finished\n",
      "Toy 5: Fitting pdf...\n",
      "------------------------------------------------------------------\n",
      "| FCN = 2.978E+05               |    Ncalls=1557 (1557 total)    |\n",
      "| EDM = 2.49E-05 (Goal: 5E-06)  |            up = 0.5            |\n",
      "------------------------------------------------------------------\n",
      "|  Valid Min.   | Valid Param.  | Above EDM | Reached call limit |\n",
      "------------------------------------------------------------------\n",
      "|     True      |     True      |   False   |       False        |\n",
      "------------------------------------------------------------------\n",
      "| Hesse failed  |   Has cov.    | Accurate  | Pos. def. | Forced |\n",
      "------------------------------------------------------------------\n",
      "|     False     |     True      |   True    |   True    | False  |\n",
      "------------------------------------------------------------------\n",
      "Function minimum: 297819.47573215386\n",
      "----------------------------------------------------------------------------------------------\n",
      "|   | Name     |   Value   | Hesse Err | Minos Err- | Minos Err+ | Limit-  | Limit+  | Fixed |\n",
      "----------------------------------------------------------------------------------------------\n",
      "| 0 | Ctt      |   1.07    |   0.18    |            |            |  -1.5   |   1.5   |       |\n",
      "| 1 | omega_p  |   -6.3    |    1.0    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 2 | phi_s    |   15.4    |    1.2    |            |            | 14.8182 | 23.5818 |       |\n",
      "| 3 | p4040_s  |   1.39    |   0.17    |            |            |0.00501244| 2.01499 |       |\n",
      "| 4 | rho_p    |   0.24    |   0.34    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 5 | Dbar_p   |   -5.4    |    0.3    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 6 | p3770_s  |   2.33    |   0.24    |            |            |0.918861 | 4.08114 |       |\n",
      "| 7 | phi_p    |   6.10    |   0.27    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 8 | psi2s_p  |  -2.080   |   0.031   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 9 | DDstar_s |  -0.300   |   0.026   |            |            |  -0.3   |   0.3   |       |\n",
      "| 10| bplus_2  |   -0.16   |    0.09   |            |            |   -2    |    2    |       |\n",
      "| 11| Dbar_s   |  -0.300   |   0.016   |            |            |  -0.3   |   0.3   |       |\n",
      "| 12| DDstar_p |   -5.4    |    0.5    |            |            |-6.28319 | 6.28319 |       |\n",
      "| 13| p4415_p  |   -2.31   |    0.16   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 14| omega_s  |    7.2    |    1.2    |            |            | 4.19232 | 9.40768 |       |\n",
      "| 15| bplus_1  |   0.90    |   0.07    |            |            |   -2    |    2    |       |\n",
      "| 16| p4040_p  |   -2.87   |    0.13   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 17| p4160_s  |   2.10    |   0.17    |            |            | 0.71676 | 3.68324 |       |\n",
      "| 18| p3770_p  |   -2.10   |    0.12   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 19| jpsi_p   |   4.591   |   0.025   |            |            |-6.28319 | 6.28319 |       |\n",
      "| 20| bplus_0  |   -0.45   |    0.03   |            |            |   -2    |    2    |       |\n",
      "| 21| rho_s    |    1.0    |    0.4    |            |            |0.0253049| 2.0747  |       |\n",
      "| 22| p4415_s  |   1.38    |   0.19    |            |            |0.126447 | 2.35355 |       |\n",
      "| 23| p4160_p  |   4.28    |   0.10    |            |            |-6.28319 | 6.28319 |       |\n",
      "----------------------------------------------------------------------------------------------\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|          |      Ctt  omega_p    phi_s  p4040_s    rho_p   Dbar_p  p3770_s    phi_p  psi2s_p DDstar_s  bplus_2   Dbar_s DDstar_p  p4415_p  omega_s  bplus_1  p4040_p  p4160_s  p3770_p   jpsi_p  bplus_0    rho_s  p4415_s  p4160_p |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "|      Ctt |    1.000    0.012   -0.033    0.262    0.086   -0.343   -0.022   -0.008   -0.410    0.010   -0.591    0.000   -0.126    0.105    0.017    0.070    0.019    0.200   -0.262    0.446    0.062    0.020    0.105   -0.117 |\n",
      "|  omega_p |    0.012    1.000    0.053   -0.001    0.145   -0.004    0.003    0.024   -0.003   -0.002   -0.025   -0.000   -0.025    0.014    0.642   -0.062    0.008    0.004   -0.000   -0.006    0.073   -0.076   -0.007    0.012 |\n",
      "|    phi_s |   -0.033    0.053    1.000   -0.001    0.139   -0.011   -0.017    0.738    0.016   -0.002    0.091   -0.002   -0.017   -0.015    0.070   -0.081   -0.012   -0.022    0.008    0.016    0.085    0.062    0.001   -0.001 |\n",
      "|  p4040_s |    0.262   -0.001   -0.001    1.000   -0.002    0.060    0.147   -0.001   -0.195    0.001    0.018    0.013   -0.125   -0.157   -0.002   -0.019   -0.176   -0.176   -0.171    0.054    0.004   -0.000    0.089   -0.368 |\n",
      "|    rho_p |    0.086    0.145    0.139   -0.002    1.000   -0.006    0.020    0.088   -0.019   -0.008   -0.166    0.000   -0.092    0.067    0.490   -0.267    0.039    0.028   -0.005   -0.003    0.311    0.071   -0.035    0.050 |\n",
      "|   Dbar_p |   -0.343   -0.004   -0.011    0.060   -0.006    1.000   -0.062   -0.003    0.007   -0.008   -0.026    0.014   -0.210   -0.124   -0.013   -0.023   -0.191    0.002   -0.105   -0.222   -0.030    0.029    0.051   -0.201 |\n",
      "|  p3770_s |   -0.022    0.003   -0.017    0.147    0.020   -0.062    1.000    0.001    0.260   -0.000   -0.006   -0.015    0.106   -0.023   -0.005   -0.017    0.042    0.084   -0.141   -0.089   -0.039    0.040    0.077   -0.072 |\n",
      "|    phi_p |   -0.008    0.024    0.738   -0.001    0.088   -0.003    0.001    1.000   -0.002   -0.000    0.002    0.000   -0.005    0.001    0.007   -0.007    0.001   -0.001   -0.000   -0.020    0.005    0.074   -0.000    0.001 |\n",
      "|  psi2s_p |   -0.410   -0.003    0.016   -0.195   -0.019    0.007    0.260   -0.002    1.000    0.019    0.125    0.026   -0.063    0.118    0.007    0.007    0.177   -0.061    0.243   -0.057    0.033   -0.045   -0.132    0.251 |\n",
      "| DDstar_s |    0.010   -0.002   -0.002    0.001   -0.008   -0.008   -0.000   -0.000    0.019    1.000    0.019   -0.001    0.043   -0.005   -0.006   -0.021   -0.001   -0.003    0.015    0.031   -0.012    0.010    0.005    0.011 |\n",
      "|  bplus_2 |   -0.591   -0.025    0.091    0.018   -0.166   -0.026   -0.006    0.002    0.125    0.019    1.000   -0.006    0.270   -0.265   -0.002   -0.370   -0.102   -0.082    0.077   -0.219    0.158   -0.164    0.189   -0.165 |\n",
      "|   Dbar_s |    0.000   -0.000   -0.002    0.013    0.000    0.014   -0.015    0.000    0.026   -0.001   -0.006    1.000   -0.003   -0.003   -0.002   -0.005   -0.005    0.002    0.005    0.028   -0.004    0.005    0.010   -0.004 |\n",
      "| DDstar_p |   -0.126   -0.025   -0.017   -0.125   -0.092   -0.210    0.106   -0.005   -0.063    0.043    0.270   -0.003    1.000   -0.338   -0.062   -0.177   -0.222   -0.098    0.073   -0.339   -0.127    0.084    0.009   -0.212 |\n",
      "|  p4415_p |    0.105    0.014   -0.015   -0.157    0.067   -0.124   -0.023    0.001    0.118   -0.005   -0.265   -0.003   -0.338    1.000    0.022    0.128    0.152   -0.085    0.086    0.170    0.011    0.004   -0.158    0.368 |\n",
      "|  omega_s |    0.017    0.642    0.070   -0.002    0.490   -0.013   -0.005    0.007    0.007   -0.006   -0.002   -0.002   -0.062    0.022    1.000   -0.185    0.010   -0.004    0.005    0.023    0.214   -0.208   -0.016    0.025 |\n",
      "|  bplus_1 |    0.070   -0.062   -0.081   -0.019   -0.267   -0.023   -0.017   -0.007    0.007   -0.021   -0.370   -0.005   -0.177    0.128   -0.185    1.000    0.091    0.008   -0.016    0.015   -0.924    0.227   -0.076    0.102 |\n",
      "|  p4040_p |    0.019    0.008   -0.012   -0.176    0.039   -0.191    0.042    0.001    0.177   -0.001   -0.102   -0.005   -0.222    0.152    0.010    0.091    1.000    0.341    0.187    0.070   -0.020    0.010    0.105    0.004 |\n",
      "|  p4160_s |    0.200    0.004   -0.022   -0.176    0.028    0.002    0.084   -0.001   -0.061   -0.003   -0.082    0.002   -0.098   -0.085   -0.004    0.008    0.341    1.000   -0.010    0.013   -0.035    0.040    0.331   -0.117 |\n",
      "|  p3770_p |   -0.262   -0.000    0.008   -0.171   -0.005   -0.105   -0.141   -0.000    0.243    0.015    0.077    0.005    0.073    0.086    0.005   -0.016    0.187   -0.010    1.000   -0.048    0.024   -0.022   -0.090    0.191 |\n",
      "|   jpsi_p |    0.446   -0.006    0.016    0.054   -0.003   -0.222   -0.089   -0.020   -0.057    0.031   -0.219    0.028   -0.339    0.170    0.023    0.015    0.070    0.013   -0.048    1.000    0.184   -0.115   -0.054    0.117 |\n",
      "|  bplus_0 |    0.062    0.073    0.085    0.004    0.311   -0.030   -0.039    0.005    0.033   -0.012    0.158   -0.004   -0.127    0.011    0.214   -0.924   -0.020   -0.035    0.024    0.184    1.000   -0.275   -0.026    0.038 |\n",
      "|    rho_s |    0.020   -0.076    0.062   -0.000    0.071    0.029    0.040    0.074   -0.045    0.010   -0.164    0.005    0.084    0.004   -0.208    0.227    0.010    0.040   -0.022   -0.115   -0.275    1.000    0.009   -0.023 |\n",
      "|  p4415_s |    0.105   -0.007    0.001    0.089   -0.035    0.051    0.077   -0.000   -0.132    0.005    0.189    0.010    0.009   -0.158   -0.016   -0.076    0.105    0.331   -0.090   -0.054   -0.026    0.009    1.000   -0.067 |\n",
      "|  p4160_p |   -0.117    0.012   -0.001   -0.368    0.050   -0.201   -0.072    0.001    0.251    0.011   -0.165   -0.004   -0.212    0.368    0.025    0.102    0.004   -0.117    0.191    0.117    0.038   -0.023   -0.067    1.000 |\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
      "Hesse errors: OrderedDict([(<zfit.Parameter 'Ctt' floating=True>, {'error': 0.17756068569158656}), (<zfit.Parameter 'omega_p' floating=True>, {'error': 0.974810495678307}), (<zfit.Parameter 'phi_s' floating=True>, {'error': 1.2048099140867805}), (<zfit.Parameter 'p4040_s' floating=True>, {'error': 0.1746148426443993}), (<zfit.Parameter 'rho_p' floating=True>, {'error': 0.3388255648988139}), (<zfit.Parameter 'Dbar_p' floating=True>, {'error': 0.33809702511271134}), (<zfit.Parameter 'p3770_s' floating=True>, {'error': 0.23753713150793}), (<zfit.Parameter 'phi_p' floating=True>, {'error': 0.2688131804807581}), (<zfit.Parameter 'psi2s_p' floating=True>, {'error': 0.031223719302520436}), (<zfit.Parameter 'DDstar_s' floating=True>, {'error': 0.025508419572656593}), (<zfit.Parameter 'bplus_2' floating=True>, {'error': 0.09027231656247292}), (<zfit.Parameter 'Dbar_s' floating=True>, {'error': 0.015661261261986303}), (<zfit.Parameter 'DDstar_p' floating=True>, {'error': 0.46545022942963765}), (<zfit.Parameter 'p4415_p' floating=True>, {'error': 0.1639990663169446}), (<zfit.Parameter 'omega_s' floating=True>, {'error': 1.248215198422895}), (<zfit.Parameter 'bplus_1' floating=True>, {'error': 0.06888293344375596}), (<zfit.Parameter 'p4040_p' floating=True>, {'error': 0.13379105407082625}), (<zfit.Parameter 'p4160_s' floating=True>, {'error': 0.17211479906979676}), (<zfit.Parameter 'p3770_p' floating=True>, {'error': 0.11563312352469346}), (<zfit.Parameter 'jpsi_p' floating=True>, {'error': 0.025354427208602814}), (<zfit.Parameter 'bplus_0' floating=True>, {'error': 0.03394709032162613}), (<zfit.Parameter 'rho_s' floating=True>, {'error': 0.354306489482555}), (<zfit.Parameter 'p4415_s' floating=True>, {'error': 0.19367485659297867}), (<zfit.Parameter 'p4160_p' floating=True>, {'error': 0.1046625681834934})])\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Toy 6/10\n",
      "Time taken: 13 min, 46 s\n",
      "Projected time left: 9 min, 8 s\n",
      "Toy 6: Generating data...\n",
      "Toy 6: Data generation finished\n",
      "Toy 6: Loading data...\n",
      "Toy 6: Loading data finished\n",
      "Toy 6: Fitting pdf...\n"
     ]
    }
   ],
>>>>>>> 8d1244c7112250ba6836e4821666773a87d7b683
   "source": [
    "# # zfit.run.numeric_checks = False   \n",
    "\n",
    "# fitting_range = 'cut'\n",
    "# total_BR = 1.7e-10 + 4.9e-10 + 2.5e-9 + 6.02e-5 + 4.97e-6 + 1.38e-9 + 4.2e-10 + 2.6e-9 + 6.1e-10 + 4.37e-7\n",
    "# cut_BR = 1.0 - (6.02e-5 + 4.97e-6)/total_BR\n",
    "\n",
    "# Ctt_list = []\n",
    "# Ctt_error_list = []\n",
    "\n",
<<<<<<< HEAD
    "# nr_of_toys = 1\n",
    "# if fitting_range == 'cut':\n",
    "#     nevents = int(pdg[\"number_of_decays\"]*cut_BR)\n",
    "# else:\n",
    "#     nevents = int(pdg[\"number_of_decays\"])\n",
    "# # nevents = pdg[\"number_of_decays\"]\n",
    "# event_stack = 1000000\n",
    "# # nevents *= 41\n",
    "# # zfit.settings.set_verbosity(10)\n",
    "# calls = int(nevents/event_stack + 1)\n",
=======
    "nr_of_toys = 25\n",
    "if fitting_range == 'cut':\n",
    "    nevents = int(pdg[\"number_of_decays\"]*cut_BR)\n",
    "else:\n",
    "    nevents = int(pdg[\"number_of_decays\"])\n",
    "# nevents = pdg[\"number_of_decays\"]\n",
    "event_stack = 1000000\n",
    "nevents *= 41\n",
    "# zfit.settings.set_verbosity(10)\n",
    "calls = int(nevents/event_stack + 1)\n",
>>>>>>> 8d1244c7112250ba6836e4821666773a87d7b683
    "\n",
    "# total_samp = []\n",
    "\n",
    "# start = time.time()\n",
    "\n",
    "# sampler = total_f.create_sampler(n=event_stack)\n",
    "\n",
    "# for toy in range(nr_of_toys):\n",
    "    \n",
    "#     ### Generate data\n",
    "    \n",
    "# #     clear_output(wait=True)\n",
    "    \n",
    "#     print(\"Toy {}: Generating data...\".format(toy))\n",
    "    \n",
    "#     dirName = 'data/zfit_toys/toy_{0}'.format(toy)\n",
    "    \n",
    "#     if not os.path.exists(dirName):\n",
    "#         os.mkdir(dirName)\n",
    "#         print(\"Directory \" , dirName ,  \" Created \")\n",
    "    \n",
    "#     reset_param_values()\n",
    "    \n",
    "#     if fitting_range == 'cut':\n",
    "        \n",
    "#         sampler.resample(n=nevents)\n",
    "#         s = sampler.unstack_x()\n",
    "#         sam = zfit.run(s)\n",
    "#         calls = 0\n",
    "#         c = 1\n",
    "        \n",
    "#     else:    \n",
    "#         for call in range(calls):\n",
    "\n",
    "#             sampler.resample(n=event_stack)\n",
    "#             s = sampler.unstack_x()\n",
    "#             sam = zfit.run(s)\n",
    "\n",
    "#             c = call + 1\n",
    "\n",
    "#             with open(\"data/zfit_toys/toy_{0}/{1}.pkl\".format(toy, call), \"wb\") as f:\n",
    "#                 pkl.dump(sam, f, pkl.HIGHEST_PROTOCOL)\n",
    "            \n",
    "#     print(\"Toy {}: Data generation finished\".format(toy))\n",
    "        \n",
    "#     ### Load data\n",
    "    \n",
    "#     print(\"Toy {}: Loading data...\".format(toy))\n",
    "    \n",
    "#     if fitting_range == 'cut':\n",
    "        \n",
    "#         total_samp = sam\n",
    "    \n",
    "#     else:\n",
    "                \n",
    "#         for call in range(calls):\n",
    "#             with open(r\"data/zfit_toys/toy_0/{}.pkl\".format(call), \"rb\") as input_file:\n",
    "#                 sam = pkl.load(input_file)\n",
    "#             total_samp = np.append(total_samp, sam)\n",
    "\n",
    "#         total_samp = total_samp.astype('float64')\n",
    "    \n",
    "#     if fitting_range == 'full':\n",
    "\n",
    "#         data = zfit.data.Data.from_numpy(array=total_samp[:int(nevents)], obs=obs)\n",
    "    \n",
    "#         print(\"Toy {}: Loading data finished\".format(toy))\n",
    "\n",
    "#         ### Fit data\n",
    "\n",
    "#         print(\"Toy {}: Fitting pdf...\".format(toy))\n",
    "\n",
    "#         for param in total_f.get_dependents():\n",
    "#             param.randomize()\n",
    "\n",
    "#         nll = zfit.loss.UnbinnedNLL(model=total_f, data=data, fit_range = (x_min, x_max), constraints = constraints)\n",
    "\n",
    "#         minimizer = zfit.minimize.MinuitMinimizer(verbosity = 5)\n",
    "#         # minimizer._use_tfgrad = False\n",
    "#         result = minimizer.minimize(nll)\n",
    "\n",
    "#         print(\"Toy {}: Fitting finished\".format(toy))\n",
    "\n",
    "#         print(\"Function minimum:\", result.fmin)\n",
    "#         print(\"Hesse errors:\", result.hesse())\n",
    "\n",
    "#         params = result.params\n",
    "#         Ctt_list.append(params[Ctt]['value'])\n",
    "#         Ctt_error_list.append(params[Ctt]['minuit_hesse']['error'])\n",
    "\n",
    "#         #plotting the result\n",
    "\n",
    "#         plotdirName = 'data/plots'.format(toy)\n",
    "\n",
    "#         if not os.path.exists(plotdirName):\n",
    "#             os.mkdir(plotdirName)\n",
    "# #             print(\"Directory \" , dirName ,  \" Created \")\n",
    "        \n",
    "#         probs = total_f.pdf(test_q, norm_range=False)\n",
    "#         calcs_test = zfit.run(probs)\n",
    "#         plt.clf()\n",
    "#         plt.plot(test_q, calcs_test, label = 'pdf')\n",
    "#         plt.legend()\n",
    "#         plt.ylim(0.0, 6e-6)\n",
    "#         plt.savefig(plotdirName + '/toy_fit_full_range{}.png'.format(toy))\n",
    "\n",
    "#         print(\"Toy {0}/{1}\".format(toy+1, nr_of_toys))\n",
    "#         print(\"Time taken: {}\".format(display_time(int(time.time() - start))))\n",
    "#         print(\"Projected time left: {}\".format(display_time(int((time.time() - start)/(c+calls*(toy))*((nr_of_toys-toy)*calls-c)))))\n",
    "    \n",
    "#     if fitting_range == 'cut':\n",
    "        \n",
    "#         _1 = np.where((total_samp >= x_min) & (total_samp <= (jpsi_mass - 60.)))\n",
    "        \n",
    "#         tot_sam_1 = total_samp[_1]\n",
    "    \n",
    "#         _2 = np.where((total_samp >= (jpsi_mass + 70.)) & (total_samp <= (psi2s_mass - 50.)))\n",
    "        \n",
    "#         tot_sam_2 = total_samp[_2]\n",
    "\n",
    "#         _3 = np.where((total_samp >= (psi2s_mass + 50.)) & (total_samp <= x_max))\n",
    "        \n",
    "#         tot_sam_3 = total_samp[_3]\n",
    "\n",
    "#         tot_sam = np.append(tot_sam_1, tot_sam_2)\n",
    "#         tot_sam = np.append(tot_sam, tot_sam_3)\n",
    "    \n",
    "#         data = zfit.data.Data.from_numpy(array=tot_sam[:int(nevents)], obs=obs_fit)\n",
    "        \n",
    "#         print(\"Toy {}: Loading data finished\".format(toy))\n",
    "        \n",
    "#         ### Fit data\n",
    "\n",
    "#         print(\"Toy {}: Fitting pdf...\".format(toy))\n",
    "\n",
    "#         for param in total_f_fit.get_dependents():\n",
    "#             param.randomize()\n",
    "\n",
    "#         nll = zfit.loss.UnbinnedNLL(model=total_f_fit, data=data, constraints = constraints)\n",
    "\n",
    "#         minimizer = zfit.minimize.MinuitMinimizer(verbosity = 5)\n",
    "#         # minimizer._use_tfgrad = False\n",
    "#         result = minimizer.minimize(nll)\n",
    "\n",
    "#         print(\"Function minimum:\", result.fmin)\n",
    "#         print(\"Hesse errors:\", result.hesse())\n",
    "\n",
    "#         params = result.params\n",
    "        \n",
    "#         if result.converged:\n",
    "#             Ctt_list.append(params[Ctt]['value'])\n",
    "#             Ctt_error_list.append(params[Ctt]['minuit_hesse']['error'])\n",
    "\n",
    "#         #plotting the result\n",
    "\n",
    "#         plotdirName = 'data/plots'.format(toy)\n",
    "\n",
    "#         if not os.path.exists(plotdirName):\n",
    "#             os.mkdir(plotdirName)\n",
    "#         #         print(\"Directory \" , dirName ,  \" Created \")\n",
    "        \n",
    "#         plt.clf()\n",
    "#         plt.hist(tot_sam, bins = int((x_max-x_min)/7.), label = 'toy data')\n",
    "#         plt.savefig(plotdirName + '/toy_histo_cut_region{}.png'.format(toy))\n",
    "\n",
    "        \n",
    "#         probs = total_f_fit.pdf(test_q, norm_range=False)\n",
    "#         calcs_test = zfit.run(probs)\n",
    "#         plt.clf()\n",
    "#         plt.plot(test_q, calcs_test, label = 'pdf')\n",
    "#         plt.axvline(x=jpsi_mass-60.,color='red', linewidth=0.7, linestyle = 'dotted')\n",
    "#         plt.axvline(x=jpsi_mass+70.,color='red', linewidth=0.7, linestyle = 'dotted')\n",
    "#         plt.axvline(x=psi2s_mass-50.,color='red', linewidth=0.7, linestyle = 'dotted')\n",
    "#         plt.axvline(x=psi2s_mass+50.,color='red', linewidth=0.7, linestyle = 'dotted')\n",
    "#         plt.legend()\n",
    "#         plt.ylim(0.0, 1.5e-6)\n",
    "#         plt.savefig(plotdirName + '/toy_fit_cut_region{}.png'.format(toy))\n",
    "        \n",
    "#         print(\"Toy {0}/{1}\".format(toy+1, nr_of_toys))\n",
    "#         print(\"Time taken: {}\".format(display_time(int(time.time() - start))))\n",
    "#         print(\"Projected time left: {}\".format(display_time(int((time.time() - start)/(toy+1))*((nr_of_toys-toy-1)))))\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'Ctt_list' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-40-c798fff4d6cd>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[1;32mwith\u001b[0m \u001b[0mopen\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"data/results/Ctt_list.pkl\"\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;34m\"wb\"\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;32mas\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 2\u001b[1;33m     \u001b[0mpkl\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mdump\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mCtt_list\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpkl\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mHIGHEST_PROTOCOL\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m      3\u001b[0m \u001b[1;32mwith\u001b[0m \u001b[0mopen\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"data/results/Ctt_error_list.pkl\"\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;34m\"wb\"\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;32mas\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      4\u001b[0m     \u001b[0mpkl\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mdump\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mCtt_error_list\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpkl\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mHIGHEST_PROTOCOL\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mNameError\u001b[0m: name 'Ctt_list' is not defined"
     ]
    }
   ],
   "source": [
    "# with open(\"data/results/Ctt_list.pkl\", \"wb\") as f:\n",
    "#     pkl.dump(Ctt_list, f, pkl.HIGHEST_PROTOCOL)\n",
    "# with open(\"data/results/Ctt_error_list.pkl\", \"wb\") as f:\n",
    "#     pkl.dump(Ctt_error_list, f, pkl.HIGHEST_PROTOCOL)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print('{0}/{1} fits converged'.format(len(Ctt_list), nr_of_toys))\n",
    "# print('Mean Ctt value = {}'.format(np.mean(Ctt_list)))\n",
    "# print('Mean Ctt error = {}'.format(np.mean(Ctt_error_list)))\n",
    "# print('95 Sensitivy = {}'.format(((2*np.mean(Ctt_error_list))**2)*4.2/1000))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plt.hist(tot_sam, bins = int((x_max-x_min)/7.))\n",
    "\n",
    "# plt.show()\n",
    "\n",
    "# # _ = np.where((total_samp >= x_min) & (total_samp <= (jpsi_mass - 50.)))\n",
    "\n",
    "# tot_sam.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Ctt.floating = False"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
<<<<<<< HEAD
   "metadata": {},
   "outputs": [],
   "source": [
    "# zfit.run(nll.value())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# result.fmin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# BR_steps = np.linspace(0.0, 1e-3, 11)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# CLS Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Step: 0/11\n",
      "Current Ctt: 0.0\n",
      "Ctt floating: True\n",
      "Toy 0/1\n"
     ]
    },
    {
     "ename": "RuntimeError",
     "evalue": "exception was raised in user function\nUser function arguments:\n     p4415_s = +0.738610\n     psi2s_p = +1.596704\n     p4415_p = +5.534902\n       phi_s = +14.820266\n    DDstar_s = +0.055634\n     p3770_s = +3.056326\n     p4160_s = +3.013799\n     omega_s = +8.221504\n     p3770_p = +3.227282\n       phi_p = -1.956867\n         Ctt = -0.074235\n    DDstar_p = -3.658252\n       rho_p = +3.220475\n      Dbar_p = +1.482929\n     p4040_p = -1.211697\n       rho_s = +0.426552\n      jpsi_p = -4.776017\n     omega_p = +2.899561\n     p4160_p = -1.589904\n     p4040_s = +1.161008\n      Dbar_s = +0.130791\nOriginal python exception in user function:\nKeyboardInterrupt: \n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 101, in func\n    loss_evaluated = self.sess.run(loss_val)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 929, in run\n    run_metadata_ptr)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1152, in _run\n    feed_dict_tensor, options, run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1328, in _do_run\n    run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1334, in _do_call\n    return fn(*args)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1319, in _run_fn\n    options, feed_dict, fetch_list, target_list, run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1407, in _call_tf_sessionrun\n    run_metadata)\n",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mRuntimeError\u001b[0m                              Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-47-feca68d05694>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m     76\u001b[0m             \u001b[0mminimizer\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mzfit\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mminimize\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mMinuitMinimizer\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mverbosity\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;36m5\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     77\u001b[0m             \u001b[1;31m# minimizer._use_tfgrad = False\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 78\u001b[1;33m             \u001b[0mresult\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mminimizer\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mminimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mnll\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     79\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     80\u001b[0m     \u001b[1;31m#         print(\"Function minimum:\", result.fmin)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\u001b[0m in \u001b[0;36mminimize\u001b[1;34m(self, loss, params)\u001b[0m\n\u001b[0;32m    205\u001b[0m             \u001b[0mtuple\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mstack\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0menter_context\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mparam\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mset_sess\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0msess\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0mparam\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    206\u001b[0m             \u001b[1;32mtry\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 207\u001b[1;33m                 \u001b[1;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_hook_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    208\u001b[0m             \u001b[1;32mexcept\u001b[0m \u001b[1;33m(\u001b[0m\u001b[0mFailMinimizeNaN\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mRuntimeError\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;32mas\u001b[0m \u001b[0merror\u001b[0m\u001b[1;33m:\u001b[0m  \u001b[1;31m# iminuit raises RuntimeError if user raises Error\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    209\u001b[0m                 \u001b[0mfail_result\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mstrategy\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mfit_result\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\u001b[0m in \u001b[0;36m_hook_minimize\u001b[1;34m(self, loss, params)\u001b[0m\n\u001b[0;32m    214\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    215\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0m_hook_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 216\u001b[1;33m         \u001b[1;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_call_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    217\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    218\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0m_call_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\u001b[0m in \u001b[0;36m_call_minimize\u001b[1;34m(self, loss, params)\u001b[0m\n\u001b[0;32m    218\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0m_call_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    219\u001b[0m         \u001b[1;32mtry\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 220\u001b[1;33m             \u001b[1;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_minimize\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mloss\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mparams\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mparams\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    221\u001b[0m         \u001b[1;32mexcept\u001b[0m \u001b[0mNotImplementedError\u001b[0m \u001b[1;32mas\u001b[0m \u001b[0merror\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    222\u001b[0m             \u001b[1;32mtry\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\u001b[0m in \u001b[0;36m_minimize\u001b[1;34m(self, loss, params)\u001b[0m\n\u001b[0;32m    136\u001b[0m             minimizer_setter)\n\u001b[0;32m    137\u001b[0m         \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_minuit_minimizer\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mminimizer\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 138\u001b[1;33m         \u001b[0mresult\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mminimizer\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mmigrad\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m**\u001b[0m\u001b[0mminimize_options\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    139\u001b[0m         \u001b[0mparams_result\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;33m[\u001b[0m\u001b[0mp_dict\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0mp_dict\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mresult\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    140\u001b[0m         \u001b[0mresult_vals\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;33m[\u001b[0m\u001b[0mres\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;34m\"value\"\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0mres\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mparams_result\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32miminuit\\_libiminuit.pyx\u001b[0m in \u001b[0;36miminuit._libiminuit.Minuit.migrad\u001b[1;34m()\u001b[0m\n",
      "\u001b[1;31mRuntimeError\u001b[0m: exception was raised in user function\nUser function arguments:\n     p4415_s = +0.738610\n     psi2s_p = +1.596704\n     p4415_p = +5.534902\n       phi_s = +14.820266\n    DDstar_s = +0.055634\n     p3770_s = +3.056326\n     p4160_s = +3.013799\n     omega_s = +8.221504\n     p3770_p = +3.227282\n       phi_p = -1.956867\n         Ctt = -0.074235\n    DDstar_p = -3.658252\n       rho_p = +3.220475\n      Dbar_p = +1.482929\n     p4040_p = -1.211697\n       rho_s = +0.426552\n      jpsi_p = -4.776017\n     omega_p = +2.899561\n     p4160_p = -1.589904\n     p4040_s = +1.161008\n      Dbar_s = +0.130791\nOriginal python exception in user function:\nKeyboardInterrupt: \n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 101, in func\n    loss_evaluated = self.sess.run(loss_val)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 929, in run\n    run_metadata_ptr)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1152, in _run\n    feed_dict_tensor, options, run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1328, in _do_run\n    run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1334, in _do_call\n    return fn(*args)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1319, in _run_fn\n    options, feed_dict, fetch_list, target_list, run_metadata)\n  File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\", line 1407, in _call_tf_sessionrun\n    run_metadata)\n"
     ]
    }
   ],
=======
   "metadata": {},
   "outputs": [],
>>>>>>> 8d1244c7112250ba6836e4821666773a87d7b683
   "source": [
    "# zfit.run.numeric_checks = False   \n",
    "\n",
    "fitting_range = 'cut'\n",
    "total_BR = 1.7e-10 + 4.9e-10 + 2.5e-9 + 6.02e-5 + 4.97e-6 + 1.38e-9 + 4.2e-10 + 2.6e-9 + 6.1e-10 + 4.37e-7\n",
    "cut_BR = 1.0 - (6.02e-5 + 4.97e-6)/total_BR\n",
    "\n",
    "Ctt_list = []\n",
    "Ctt_error_list = []\n",
    "\n",
    "nr_of_toys = 1\n",
    "nevents = int(pdg[\"number_of_decays\"]*cut_BR)\n",
    "# nevents = pdg[\"number_of_decays\"]\n",
    "event_stack = 1000000\n",
    "# nevents *= 41\n",
    "# zfit.settings.set_verbosity(10)\n",
    "\n",
    "mi = 0.0\n",
    "ma = 1e-3\n",
    "ste = 11\n",
    "\n",
    "BR_steps = np.linspace(mi, ma, ste)\n",
    "\n",
    "Ctt_steps = np.sqrt(BR_steps/4.2*1000)\n",
    "\n",
    "total_samp = []\n",
    "\n",
    "start = time.time()\n",
    "\n",
    "Nll_list = []\n",
    "\n",
    "sampler = total_f.create_sampler(n=nevents)\n",
    "\n",
    "__ = -1\n",
    "\n",
    "for Ctt_step in Ctt_steps:\n",
    "    \n",
    "    __ += 1\n",
    "    \n",
    "    for floaty in [True, False]:\n",
    "\n",
    "        Ctt_floating = floaty\n",
    "        \n",
    "        Nll_list.append([])\n",
    "\n",
    "        while len(Nll_list[-1]) <= nr_of_toys:\n",
    "            \n",
    "            print('Step: {0}/{1}'.format(__, ste))\n",
    "            \n",
    "            print('Current Ctt: {0}'.format(Ctt_step))\n",
    "            print('Ctt floating: {0}'.format(floaty))\n",
    "            \n",
    "            print('Toy {0}/{1}'.format(len(Nll_list[-1]), nr_of_toys))\n",
    "            \n",
    "            reset_param_values()\n",
    "            \n",
    "            if floaty:\n",
    "                Ctt.set_value(Ctt_step)\n",
    "            else:\n",
    "                Ctt.set_value(0.0)\n",
    "\n",
    "            sampler.resample(n=nevents)\n",
    "            s = sampler.unstack_x()\n",
    "            total_samp = zfit.run(s)\n",
    "            calls = 0\n",
    "            c = 1\n",
    "\n",
    "            data = zfit.data.Data.from_numpy(array=total_samp[:int(nevents)], obs=obs_fit)\n",
    "\n",
    "            ### Fit data\n",
    "\n",
    "            for param in total_f_fit.get_dependents():\n",
    "                param.randomize()\n",
    "\n",
    "            nll = zfit.loss.UnbinnedNLL(model=total_f_fit, data=data, constraints = constraints)\n",
    "\n",
    "            minimizer = zfit.minimize.MinuitMinimizer(verbosity = 5)\n",
    "            # minimizer._use_tfgrad = False\n",
    "            result = minimizer.minimize(nll)\n",
    "\n",
    "    #         print(\"Function minimum:\", result.fmin)\n",
    "    #         print(\"Hesse errors:\", result.hesse())\n",
    "\n",
    "            params = result.params\n",
    "\n",
    "            if result.converged:\n",
    "                Nll_list[-1].append(result.fmin)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dirName = 'data/CLs'\n",
    "\n",
    "CLs_values = []\n",
    "\n",
    "for i in range(len(Nll_list)/2):\n",
    "    CLs_values.append([])\n",
    "    for j in range(nr_of_toys):\n",
    "        CLs_values[i].append(Nll_list[i][j]-Nll_list[i+1][j])\n",
    "\n",
    "if not os.path.exists(dirName):\n",
    "    os.mkdir(dirName)\n",
    "    print(\"Directory \" , dirName ,  \" Created \")\n",
    "\n",
    "with open(\"'data/CLs/CLs_Nll_list.pkl\", \"wb\") as f:\n",
    "    pkl.dump(Nll_list, f, pkl.HIGHEST_PROTOCOL)\n",
    "\n",
    "with open(\"'data/CLs/CLs_list.pkl\", \"wb\") as f:\n",
    "    pkl.dump(CLs_values, f, pkl.HIGHEST_PROTOCOL)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " 1.0E-03\n"
     ]
    }
   ],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}