{
"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 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 = 1000000\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": 3,
"metadata": {},
"outputs": [],
"source": [
"def formfactor( q2, subscript): #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",
" b0 = ztf.constant(pdg[\"b0\"])\n",
" bplus = ztf.constant(pdg[\"bplus\"])\n",
" bT = ztf.constant(pdg[\"bT\"])\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",
"\n",
" for i in range(N):\n",
" _sum += b0[i]*(tf.pow(z,i))\n",
"\n",
" return tf.complex(prefactor * _sum, ztf.constant(0.0))\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",
" b = bT\n",
" else:\n",
" b = bplus\n",
"\n",
" for i in range(N):\n",
" _sum += b[i] * (tf.pow(z, i) - ((-1)**(i-N)) * (i/N) * tf.pow(z, N))\n",
"\n",
" return tf.complex(prefactor * _sum, ztf.constant(0.0))\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",
"def bifur_gauss(q, mean, sigma_L, sigma_R, scale):\n",
"\n",
" _exp = tf.where(q < mean, ztf.exp(- tf.pow((q-mean),2) / (2 * sigma_L**2)), ztf.exp(- tf.pow((q-mean),2) / (2 * sigma_R**2)))\n",
"\n",
" #Scale so the total area under curve is 1 and the top of the cusp is continuous\n",
"\n",
" dgamma = scale*_exp/(ztf.sqrt(2*np.pi))*2*(sigma_L*sigma_R)/(sigma_L+sigma_R)\n",
"\n",
" com = ztf.complex(dgamma, ztf.constant(0.0))\n",
"\n",
" return com\n",
"\n",
"def axiv_nonres(q):\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 = ztf.sqrt(tf.abs(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. * kabs**2. * beta**2. *tf.abs(tf.complex(C10eff, ztf.constant(0.0))*formfactor(q2, \"+\"))**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(tf.complex(C10eff, ztf.constant(0.0)) * formfactor(q2, \"0\")), 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 *ztf.sqrt(q2)\n",
"\n",
"def vec(q, funcs):\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 = ztf.sqrt(tf.abs(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 = kabs**2 * (1. - 1./3. * beta**2)\n",
"\n",
" abs_bracket = tf.abs(c9eff(q, funcs) * formfactor(q2, \"+\") + tf.complex(2.0 * C7eff * (mb + ms)/(mB + mK), ztf.constant(0.0)) * formfactor(q2, \"T\"))**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 * ztf.sqrt(q2)\n",
"\n",
"def c9eff(q, funcs):\n",
"\n",
" C9eff_nr = tf.complex(ztf.constant(pdg['C9eff']), ztf.constant(0.0))\n",
"\n",
" c9 = C9eff_nr\n",
"\n",
" c9 = c9 + 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(-q)))\n",
" \n",
" big_bracket = tf.where(y > ztf.const(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) - 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) - 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(zfit.pdf.ZPDF):\n",
" _N_OBS = 1 # dimension, can be omitted\n",
" _PARAMS = ['jpsi_mass', 'jpsi_scale', 'jpsi_phase', 'jpsi_width',\n",
" 'psi2s_mass', 'psi2s_scale', 'psi2s_phase', 'psi2s_width'#,\n",
" #'cusp_mass', 'sigma_L', 'sigma_R', 'cusp_scale'\n",
" ] # the name of the parameters\n",
"\n",
" def _unnormalized_pdf(self, x):\n",
" \n",
" x = x.unstack_x()\n",
"\n",
" def jpsi_res(q):\n",
" return resonance(q, _mass = self.params['jpsi_mass'], scale = self.params['jpsi_scale'], phase = self.params['jpsi_phase'], width = self.params['jpsi_width'])\n",
"\n",
" def psi2s_res(q):\n",
" return resonance(q, _mass = self.params['psi2s_mass'], scale = self.params['psi2s_scale'], phase = self.params['psi2s_phase'], width = self.params['psi2s_width'])\n",
"\n",
" def cusp(q):\n",
" return bifur_gauss(q, mean = self.params['cusp_mass'], sigma_L = self.params['sigma_L'], sigma_R = self.params['sigma_R'], scale = self.params['cusp_scale'])\n",
"\n",
" funcs = jpsi_res(x) + psi2s_res(x) #+ cusp(x)\n",
"\n",
" vec_f = vec(x, funcs)\n",
"\n",
" axiv_nr = axiv_nonres(x)\n",
"\n",
" tot = vec_f + axiv_nr\n",
"\n",
" return tot"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Load data"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"x_min = 2*pdg['muon_M']\n",
"x_max = (pdg[\"Bplus_M\"]-pdg[\"Ks_M\"]-0.1)\n",
"\n",
"obs = zfit.Space('q', limits = (x_min, x_max))\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 parameters"
]
},
{
"cell_type": "code",
"execution_count": 7,
"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": [
"#jpsi\n",
"\n",
"jpsi_mass, jpsi_width, jpsi_phase, jpsi_scale = pdg[\"jpsi\"]\n",
"# jpsi_scale *= pdg[\"factor_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))\n",
"jpsi_s = zfit.Parameter(\"jpsi_s\", ztf.constant(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))\n",
"psi2s_s = zfit.Parameter(\"psi2s_s\", ztf.constant(psi2s_scale))\n",
"\n",
"#cusp\n",
"\n",
"# cusp_mass, sigma_R, sigma_L, cusp_scale = 3550, 3e-7, 200, 0\n",
"\n",
"# cusp_m = zfit.Parameter(\"cusp_m\", ztf.constant(cusp_mass), floating = False)\n",
"# sig_L = zfit.Parameter(\"sig_L\", ztf.constant(sigma_L), floating = False)\n",
"# sig_R = zfit.Parameter(\"sig_R\", ztf.constant(sigma_R), floating = False)\n",
"# cusp_s = zfit.Parameter(\"cusp_s\", ztf.constant(cusp_scale), floating = False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Setup pdf"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"total_f = total_pdf(obs=obs, 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",
" #cusp_mass = cusp_m, sigma_L = sig_L, sigma_R = sig_R, cusp_scale = cusp_s) \n",
" \n",
"# print(total_pdf.obs)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Test if graphs actually work and compute values"
]
},
{
"cell_type": "code",
"execution_count": 9,
"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, 2000000)\n",
"\n",
"probs = total_f.pdf(test_q)\n",
"\n",
"calcs_test = zfit.run(probs)\n",
"res_y = zfit.run(jpsi_res(test_q))\n",
"f0_y = zfit.run(formfactor(test_q,\"0\"))\n",
"fplus_y = zfit.run(formfactor(test_q,\"+\"))\n",
"fT_y = zfit.run(formfactor(test_q,\"T\"))"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAD8CAYAAABthzNFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3X90XOV95/H3VyPNjH7/lmxLsiWwYzA/Awq0Nd24lKSmkNBtaQIN2zbQpUlLusmes7tku21O9+wuJ7tpTncbTsAbKCGbQNikbVxCQrIljpPWAYuYEBvbIBsby7ItWZL1WyON9Owfc0eWpZEsazSaO1ef1zk6mvt4NPOdx9L9zPPc594x5xwiIiKz5WW7ABER8ScFhIiIpKSAEBGRlBQQIiKSkgJCRERSUkCIiEhKCggREUlJASEiIikpIEREJKX8bBewkJqaGtfc3JztMkREcsqrr7561jlXm+7j+DogmpubaWtry3YZIiI5xcyOL8fjaIpJRERSWrGAMLPLzOwJM/vGSj2niIgsXVoBYWZPmlmXme2f1b7dzA6bWbuZPQzgnDvqnHsgnecTEZGVk+4xiKeALwBPJxvMLAQ8CrwP6AD2mtlO59wbaT6XiEhaJiYm6OjoYGxsLNulLItoNEpjYyMFBQUZefy0AsI5t9vMmmc13wS0O+eOApjZs8BdgAJCRLKqo6OD0tJSmpubMbNsl5MW5xw9PT10dHTQ0tKSkefIxDGIBuDEjO0OoMHMqs3sMeDdZvbp+X7YzB40szYza+vu7s5AeSKyWo2NjVFdXZ3z4QBgZlRXV2d0NJSJZa6pet4553qAj13sh51zO4AdAK2trfq4OxFZVkEIh6RMv5ZMjCA6gKYZ241AZwaeR0RS6Dw3yg8OdWW7DAmATATEXmCTmbWYWRi4B9iZgecRkRTu/Osf89Gn9ma7DEnDrl27uPPOOwGIxWLcdtttXH/99Xz9619f0TrSmmIys2eAbUCNmXUAn3HOPWFmDwEvAiHgSefcgbQrFZFF6R0ez3YJsoz27dvHxMQEr7322oo/d7qrmO6dp/0F4IV0HltEJIiOHTvG9u3bufnmm9m3bx/vete7ePrpp9m9ezef/OQnqamp4YYbbgCgq6uL++67j+7ubq6//nq++c1vcvnll69Yrb6+FpOISKb8xT8c4I3OgWV9zC3ryvjMB6666P0OHz7ME088wdatW7n//vv5/Oc/z+OPP85LL73Exo0b+fCHPwxAXV0dX/rSl/jc5z7H888/v6y1LoauxSQissKamprYunUrAPfddx9tbW20tLSwadMmzIz77rsvyxUmaAQhIqvSYt7pZ8rs5an9/f2+XH6rEYSIyAp755132LNnDwDPPPMMt912G2+//TZHjhyZbvMDBYSIyAq78sor+fKXv8y1115Lb28vn/rUp9ixYwd33HEHt9xyCxs2bMh2iYCmmEREVlxeXh6PPfbYBW3bt2/n0KFDc+67bds2tm3btkKVXUgjCBERSUkBISKygpqbm9m/f//F7+gDCggRWVWcC841QDP9WhQQIrJqRKNRenp6AhESyc+DiEajGXsOHaQWkVWjsbGRjo4OgvJZM8lPlMsUBYSIrBoFBQUZ+/S1INIUk4iIpKSAEBGRlBQQIiKSkgJCRERSUkCIiEhKCggREUlJASEiIikpIEREJCUFhIiIpKSAEBGRlBQQIiKSkgJCJKCCcMVSyS4FhIiIpKSAEBGRlBQQIgHlxxmmfe/00fzwt3nl7d5slyKLsGIBYWaXmdkTZvaNlXpOEfGXf2o/C8AP3+zKciWyGIsKCDN70sy6zGz/rPbtZnbYzNrN7OGFHsM5d9Q590A6xYpIMPhxdCNzLfYT5Z4CvgA8nWwwsxDwKPA+oAPYa2Y7gRDwyKyfv985p7cMIitI+2BJ16ICwjm328yaZzXfBLQ7544CmNmzwF3OuUeAO5ezSBERWXnpHINoAE7M2O7w2lIys2ozewx4t5l9eoH7PWhmbWbWFpQPFheRBDMDNLrJFYudYkrFUrTN+//unOsBPnaxB3XO7QB2ALS2tur3SGSJEifKpfozFVmcdEYQHUDTjO1GoDO9ckRExC/SCYi9wCYzazGzMHAPsHN5yhIRkWxb7DLXZ4A9wGYz6zCzB5xzceAh4EXgIPCcc+5A5koVkUvh5/lZTXzlhsWuYrp3nvYXgBeWtSIRCTw/h5ecp0ttiIhISgoIkYDy89nKmmLKDQoIEVkx+oyK3KKAEAko5+OZftMQIicoIERkxWgAkVsUECIikpICQiSg/PhuXVNLuUUBISIrxo+hJfNTQIjIijMtdM0JCggRWTEaQOQWBYRIQPl5OkfHInKDAkIkoCZ9nBA+Lk1mUECIBNSU9sKSJgWESEBNTfk3IDTFlBsUECIB5cd80KAmtyggRALKz1NMGkDkBgWESED5cYrJzxcQlLkUECIB5cN8kByjgBAJKD9OMekM6tyigBAJqEkfDiE0xZRbFBAiAeXDAYTkGAWESED5eYrJf5VJKgoIkYDy46U28rxDED4sTVJQQIgElPPhXjjPSwg/jm5kLgWESEBNTmW7grmSl9jw4fFzSUEBIRJQfnyXPn0Mwoe1yVwKCJGA8mVAJI9BZLcMWaQVCwgzu9LMHjOzb5jZx1fqeUVWqykfTjElD1L78TIgMteiAsLMnjSzLjPbP6t9u5kdNrN2M3t4ocdwzh10zn0M+BDQuvSSRWQx/DiCkNyy2BHEU8D2mQ1mFgIeBW4HtgD3mtkWM7vGzJ6f9VXn/cwHgR8D/7hsr0BEUvJzQOjzIHJD/mLu5JzbbWbNs5pvAtqdc0cBzOxZ4C7n3CPAnfM8zk5gp5l9G/jaUosWkYvzY0D4sCRZwKICYh4NwIkZ2x3AzfPd2cy2Ab8JRIAXFrjfg8CDAOvXr0+jPJHVzY/T/MmSTEOInJBOQKT6H573V9I5twvYdbEHdc7tAHYAtLa2+vBXXCQ3xCf9++ejeMgN6axi6gCaZmw3Ap3plSMiy2Xch2fKaYopt6QTEHuBTWbWYmZh4B5g5/KUJSLpmoj7LyCmaQiRExa7zPUZYA+w2cw6zOwB51wceAh4ETgIPOecO5C5UkXkUvhyBKFT5HLKYlcx3TtP+wsscMBZRLJnwocBkaRPlssNutSGSEDFfDjFpGMQuUUBIRJQ4z4MiCStcs0NCgiRAJl5lVQ/TjElPydb+ZAbFBAiATIx49wHP44gJv149p7MSwEhEiCx+OT0bT+OIJKX/1BO5AYFhEiAzBw1+HkEMenHa5HLHAoIkQCZee5DzIcjiElvBDHh48uAyHkKCJEAiU34ewQxEU8Egx+X4MpcCgiRABkej0/fHh2fXOCe2TE+majJj+ElcykgRAJkOHY+FIZi8QXumR3JYPDjZUBkLgWESIAMzwiFET+OIJIBEfdfbTKXAkIkQJKjhurisD9HEN7IQQepc4MCQiRAkiOI2tLIBaMJvzg/gtAUUy5QQIgESHLUUF8W9WdAeCMHBURuUECIBEjyIHVdaYRhXx6DSNTkx3M0ZC4FhEiADI5NUFgQorywwJ8jCE0x5RQFhEiA9A6PU1UcpiSaz8j4JHGfvVNPHqTWKqbcoIAQCZCe4XGqS8JUFYcB6BuZyHJFFxrzzvQem/BXcElqCgiRAEmOIJIB0Ts8nuWKLjTiTXv5cQmuzKWAEAmQZEBUF0cA6BmKZbmiCw3NCAinzx/1PQWESEA45+gZjlFdHKa6JDGC6PHRCMI5x/D4JKE8Y3LKaZopByggRAJiYDTO2MQU9WVRX04xxeJTTE451pRFARiM+ev4iMylgBAJiBN9IwA0VBRSWRTGzF9TTMllt3VliemvwTEdh/A7BYRIQJw8NwpAY2URoTyjrjTCqf6xLFd1XvIkvuQIYkgB4XsKCJGAONmXCIiGykIgERTJUYUfJKeU1pYn6hsY0xST3ykgRALind4RisMhKosKAGiqLKTDCw0/SB4PeVd9CQBnfTT9JakpIEQC4tDpATbVl2JmADRVFXGqf8w3Z1MnA+KKtWUAdA8qIPxuxQLCzLaZ2Y/M7DEz27ZSzyuyGjjnOHx6kCvWlE63NVYWMjnl6Dznj+MQPUOJgNhQVURhQYiuAQWE3y0qIMzsSTPrMrP9s9q3m9lhM2s3s4cv8jAOGAKiQMfSyhWRVLoGY/SNTLB5RkBsrEvcPnxmMFtlXaB3eJxQnlFeWEBdWYQujSB8b7EjiKeA7TMbzCwEPArcDmwB7jWzLWZ2jZk9P+urDviRc+524D8Af7F8L0FE2o71AXB9U8V02xVrSjGDg6cGslXWBboHY1QVh8nLM+pLo5z20QorSS1/MXdyzu02s+ZZzTcB7c65owBm9ixwl3PuEeDOBR6uD4jM949m9iDwIMD69esXU57IqvfK2z0UhUNc3VA+3VYcyae5upg3Ov0RECf6RmjyVlhdVlvM9984k+WK5GLSOQbRAJyYsd3htaVkZr9pZo8DXwG+MN/9nHM7nHOtzrnW2traNMoTWR2cc/y4/Sw3bqikIHThn/RV68p47cQ5X1z36J3eEZqqigDYWFdCz/C4r870lrnSCQhL0Tbvb6Fz7m+dc3/onPuwc25XGs8rIjO8eWaII93DvP+qNXP+7Rcvr+b0wBhHuoezUNl54/EpOs+N0lSZCIhN9d7xkdP+OD4iqaUTEB1A04ztRqAzvXJE5FL93b6T5BlsTxEQv7wxMQr/8VvdK13WBd48M8iUgyvWJoLhusZyzODlt3uyWpcsLJ2A2AtsMrMWMwsD9wA7l6csEVmM4Vicr718nO1Xr6G2dO6hvfXVRWysK+H5109lobrzXu/oB+DahsRB9IqiMFevK2f3m9kNLlnYYpe5PgPsATabWYeZPeCciwMPAS8CB4HnnHMHMleqiMz2xV1HGBiL869/+bJ573P3jY20He+jvWtoBSu70D+1n6W2NEJTVeF02/ar1/DTd85ltS5Z2KICwjl3r3NurXOuwDnX6Jx7wmt/wTn3Lufc5c65/5rZUkVkptdOnOPx3Uf4l+9u4N3rK+e932/d0Ei0II+/fumtFazuvKFYnF2Hu7jtyrrps7wBPtTaRCQ/j7/83mFfHESXuXSpDZEc9NaZQR58uo36sih/dueWBe9bWxrh/q0tfOu1Tn6YhSmdp/ccY3h8kg+1Nl3QXlsa4U9+dRPf2X+az33vMJNTCgm/UUCI5BDnHM+/3sndj+1hysGTv/+e6Q8HWsgnbt3E5vpSPvG1n/Lq8b4VqDThx2+d5a/+31vcdmV9ylHOx997OR9qbeTRHxzh9v+5my/96Cj7T/YzNjG5YjXK/MzPQ7vW1lbX1taW7TJEsi4+OcXut7p5bNdRXjnWy3WN5Xzhd26YPq9gMU70jnDfEy9zsm+Uj25t5ve3ttBQUXjxH7xEvcPjvHaij7/f18k/vN7J5vpSvvoHN1Ndkvr8WOcc3/75Kb646wgHZpzUV18Woao4Qlk0n9JoAcWREAWhPApCeYRDlridn0dBXuJ2OD9v+ns4P4/wrLaCkBHxtiP5IUqi+ZRG8ykJ55OXl2rVfu4ys1edc61pP44CQsR/nHOcPDfKq8f72HOkh++9cYbe4XHWlkf52Hsv5yM3ryc/dOkTAH3D4/yXbx/k7/Z1MOUSl96+pqGC5uoiaksjlETzKY7kk59nzNw1xKemGJuYYnR8ktGJScYmJhkdn2QoFqdneJyeoRi9w+Oc6h+bvsZSaTSfe29azydu3UhptGBR9Z08N0rbsV6O94zwTu8I50YmGBibYGB0gtGJSeKTjvHJKSYmp5iITzHhbafDDErCibAojRZQXlRAXWmEutIodWUR6koj1JdFWV9VxLqKQkI5ECYKCJEAGJuYpGsgxom+EY50D3GkK3HS25tnBqd3tCWRfH7lijruuGYNt15RTzg//ZnhE70j/MPrnew50sObZwY5s8Qrq4ZDeVSXhKkuCVNVnNiZbq4vZcu6MlqbK4nkh9Ku9WKcc8SnHOPxRHCMx6e8EDnfFpvxbxOTibAbik0wOBZnYCzO4Fji9uDYBH0jE3QPxugaGGN4/MKprnB+Hs3VRbTUFLN5TRnXN5VzbWMFNfOMjrJFASHiU5NTjv7RCfpGxunzLifRNzJO10CMUwNjnOkf41T/GKcHxuZcaqIkks/ltcVcXlvCdU0V3LihkivWlC5ptHApRscn6R0ZZzgWZygWZ8o7YJxcdJRnRlE4n8KCENFwHtGCEIUFoTmX9giaoVicroHE/9XxnhHePjvM0e5hjp4d4tjZYZLH1RsrC7llYw3bNtfySxtrKFvkiClTFBAiK2DK29n3ejv7vpGJxE5/enuc3mEvDLy2c6MTzPdnVVUcZk1ZlLXlUerLo6wti7KmPEpDRSGX15VQVxq5YCmo+NdwLM7+k/283tFP2/Fe/rm9h8FYnPw849Yr6vitGxv5lc11yzLiu1QKCJFLNDXznX1yx57cyc8TAP2jE8y3+jKcn0dVUZjK4jBVxQVUFIXPbxcVUFkcprIoTFVxmIqiAmpKIkQLMj/lItkxMTnFvnfO8b0Dp/n71zo5OxRjTVmUP3zvZdzznvUUhlfu/14BIatayp19ckpnZJxzwxMXvMvvG5ng3Mj4/Dv7UB6VxQXTO/TEzr1gRgCEZwRAAVXFYQoLQnq3LynFJ6fYdbibHT86yitv91JfFuHP7tzCHdesXZHfGQWEBMbklGNgvmmc6Xn8xA6+d2Scc5ews68oKqC6ODK9PTsAkttFYe3sJTNePtrDf37+DQ50DnDblfX85W9fR3lRZo9RKCDEt+KTU/SNTNAzHKNnaJyz3hLInqFxeoZjnB06f+D2YnP2M3f253fuF25XFCXe0Vd67/aLtbMXn4lPTvE3/3SMz373EOsqCnny99/DxrqSjD2fAkJWjHOOgbE4PUOx6TXvZ4cSO/ze4RhnvbZEACR2/Kl+rUJ5RmVRmOrixLLI2dM4lUV6Zy/B9urxPv7wK68C8NU/uPmCzxBfTgoISdvE5BRnh2J0DcToGowl1n4PJk506hqI0e3dPjsUY2Iy9e9JWTSfmpJIYi18sfe9JELNjO3k7fLCgsCdsSpyqdq7hvid//0TQnnGtx7aSl1pdNmfY7kCYlGfSS25JT45xZnBGKf7Rzndf+FOv2twjG4vDHrneadfVRymrjRCbWmEy+tKqC2NUDsrBGpKIlQWhbOyhE8kl22sK+FvPvoe7v7iHv74qz/l2Qd/0bdnZysgcszE5BRnBsY43T9GZ/8Yp/tHOdU/xqlzY5waSGx3D8bmHMDNzzNqSxNnujZWFnHDhkpqSyLepQSiiUsLlEWoLo5opy+SYVetK+e//ebVfOrrP+Opfz7GA7e0ZLuklBQQPjMyHqejb5SOvhFO9Ca+d/SN0nlulM7+Mc4Oxea86y8Oh1hbUcja8iib62tZU564vaY8ypqyKPVlUSo0vSPiK79xfQM7X+vkcy8e5gPXrc3IVFO6FBArbDw+RUdf4kJkHX2jnPACoMPb7pl16YVIfh6NlYWsqyjkijVlrK2Iejv/RAisLY8u+kJoIuIfZsaff+Aqbvv8D3n8h0cv+rke2aCAyIDJKcfJvlHe7hnm2Nlh3va+jvUM09E3esEHo4RDeTRUFtJYWcj715XTWFlIU1VR4ntlETUlYa3iEQmolppifuP6Bv7PT47zJ7duyvj5EZdKAZGG0fFJ2ruGePPMIG92DXKka4i3zw7zTu/IBat+isMhmmuKuaahnA9et47m6mI2VBfRWFlEXWlEUz8iq9hHtzbzzZ928K2fneR3f7E52+VcQAGxCOPxqUQInBnkzTNDtHclvp/oG5k+HhAO5dFSU8ymulLet2UNLTVFtNSU0FxTRG2JLsAmIqld3VDO1Q1lPNd2QgHhd6Pjkxw8PcCBk/3sPznAgVP9HD49OD0iKAgZLTXFXNNYzm/d0Mi76kvYVF9Kc3VRxi/JLCLB9IFr1/HIdw5xqn+UteXL/yl/S7XqA+LMwBhtx/rYe6yXtuO9vNE5ML1EtLKogKsbynnglsu4al0ZV64tZUN1ceCvgS8iK+vWK+p45DuHeOlQFx+5eUO2y5m26gLCOcerx/v4/sEz/OBQF2+eGQIgWpDHu5sq+aNtG7mmsZyrG8pZVx7V1JCIZNzGuhIaKgr58VtnFRDZMDIe5yt7jvPMK+9wrGeE/DzjppYq7r6xkZtaqrlqXZlGBiKSFWbGjRsq2XusN9ulXGBVBMR395/mP/39fs4OxbippYpP3LqJ911Vn/WPBRQRSbq+qYKdP+vkdP8Ya8r9cdJc4APii7uO8NnvHuKahnIe/1c3cOOGqmyXJCIyx3VNFQD8/GT/6gsIM/tl4CPec25xzv1Spp/zxQOn+ex3D/HB69bxud++TtcYEhHfSn4+xJHuId5HfZarSVjUHtPMnjSzLjPbP6t9u5kdNrN2M3t4ocdwzv3IOfcx4Hngy0sveXHG41N85lsH2LK2jL/8kMJBRPytvLCAmpIwR7uHsl3KtMXuNZ8Cts9sMLMQ8ChwO7AFuNfMtpjZNWb2/Kyvuhk/+jvAM8tQ+4JePHCa0wNj/Lvtm3XwWURywmW1JRztHs52GdMWNcXknNttZs2zmm8C2p1zRwHM7FngLufcI8CdqR7HzNYD/c65gSVXvEj/ePAMNSUR3rupNtNPJSKyLDZUFbH7re5slzEtnbfWDcCJGdsdXttCHgD+ZqE7mNmDZtZmZm3d3UvvqJ8c7eUXLqvSdY5EJGesKY/SPRi74IKe2ZROQKTa8y74qpxzn3HO/fNF7rPDOdfqnGutrV3au/+hWJzTA2NsWVe2pJ8XEcmGurIoUw56hmLZLgVILyA6gKYZ241AZ3rlLI93ekYA2FBVnOVKREQWr740AsCZgdwPiL3AJjNrMbMwcA+wc3nKSs+p/lEAGir9c9ErEZGLqS9LnP9wZmAsy5UkLHaZ6zPAHmCzmXWY2QPOuTjwEPAicBB4zjl3IHOlLt5QLA5AWTTw5wGKSIBUFYcBODc6keVKEha7iuneedpfAF5Y1oqWweBYIiBKFBAikkOSl/8Z8ElABPIEgeQIojSiay2JSO5IvqkdGFNAZMzQWJw8S1zCW0QkV4TyjNJIPgOj8WyXAgQ0IOJTjvxQnj7LQURyTllhgUYQIiIyV2k0X8cgMsktfL6eiIhvFYZDjE5MZrsMIKABAalP8xYR8btIfh6xialslwEEOCBERHJRtCBELK4RROZohklEclQ0P8SYRhCZpQVMIpKLIgV5jGkEISIis0XzQzoGkUmaYRKRXKURxAowrWMSkRwULQgxpmWuIiIyW0HIiE/6Yx4kkAHhnD86V0TkUoXMiAfgI0d9TauYRCQXhfISu+UpH4REYANCRCQXhby9sh9GEYEMCM0wiUiumh5B+GBHFsiAAF2LSURyU35eYu+lEYSIiFwgzwuISR+sZApkQGS/W0VEliY5gpjUFFPm6NPkRCQX5U1PMWX/chuBDQgRkVyUHEH4IB+CGRA+GJmJiCxJSCOIzNMEk4jkopA3PT6pVUwiIjJTcgShgMgQp3VMIpKjkutr/LAXC2RAAJpjEpGc5odjqSsWEGa2xcyeM7MvmtndK/W8IiK55PwS/ewnxKICwsyeNLMuM9s/q327mR02s3Yze/giD3M78NfOuY8Dv7vEehfFD8krIrIU0/Hgg/1Y/iLv9xTwBeDpZIOZhYBHgfcBHcBeM9sJhIBHZv38/cBXgM+Y2QeB6vTKvjjNMIlILvLTOb6LCgjn3G4za57VfBPQ7pw7CmBmzwJ3OeceAe6c56H+2AuWv11auSIiq4MPBhCLHkGk0gCcmLHdAdw83529gPmPQDHwPxa434PAgwDr169PozwRkdxj3vxHLk0xpZJqIDTvS3LOHcPb8S/EObcD2AHQ2tq65C7StZhEJBedX+aa/YRIZxVTB9A0Y7sR6EyvHBGR1c1PB6nTCYi9wCYzazGzMHAPsHN5ykqP80PPiogswfQIwge7scUuc30G2ANsNrMOM3vAORcHHgJeBA4CzznnDmSu1EujGSYRyU3eMQgfTDEtdhXTvfO0vwC8sKwViYisYn56cxvIS21kP3dFRNKTM1NMuchHISwismh+2ncFNiBERHJRcom+RhAZ4oeOFRFZivOX6sv+jiyQAQE6UU5EclPOLXMVEZGVoQ8MyjA/DM1ERJbCfHSYOpABAf5aCSAicqn8cEWIwAaEiEhO0hRTZvkgeEVEliQoF+vzNS1iEpFclHOfSS0iIiujtiTCHdespaIonO1S0vrAIN/Kfu6KiCzNlnVlPPqRG7JdBhDoEYTmmERE0hHggBARkXQEMiD8cPRfRCTXBTIgQKuYRETSFdiAEBGR9AQ0IDTHJCKSroAGhNYwiYikK7ABISIi6TE/XDFwPmbWDQwDZ7Ndi0/UoL4A9cNM6osE9cN5NUCxc6423QfydUAAmFmbc64123X4gfoiQf1wnvoiQf1w3nL2haaYREQkJQWEiIiklAsBsSPbBfiI+iJB/XCe+iJB/XDesvWF749BiIhIduTCCEJERLLA1wFhZtvN7LCZtZvZw9muZ7mZ2ZNm1mVm+2e0VZnZ983sLe97pdduZva/vL543cxumPEzv+fd/y0z+71svJZ0mFmTmf3AzA6a2QEz+zde+2rsi6iZvWJmP/P64i+89hYze9l7XV83s7DXHvG2271/b57xWJ/22g+b2a9l5xWlx8xCZrbPzJ73tldrPxwzs5+b2Wtm1ua1Zf7vwznnyy8gBBwBLgPCwM+ALdmua5lf478AbgD2z2j778DD3u2Hgc96t38d+A6Jk8R/AXjZa68CjnrfK73bldl+bZfYD2uBG7zbpcCbwJZV2hcGlHi3C4CXvdf4HHCP1/4Y8HHv9h8Bj3m37wG+7t3e4v3NRIAW728plO3Xt4T++LfA14Dnve3V2g/HgJpZbRn/+/DzCOImoN05d9Q5Nw48C9yV5ZqWlXNuN9A7q/ku4Mve7S8DvzGj/WmX8BOgwszWAr8GfN851+uc6wO+D2zPfPXLxzl3yjn3U+/2IHAQaGB19oVzzg15mwXelwNuBb7htc/ui2QffQP4VUv0sX02AAACo0lEQVR8qPFdwLPOuZhz7m2gncTfVM4ws0bgDuBL3raxCvthARn/+/BzQDQAJ2Zsd3htQVfvnDsFiR0nUOe1z9cfgeonb2rg3STeOa/KvvCmVV4Dukj8ER8Bzjnn4t5dZr6u6dfs/Xs/UE0w+uKvgH8PTHnb1azOfoDEm4TvmdmrZvag15bxvw8/fyZ1quvtreYlV/P1R2D6ycxKgG8Cn3TODdj8H+oR6L5wzk0C15tZBfB3wJWp7uZ9D2RfmNmdQJdz7lUz25ZsTnHXQPfDDFudc51mVgd838wOLXDfZesLP48gOoCmGduNQGeWallJZ7zhIN73Lq99vv4IRD+ZWQGJcPiqc+5vveZV2RdJzrlzwC4S88gVZpZ8QzfzdU2/Zu/fy0lMW+Z6X2wFPmhmx0hML99KYkSx2voBAOdcp/e9i8SbhptYgb8PPwfEXmCTt2ohTOLA084s17QSdgLJ1QW/B3xrRvvveisUfgHo94aVLwLvN7NKbxXD+722nOHNFT8BHHTOfX7GP63Gvqj1Rg6YWSFwG4ljMj8A7vbuNrsvkn10N/CSSxyR3Anc463uaQE2Aa+szKtIn3Pu0865RudcM4m//Zeccx9hlfUDgJkVm1lp8jaJ3+v9rMTfR7aPzl/kyP2vk1jRcgT402zXk4HX9wxwCpggke4PkJg3/UfgLe97lXdfAx71+uLnQOuMx7mfxMG3duCj2X5dS+iHW0gMdV8HXvO+fn2V9sW1wD6vL/YDf+61X0Zix9YO/F8g4rVHve12798vm/FYf+r10WHg9my/tjT6ZBvnVzGtun7wXvPPvK8DyX3hSvx96ExqERFJyc9TTCIikkUKCBERSUkBISIiKSkgREQkJQWEiIikpIAQEZGUFBAiIpKSAkJERFL6/1b5DeLP7ev5AAAAAElFTkSuQmCC\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, 6e-6)\n",
"plt.yscale('log')\n",
"# plt.xlim(3080, 3110)\n",
"plt.savefig('test.png')\n",
"# print(jpsi_width)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Adjust scaling of different parts"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"# total_f.update_integration_options(draws_per_dim=20000000, mc_sampler=None)\n",
"# inte = total_f.integrate(limits = (3080, 3112), norm_range=False)\n",
"# inte_fl = zfit.run(inte)\n",
"# print(inte_fl)\n",
"# print(pdg[\"jpsi_BR\"]/pdg[\"NR_BR\"], inte_fl*pdg[\"psi2s_auc\"]/pdg[\"NR_auc\"])"
]
},
{
"cell_type": "code",
"execution_count": 12,
"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",
"\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() - _))))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Tensorflow scaling"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [],
"source": [
"# def scaling_func(x):\n",
"\n",
"# funcs = resonance(x, _mass = ztf.constant(jpsi_mass), scale = ztf.constant(jpsi_scale), phase = ztf.constant(jpsi_phase), width = ztf.constant(jpsi_width)) + resonance(x, _mass = ztf.constant(psi2s_mass), scale = ztf.constant(psi2s_scale), phase = ztf.constant(psi2s_phase), width = ztf.constant(psi2s_width))\n",
"\n",
"# vec_f = vec(x, funcs)\n",
"\n",
"# axiv_nr = axiv_nonres(x)\n",
"\n",
"# tot = vec_f + axiv_nr\n",
"\n",
"# return tot\n",
"\n",
"\n",
"# def s_func(x):\n",
" \n",
"# q = ztf.constant(x)\n",
" \n",
"# return zfit.run(scaling_func(q))\n",
" \n",
"\n",
"# print(integrate.quad(s_func, x_min, x_max, limit = 50))"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [],
"source": [
"# factor_jpsi = pdg[\"NR_auc\"]*pdg[\"jpsi_BR\"]/(pdg[\"NR_BR\"]*pdg[\"jpsi_auc\"])\n",
"# factor_jpsi = pdg[\"NR_auc\"]*pdg[\"jpsi_BR\"]/(pdg[\"NR_BR\"]*inte_fl)\n",
"# print(np.sqrt(factor_jpsi)*jpsi_scale)\n",
"# print(np.sqrt(factor_jpsi))\n",
"# # print(psi2s_scale)\n",
"# factor_psi2s = pdg[\"NR_auc\"]*pdg[\"psi2s_BR\"]/(pdg[\"NR_BR\"]*pdg[\"psi2s_auc\"])\n",
"# factor_psi2s = pdg[\"NR_auc\"]*pdg[\"psi2s_BR\"]/(pdg[\"NR_BR\"]*inte_fl)\n",
"# print(np.sqrt(factor_psi2s)*psi2s_scale)\n",
"# print(np.sqrt(factor_psi2s))"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [],
"source": [
"# def _t_f(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",
"# funcs = psi2s_res(xq) + jpsi_res(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 t_f(x):\n",
"# _ = np.array(x)\n",
"# probs = zfit.run(_t_f(_))\n",
"# return probs"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [],
"source": [
"# print(36000*(1+ pdg[\"jpsi_BR\"]/pdg[\"NR_BR\"] + pdg[\"psi2s_BR\"]/pdg[\"NR_BR\"]))"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [],
"source": [
"# start = time.time()\n",
"\n",
"# result, err = integrate.quad(lambda x: t_f(x), x_min, x_max, limit = 5)\n",
"# print(result, \"{0:.2f} %\".format(err/result))\n",
"# print(\"Time:\", time.time()-start)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Sampling\n",
"## One sample\n",
"! total_f.sample() always returns the same set !"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [],
"source": [
"# nevents = int(pdg[\"number_of_decays\"])\n",
"# event_stack = 5000\n",
"\n",
"# calls = int(nevents/event_stack + 1)\n",
"\n",
"# total_samp = []\n",
"\n",
"# start = time.time()\n",
"\n",
"# samp = total_f.sample(n=event_stack)\n",
"# s = samp.unstack_x()\n",
"\n",
"# for call in range(calls):\n",
"\n",
"# sam = zfit.run(s)\n",
"# clear_output(wait=True)\n",
" \n",
"# # if call != 0:\n",
"# # print(np.sum(_last_sam-sam))\n",
" \n",
"# # _last_sam = sam\n",
" \n",
"# c = call + 1 \n",
"# print(\"{0}/{1}\".format(c, calls))\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-c)))))\n",
" \n",
"# with open(\"data/zfit_toys/toy_1/{}.pkl\".format(call), \"wb\") as f:\n",
"# pkl.dump(sam, f, pkl.HIGHEST_PROTOCOL)"
]
},
{
"cell_type": "code",
"execution_count": 19,
"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_1/{}.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",
"# print(total_samp[:nevents].shape)"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [],
"source": [
"# bins = int((x_max-x_min)/7)\n",
"\n",
"# # calcs = zfit.run(total_test_tf(samp))\n",
"\n",
"# plt.hist(total_samp[:event_stack], bins = bins, range = (x_min,x_max))\n",
"\n",
"# # plt.plot(sam, calcs, '.')\n",
"# # plt.plot(test_q, calcs_test)\n",
"# plt.ylim(0, 20)\n",
"# # plt.xlim(3000, 3750)\n",
"\n",
"# plt.savefig('test2.png')\n",
"# 1-(0.21+0.62)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Toys"
]
},
{
"cell_type": "code",
"execution_count": 21,
"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.1, dtype=dtype),\n",
" tf.constant(0.7, dtype=dtype),\n",
" tf.constant(0.2, dtype=dtype)]),\n",
" components_distribution=tfd.Uniform(low=[tf.constant(x_min, dtype=dtype), \n",
" tf.constant(2800, dtype=dtype),\n",
" tf.constant(3550, dtype=dtype)], \n",
" high=[tf.constant(x_max, dtype=dtype),\n",
" tf.constant(3300, dtype=dtype), \n",
" tf.constant(3900, 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)\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": 22,
"metadata": {},
"outputs": [],
"source": [
"# total_f._sample_and_weights = UniformSampleAndWeights"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"scrolled": false
},
"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"
]
},
{
"ename": "InvalidArgumentError",
"evalue": "assertion failed: [32254] [1.1806250218058597e-06 1.1806250218058597e-06 1.1806250218058597e-06...] [8.2086369510094e-06 2.2747180268017141e-06 1.6044321221245944e-06...]\n\t [[node create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert (defined at c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py:294) ]]\n\nCaused by op 'create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert', defined at:\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\runpy.py\", line 193, in _run_module_as_main\n \"__main__\", mod_spec)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\runpy.py\", line 85, in _run_code\n exec(code, run_globals)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel_launcher.py\", line 16, in <module>\n app.launch_new_instance()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\traitlets\\config\\application.py\", line 658, in launch_instance\n app.start()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelapp.py\", line 505, in start\n self.io_loop.start()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\platform\\asyncio.py\", line 148, in start\n self.asyncio_loop.run_forever()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\base_events.py\", line 539, in run_forever\n self._run_once()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\base_events.py\", line 1775, in _run_once\n handle._run()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\events.py\", line 88, in _run\n self._context.run(self._callback, *self._args)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\ioloop.py\", line 690, in <lambda>\n lambda f: self._run_callback(functools.partial(callback, future))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\ioloop.py\", line 743, in _run_callback\n ret = callback()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 781, in inner\n self.run()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 742, in run\n yielded = self.gen.send(value)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 357, in process_one\n yield gen.maybe_future(dispatch(*args))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 267, in dispatch_shell\n yield gen.maybe_future(handler(stream, idents, msg))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 534, in execute_request\n user_expressions, allow_stdin,\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\ipkernel.py\", line 294, in do_execute\n res = shell.run_cell(code, store_history=store_history, silent=silent)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\zmqshell.py\", line 536, in run_cell\n return super(ZMQInteractiveShell, self).run_cell(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 2848, in run_cell\n raw_cell, store_history, silent, shell_futures)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 2874, in _run_cell\n return runner(coro)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\async_helpers.py\", line 67, in _pseudo_sync_runner\n coro.send(None)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3049, in run_cell_async\n interactivity=interactivity, compiler=compiler, result=result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3214, in run_ast_nodes\n if (yield from self.run_code(code, result)):\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3296, in run_code\n exec(code_obj, self.user_global_ns, self.user_ns)\n File \"<ipython-input-23-3eb82b0d7abd>\", line 13, in <module>\n sampler = total_f.create_sampler(n=event_stack)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 814, in create_sampler\n limits=limits, n=n, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 835, in _create_sampler_tensor\n sample = self._single_hook_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 877, in _single_hook_sample\n return self._hook_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 491, in _hook_sample\n samples = super()._hook_sample(limits=limits, n=n, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 880, in _hook_sample\n return self._norm_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 884, in _norm_sample\n return self._limits_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 888, in _limits_sample\n return self._call_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 898, in _call_sample\n return self._fallback_sample(n=n, limits=limits)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 931, in _fallback_sample\n sample_and_weights_factory=self._sample_and_weights)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py\", line 346, in accept_reject_sample\n back_prop=False)[1] # backprop not needed here\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3556, in while_loop\n return_same_structure)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3087, in BuildLoop\n pred, body, original_loop_vars, loop_vars, shape_invariants)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3022, in _BuildLoop\n body_result = body(*packed_vars_for_body)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py\", line 294, in sample_body\n message=\"Not all weights are >= probs so the sampling \"\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\check_ops.py\", line 1023, in assert_greater_equal\n return control_flow_ops.Assert(condition, data, summarize=summarize)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\tf_should_use.py\", line 193, in wrapped\n return _add_should_use_warning(fn(*args, **kwargs))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 168, in Assert\n guarded_assert = cond(condition, no_op, true_assert, name=\"AssertGuard\")\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\deprecation.py\", line 507, in new_func\n return func(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 2108, in cond\n orig_res_f, res_f = context_f.BuildCondBranch(false_fn)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 1941, in BuildCondBranch\n original_result = fn()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 166, in true_assert\n condition, data, summarize, name=\"Assert\")\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_logging_ops.py\", line 72, in _assert\n name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\op_def_library.py\", line 788, in _apply_op_helper\n op_def=op_def)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\deprecation.py\", line 507, in new_func\n return func(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\ops.py\", line 3300, in create_op\n op_def=op_def)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\ops.py\", line 1801, in __init__\n self._traceback = tf_stack.extract_stack()\n\nInvalidArgumentError (see above for traceback): assertion failed: [32254] [1.1806250218058597e-06 1.1806250218058597e-06 1.1806250218058597e-06...] [8.2086369510094e-06 2.2747180268017141e-06 1.6044321221245944e-06...]\n\t [[node create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert (defined at c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py:294) ]]\n",
"output_type": "error",
"traceback": [
"\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[1;31mInvalidArgumentError\u001b[0m Traceback (most recent call last)",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_do_call\u001b[1;34m(self, fn, *args)\u001b[0m\n\u001b[0;32m 1333\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-> 1334\u001b[1;33m \u001b[1;32mreturn\u001b[0m \u001b[0mfn\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m*\u001b[0m\u001b[0margs\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 1335\u001b[0m \u001b[1;32mexcept\u001b[0m \u001b[0merrors\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mOpError\u001b[0m \u001b[1;32mas\u001b[0m \u001b[0me\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_run_fn\u001b[1;34m(feed_dict, fetch_list, target_list, options, run_metadata)\u001b[0m\n\u001b[0;32m 1318\u001b[0m return self._call_tf_sessionrun(\n\u001b[1;32m-> 1319\u001b[1;33m options, feed_dict, fetch_list, target_list, run_metadata)\n\u001b[0m\u001b[0;32m 1320\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_call_tf_sessionrun\u001b[1;34m(self, options, feed_dict, fetch_list, target_list, run_metadata)\u001b[0m\n\u001b[0;32m 1406\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_session\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0moptions\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mfeed_dict\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mfetch_list\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mtarget_list\u001b[0m\u001b[1;33m,\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m-> 1407\u001b[1;33m run_metadata)\n\u001b[0m\u001b[0;32m 1408\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;31mInvalidArgumentError\u001b[0m: assertion failed: [32254] [1.1806250218058597e-06 1.1806250218058597e-06 1.1806250218058597e-06...] [8.2086369510094e-06 2.2747180268017141e-06 1.6044321221245944e-06...]\n\t [[{{node create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert}}]]",
"\nDuring handling of the above exception, another exception occurred:\n",
"\u001b[1;31mInvalidArgumentError\u001b[0m Traceback (most recent call last)",
"\u001b[1;32m<ipython-input-23-3eb82b0d7abd>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m 23\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0mcall\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mcalls\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 24\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 25\u001b[1;33m \u001b[0msampler\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mresample\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mn\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mevent_stack\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 26\u001b[0m \u001b[0ms\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0msampler\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0munstack_x\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 27\u001b[0m \u001b[0msam\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mzfit\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mrun\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0ms\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\data.py\u001b[0m in \u001b[0;36mresample\u001b[1;34m(self, param_values, n)\u001b[0m\n\u001b[0;32m 637\u001b[0m \u001b[1;32mraise\u001b[0m \u001b[0mRuntimeError\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Cannot set a new `n` if not a Tensor-like object was given\"\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 638\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mn_samples\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mload\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mvalue\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mn\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0msession\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;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 639\u001b[1;33m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0msess\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mrun\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0msample_holder\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0minitializer\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 640\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_initial_resampled\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;32mTrue\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 641\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36mrun\u001b[1;34m(self, fetches, feed_dict, options, run_metadata)\u001b[0m\n\u001b[0;32m 927\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[0;32m 928\u001b[0m result = self._run(None, fetches, feed_dict, options_ptr,\n\u001b[1;32m--> 929\u001b[1;33m run_metadata_ptr)\n\u001b[0m\u001b[0;32m 930\u001b[0m \u001b[1;32mif\u001b[0m \u001b[0mrun_metadata\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 931\u001b[0m \u001b[0mproto_data\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mtf_session\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mTF_GetBuffer\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mrun_metadata_ptr\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_run\u001b[1;34m(self, handle, fetches, feed_dict, options, run_metadata)\u001b[0m\n\u001b[0;32m 1150\u001b[0m \u001b[1;32mif\u001b[0m \u001b[0mfinal_fetches\u001b[0m \u001b[1;32mor\u001b[0m \u001b[0mfinal_targets\u001b[0m \u001b[1;32mor\u001b[0m \u001b[1;33m(\u001b[0m\u001b[0mhandle\u001b[0m \u001b[1;32mand\u001b[0m \u001b[0mfeed_dict_tensor\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 1151\u001b[0m results = self._do_run(handle, final_targets, final_fetches,\n\u001b[1;32m-> 1152\u001b[1;33m feed_dict_tensor, options, run_metadata)\n\u001b[0m\u001b[0;32m 1153\u001b[0m \u001b[1;32melse\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 1154\u001b[0m \u001b[0mresults\u001b[0m \u001b[1;33m=\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;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_do_run\u001b[1;34m(self, handle, target_list, fetch_list, feed_dict, options, run_metadata)\u001b[0m\n\u001b[0;32m 1326\u001b[0m \u001b[1;32mif\u001b[0m \u001b[0mhandle\u001b[0m \u001b[1;32mis\u001b[0m \u001b[1;32mNone\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 1327\u001b[0m return self._do_call(_run_fn, feeds, fetches, targets, options,\n\u001b[1;32m-> 1328\u001b[1;33m run_metadata)\n\u001b[0m\u001b[0;32m 1329\u001b[0m \u001b[1;32melse\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 1330\u001b[0m \u001b[1;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_do_call\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0m_prun_fn\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mhandle\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mfeeds\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mfetches\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
"\u001b[1;32mc:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\client\\session.py\u001b[0m in \u001b[0;36m_do_call\u001b[1;34m(self, fn, *args)\u001b[0m\n\u001b[0;32m 1346\u001b[0m \u001b[1;32mpass\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 1347\u001b[0m \u001b[0mmessage\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0merror_interpolation\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0minterpolate\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mmessage\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_graph\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m-> 1348\u001b[1;33m \u001b[1;32mraise\u001b[0m \u001b[0mtype\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0me\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mnode_def\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mop\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mmessage\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 1349\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 1350\u001b[0m \u001b[1;32mdef\u001b[0m \u001b[0m_extend_graph\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\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;31mInvalidArgumentError\u001b[0m: assertion failed: [32254] [1.1806250218058597e-06 1.1806250218058597e-06 1.1806250218058597e-06...] [8.2086369510094e-06 2.2747180268017141e-06 1.6044321221245944e-06...]\n\t [[node create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert (defined at c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py:294) ]]\n\nCaused by op 'create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert', defined at:\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\runpy.py\", line 193, in _run_module_as_main\n \"__main__\", mod_spec)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\runpy.py\", line 85, in _run_code\n exec(code, run_globals)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel_launcher.py\", line 16, in <module>\n app.launch_new_instance()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\traitlets\\config\\application.py\", line 658, in launch_instance\n app.start()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelapp.py\", line 505, in start\n self.io_loop.start()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\platform\\asyncio.py\", line 148, in start\n self.asyncio_loop.run_forever()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\base_events.py\", line 539, in run_forever\n self._run_once()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\base_events.py\", line 1775, in _run_once\n handle._run()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\asyncio\\events.py\", line 88, in _run\n self._context.run(self._callback, *self._args)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\ioloop.py\", line 690, in <lambda>\n lambda f: self._run_callback(functools.partial(callback, future))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\ioloop.py\", line 743, in _run_callback\n ret = callback()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 781, in inner\n self.run()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 742, in run\n yielded = self.gen.send(value)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 357, in process_one\n yield gen.maybe_future(dispatch(*args))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 267, in dispatch_shell\n yield gen.maybe_future(handler(stream, idents, msg))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\kernelbase.py\", line 534, in execute_request\n user_expressions, allow_stdin,\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tornado\\gen.py\", line 209, in wrapper\n yielded = next(result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\ipkernel.py\", line 294, in do_execute\n res = shell.run_cell(code, store_history=store_history, silent=silent)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\ipykernel\\zmqshell.py\", line 536, in run_cell\n return super(ZMQInteractiveShell, self).run_cell(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 2848, in run_cell\n raw_cell, store_history, silent, shell_futures)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 2874, in _run_cell\n return runner(coro)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\async_helpers.py\", line 67, in _pseudo_sync_runner\n coro.send(None)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3049, in run_cell_async\n interactivity=interactivity, compiler=compiler, result=result)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3214, in run_ast_nodes\n if (yield from self.run_code(code, result)):\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\IPython\\core\\interactiveshell.py\", line 3296, in run_code\n exec(code_obj, self.user_global_ns, self.user_ns)\n File \"<ipython-input-23-3eb82b0d7abd>\", line 13, in <module>\n sampler = total_f.create_sampler(n=event_stack)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 814, in create_sampler\n limits=limits, n=n, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 835, in _create_sampler_tensor\n sample = self._single_hook_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 877, in _single_hook_sample\n return self._hook_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 491, in _hook_sample\n samples = super()._hook_sample(limits=limits, n=n, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 880, in _hook_sample\n return self._norm_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 884, in _norm_sample\n return self._limits_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 888, in _limits_sample\n return self._call_sample(n=n, limits=limits, name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 898, in _call_sample\n return self._fallback_sample(n=n, limits=limits)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basemodel.py\", line 931, in _fallback_sample\n sample_and_weights_factory=self._sample_and_weights)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py\", line 346, in accept_reject_sample\n back_prop=False)[1] # backprop not needed here\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3556, in while_loop\n return_same_structure)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3087, in BuildLoop\n pred, body, original_loop_vars, loop_vars, shape_invariants)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 3022, in _BuildLoop\n body_result = body(*packed_vars_for_body)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py\", line 294, in sample_body\n message=\"Not all weights are >= probs so the sampling \"\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\check_ops.py\", line 1023, in assert_greater_equal\n return control_flow_ops.Assert(condition, data, summarize=summarize)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\tf_should_use.py\", line 193, in wrapped\n return _add_should_use_warning(fn(*args, **kwargs))\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 168, in Assert\n guarded_assert = cond(condition, no_op, true_assert, name=\"AssertGuard\")\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\deprecation.py\", line 507, in new_func\n return func(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 2108, in cond\n orig_res_f, res_f = context_f.BuildCondBranch(false_fn)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 1941, in BuildCondBranch\n original_result = fn()\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\control_flow_ops.py\", line 166, in true_assert\n condition, data, summarize, name=\"Assert\")\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_logging_ops.py\", line 72, in _assert\n name=name)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\op_def_library.py\", line 788, in _apply_op_helper\n op_def=op_def)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\util\\deprecation.py\", line 507, in new_func\n return func(*args, **kwargs)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\ops.py\", line 3300, in create_op\n op_def=op_def)\n File \"c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\framework\\ops.py\", line 1801, in __init__\n self._traceback = tf_stack.extract_stack()\n\nInvalidArgumentError (see above for traceback): assertion failed: [32254] [1.1806250218058597e-06 1.1806250218058597e-06 1.1806250218058597e-06...] [8.2086369510094e-06 2.2747180268017141e-06 1.6044321221245944e-06...]\n\t [[node create_sampler/while/assert_greater_equal/Assert/AssertGuard/Assert (defined at c:\\users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\sample.py:294) ]]\n"
]
}
],
"source": [
"# zfit.run.numeric_checks = False \n",
"\n",
"nr_of_toys = 1\n",
"nevents = int(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",
" print(\"{0}/{1}\".format(c, calls))\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-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": null,
"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": null,
"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": null,
"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, 40000)\n",
"# plt.xlim(3080, 3110)\n",
"\n",
"plt.legend()\n",
"\n",
"plt.savefig('test2.png')"
]
},
{
"cell_type": "code",
"execution_count": null,
"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": null,
"metadata": {},
"outputs": [],
"source": [
"# plt.hist(sample, weights=1 / prob(sample))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Fitting"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"nll = zfit.loss.UnbinnedNLL(model=total_f, data=data2, fit_range = (x_min, x_max))\n",
"\n",
"minimizer = zfit.minimize.MinuitMinimizer()\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)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# 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": null,
"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, 5e-4)\n",
"# plt.yscale('log')\n",
"# plt.xlim(3080, 3110)\n",
"plt.savefig('test3.png')\n",
"# print(jpsi_width)"
]
},
{
"cell_type": "code",
"execution_count": null,
"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": null,
"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": null,
"metadata": {},
"outputs": [],
"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
}
