diff --git a/.ipynb_checkpoints/raremodel-nb-checkpoint.ipynb b/.ipynb_checkpoints/raremodel-nb-checkpoint.ipynb index 47c3623..529e151 100644 --- a/.ipynb_checkpoints/raremodel-nb-checkpoint.ipynb +++ b/.ipynb_checkpoints/raremodel-nb-checkpoint.ipynb @@ -9,31 +9,9 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "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" - ] - } - ], + "outputs": [], "source": [ "import os\n", "\n", @@ -64,7 +42,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -90,7 +68,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -280,7 +258,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -335,7 +313,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -354,7 +332,7 @@ " '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", - " 'DDstar_mass', 'DDstar_scale', 'DDstar_phase',\n", + " 'Dstar_mass', 'DDstar_scale', 'DDstar_phase', 'D_mass',\n", " 'tau_mass', 'C_tt']\n", "# the name of the parameters\n", "\n", @@ -404,7 +382,7 @@ " \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['DDstar_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", @@ -431,7 +409,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -459,19 +437,9 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "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" - ] - } - ], + "outputs": [], "source": [ "# formfactors\n", "\n", @@ -580,7 +548,7 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -588,16 +556,17 @@ "\n", "Dbar_phase = 0.0\n", "DDstar_phase = 0.0\n", - "DDstar_mass = pdg['D0_M']\n", - " \n", - "Dbar_mass = (pdg['D0_M']+pdg['Dst_M'])/2\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=-1.464, upper_limit=1.464)\n", "Dbar_m = zfit.Parameter(\"Dbar_m\", ztf.constant(Dbar_mass), floating = False)\n", - "Dbar_p = zfit.Parameter(\"Dbar_p\", ztf.constant(Dbar_phase), 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=-2.0, upper_limit=2.0)#, floating = False)\n", - "DDstar_m = zfit.Parameter(\"DDstar_m\", ztf.constant(DDstar_mass), floating = False)\n", - "DDstar_p = zfit.Parameter(\"DDstar_p\", ztf.constant(DDstar_phase), floating = False)" + "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)" ] }, { @@ -609,7 +578,7 @@ }, { "cell_type": "code", - "execution_count": 10, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -626,7 +595,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -639,7 +608,7 @@ " 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", - " DDstar_mass = DDstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p,\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", @@ -663,7 +632,7 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -709,22 +678,9 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ8AAAD8CAYAAACo9anUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO2dd3yc1ZX3v2equmTLsi3cZcu4gCHGoYTezZLE7IZsTEIglZdd2JJsgby7aSxkX/ImIVsghAB5CZtgCCTBCSSEBEwgAYPBNuCGZbnJTcXq0khT7vvH88xoLM+MZmRNP9/Pxx/P3Oe255E0vznnnnuuGGNQFEVRlEziyPYEFEVRlOJDxUdRFEXJOCo+iqIoSsZR8VEURVEyjoqPoiiKknFUfBRFUZSMk5T4iMhKEdkhIk0icnuM614Redy+vl5E5kZd+5JdvkNErhyrTxGZZ/ex0+7Tk8QYy0TkVRHZIiLviEjJeB6GoiiKkhnGFB8RcQL3AlcBS4DrRGTJqGqfBTqNMQuAe4C77bZLgNXAUmAlcJ+IOMfo827gHmNMI9Bp951oDBfwP8DNxpilwEWAP8XnoCiKomSQZCyfM4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FIREbt8jTFmyBizG2iy+4vZp93mErsP7D6vGWOMK4C3jTGbAYwxHcaYYPKPQFEURck0riTqzAD2R71vAc6KV8cYExCRbqDWLn9tVNsZ9utYfdYCXcaYQIz68cZYCBgReQ6owxK7b46+CRG5CbgJoLy8/IxFixYlcevFwY4jvQwHQiycVonXpcuAisXOI314XA7m1JbFrdPj87O3Y4AFdRWUepwZnN0IzW39ADTUlR9Tvv1wL/5giEllHjwuB0d6fJwyoxpJ41yO9Pho7R0CoK7Sy/SqwloBePPNN9uNMXUT0Vcy4hPrZzU6J0+8OvHKY33CJaqfaAwXcB7wfmAA+L2IvGmM+f0xFY15AHgAYMWKFWbDhg0xuitOLv7WOna39/PzL17AgqmV2Z6OkiNc9p2XWDitgvs+cUbcOut2tPKpH77Bj/7qA5wxZ1IGZzfCtd/7E163gx9/7uxjyi/59jqa2/r5+FmzmT25jP/z6+28csfKtIrkXc9s5cfr9+FxOVh12kl8fdUpaRsrG4jI3onqK5mvuS3ArKj3M4GD8erYazDVwNEEbeOVtwM1dh+jx0o0xkvGmHZjzADwLLA8iftSRqFp/pRogiGDy5H4I8LttK4HgqFMTCkmIWNwyPHfTb0uS2Q8TkdknsNpnufAcJBSt5MKr4veocDYDYqYZMTnDaDRjkLzYAUQrB1VZy1wo/36WuAFY2UsXQustiPV5gGNwOvx+rTbvGj3gd3n02OM8RywTETKbFG6ENia/CNQFCUW/mAIlyOxkyoiPqHsfXMJGpAY4uOxXcgelwOP07ruT7P4DPqDlNji0+dT8UnEmG43e33lVqwPeSfwsDFmi4jcAWwwxqwFHgIeFZEmLGtktd12i4g8gSUGAeCWcDBArD7tIW8D1ojIncBGu28SjNEpIt/BEjQDPGuMeeaEnoqiKASCBpczsfiEr6fbokiEMYZYGllii4/bKRGRTLf4+PxByjy2+Kjlk5Bk1nwwxjyL5c6KLvtK1Gsf8NE4be8C7kqmT7u8GSsabnR5ojH+Byvcetz4/X5aWlrw+Xwn0k3OUFJSwsyZM3G73Um3Ua+bEk0gZHA5x3C7OcJut+z99oSMwRnD8qnwWh9vHqdzRHwC6Z3nwHCQUo+Tcq+LroHhtI6V7yQlPsVAS0sLlZWVzJ07N6YJn08YY+jo6KClpYV58+ZlezpKnhIIJeF2c1nXs7rmE4rtdgtbZeVeJ25XZtZ8BodH3G4tnQNpHSvf0bhaG5/PR21tbd4LD1h/iLW1tQVjxSnZIRgcO+AgfD2bbrdQHLdbeG4VXlfG1nzCbrdyr5P+Id1umAgVnygKQXjCjOdeNNpNicYfCo255uN2hi2f7P3yGEPMaDeHrUhetyNjaz7haLdyXfMZExUfRVFi4g+aiLjEYyTaLcuWT4xPstNmVgMwraokY+Iz6B8Jte4fDmD0G11cVHzyjHXr1vHBD34QgKGhIS677DJOP/10Hn/88SzPTCkkgiFDMGTwOBNvyByJdstuwEEsS//GD8zlR585k3Maakf2+aQ54MDnD1JiBxwYY1lCSmw04CCP2bhxI36/n02bNk1If0bj3RSbsIUQDiiIx0i0WzYtn9huN7fTwQULrUwwHldm1nwGhoOU2W43gP6hQOS1cixq+eQQe/bsYdGiRdx4440sW7aMa6+9loGBAX7zm9+waNEizjvvPH72s58B0NrayvXXX8+mTZs4/fTT2bVrV5ZnrxQS4Q9pz1ih1q5cCbVOXGfE8kmf+BhjLLebx0mF17IYdd0nPirJMfj6L7ew9WDPhPa55KQqvvqhpWPW27FjBw899BDnnnsun/nMZ/jOd77D97//fV544QUWLFjAxz72MQCmTp3Kgw8+yLe+9S1+9atfTehcFcVvi4l7DPEJh2JnP9otsfqUuC0xGEqj+AwFQhhjjVXuCVs+6naLh1o+OcasWbM499xzAbj++uvZsGED8+bNo7GxERHh+uuvT9vYujaqhIm43cayfJw5YPnE2ecTTaktPoP+9InBoL2+E85wAGr5JEItnxgkY6Gki9F/RN3d3QUVAq7kB2H31FjRbk6H4JDsRrvFS68TjddtiaQvneJj910yas1HiY1aPjnGvn37ePXVVwF47LHHuOyyy9i9e3dkTeexxx7L5vSUIiHsRvMkcb6Ty+nIststdsBBNGHLJ53iE+67zBMlPsMqPvFQ8ckxFi9ezCOPPMKyZcs4evQoX/jCF3jggQe4+uqrOe+885gzZ07axla3mxIm2YADALdDsh5wMEYihsiaz2AaQ5/Dlo/X5aSyRN1uY6FutxzD4XBw//33H1O2cuVKtm/fflzdiy66iIsuuihDM1OKiXACzrHWfMCKeMv2eT5juabdTgcuh+ALpN/yKfWo2y0Z1PJRFOU4hiP7fJJwuzmy73aLldV6NCVuJ4PD6Zunz2/1XeJyUOYOh1prtFs8VHxyiLlz5/Luu+9mbXzdZKqEGYl2G/tD3etypDWEeSziJRYdTYnbmVbLJ+zSK/U4cTiEco9TLZ8EqPhEUUh5mArpXpTMk8qaj9ftSOvmzbEIhcZ2uwGUuB340rjmExa28PpSudel4pMAFR+bkpISOjo6CuJDO3yeT0lJSbanouQpI6HWSYiPy5lVyydeVuvRZMzyscVHTzNNjAYc2MycOZOWlhba2tqyPZUJIXySaSoUgO4qE0Sym0zBCsfOB7dbqduZ1mg3n/0MwnuK1PJJjIqPjdvt1lM/FcUmnKXaM0ZiUbDXfNK4f2YsQmbk7J5ElLgdkaCAdOAbZfnogXKJUbeboijH4Q+E13wSH6kA2Q84CBpDMklAStzOtKbX8fmPXfNRt1tiVHwURTmOZI9UAGvNJ5sBByaJxKJgr/mkOb2OyyERV2W5faCcEhsVH0VRjiOVNR+v28FQGhfyxyLZfT6laRYfnz8UcbmBrvmMhYqPoijHMZzkkQoAXmd+BByUeZz0pzm9jjdKfNTtlhgVHyWCRrspYYYDqe3zyZb4GGMwZuwjFQAqS1z0+dInBkP+IKWekedV7nHh84eymnool1HxURTlOFLLcJC9NZ/wF6Zk1nwqvG4G/cG0icGgP0iJK9rtZr3WiLfYqPgoinIc/mAIEeu8nrGwot2y8wEbstUnGbdbujNN++wjtMNEDpTToIOYqPgoETS3mxJmOBjC43Qk5c4KbzLNRnaQUNjySUJ9Kmzx6U2T6+14y0czWyciKfERkZUiskNEmkTk9hjXvSLyuH19vYjMjbr2Jbt8h4hcOVafIjLP7mOn3acn0RgiMldEBkVkk/3v2PMIFEVJmeFAKKlgA7AsH2PAn4UzfcKWTzL7fCq96RUfnz9ESSzLR8UnJmP+domIE7gXuApYAlwnIktGVfss0GmMWQDcA9xtt10CrAaWAiuB+0TEOUafdwP3GGMagU6777hj2Owyxpxu/7s5pSegKMpxDAVClLiTFR+n3Sbzrrew+CQTal1Z4gbS7HaLemZq+SQmmd+uM4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FKx7PVVwBpjzJAxZjfQZPcXs0+7zSV2H9h9XjPGGMoEodFuShifPxgRlbEI5zLLRtBBKJWAg8iajz8tc/H5g5HsBhAdcKDiE4tkxGcGsD/qfYtdFrOOMSYAdAO1CdrGK68Fuuw+Ro8VbwyAeSKyUUReEpHzk7gnRVESMBQIRURlLMLh2NkIt07J7ZbmNZ/+4SBlMd1uGu0Wi2QSi8b6sY7+jhyvTrzyWL/VieonGuMQMNsY0yEiZwC/EJGlxpieYyYochNwE8Ds2bNjdKUoSpihUYvniQiLVDbEx9hDJmP5pHvNp38oQLln5CNV3W6JSearTQswK+r9TOBgvDoi4gKqgaMJ2sYrbwdq7D5GjxVzDNul1wFgjHkT2AUsHH0TxpgHjDErjDEr6urqkrjt4kO9bkqYVCyfXFjzSSbUuiKNodahkGFgOBgRHNCAg7FI5rfrDaDRjkLzYAUQrB1VZy1wo/36WuAFY8VdrgVW25Fq84BG4PV4fdptXrT7wO7z6URjiEidHcCAiDTYYzQn/wgUXThTRuNLxfJxZXPNxxafJNSn1O3E7RS6Byd+zWfAzhlXESU+XpcDp0PU8onDmG43Y0xARG4FngOcwMPGmC0icgewwRizFngIeFREmrAsntV22y0i8gSwFQgAtxhjggCx+rSHvA1YIyJ3Ahvtvok3BnABcIeIBIAgcLMx5uj4H0nxoRaPMpqhQIjy8uSO+xqxfDIvPsHIms/Y4iMiTCrzcLRveMLnERaYMu+IYIuI5ndLQFK/XcaYZ4FnR5V9Jeq1D/honLZ3AXcl06dd3owVDTe6POYYxpingKfGvAllTArhCHFlYkjJ8rHdc+nMGB2PkfQ6ydWfXO6ho3/ixScsMNGWD0BNmZuugfRE1+U7muFAUbebchyprPmEjxEYSGPG6HgEQ8nv8wFLfI72D034PAbsiLbogAOAmjIPnQMTL3aFgIqPom435ThSsXzC4cWD2RSfJE0fS3zSZ/lEu90AJpW5VXzioOKjRFARUsKkYvmU2d/2s2n5uJLIvg1Qmybx6Y/jdptc5qGzX91usVDxUdTtphzH6N36iQhnch7IQvbmQCgcap2s5eOlxxeIHBkxUYSPyy4/bs3HQ5daPjFR8VHU4lGOwRhjWT6uZC2f7LndwqHWLkdyc51cbuV365xg6ydewMGkMjf9w8GsHjOeq6j4KBE02E0B6zgFY0ja8nE7HbidEtnrkkkCwfCaT3L1ayu8ALSPI9y6e9DP81uPEAod/4cSPiH1uGi3cg+ARrzFQMVHUbebcgw+v+WSStbyASvibSAL+1lGAg6Sm+u0Kkt8jvT4Uh7rrme28vkfbeCpt1qOu9Y54MfjdByT2w2sNR/rurreRqPio6jbTTmGsIvIm6TlA9ZaR1YCDiJut+S+QtVXlwJwsHsw5bFebe4A4NfvHj7uWvfgMNVl7uM2u04qC7v51PIZjYqPEoXKkAJDtuVTkorl43Fmxe0WDFlzTSa9DsDUSi8OgUNdqVk+3YN+9h+1BGvjvs7jNmR3DfipKXUf165GLZ+4qPgo6nZTjmE8lk+Zx5mlfT7W/8laPi6ng2lVJSlbPvuPDgBwfuMUOgf8tHQe275rwE9N2fHiM9le80lHeHe+o+KjqL2jHINvHJZPmduVpVBra67JbjIFqK8uSdnyae216l+xZBoAm1u6jrneNeiPWDnR1FZYZW29E59VId9R8VEiaLSbAiMbJkfvWUlEadYsn9QyHADU15RyKEXLp7XHEo/zGuvwOB28c6D7mOtdA8Mx3W5up4MpFR5aVXyOQ8VHUbebcgzhwIHRkVuJKPM4s5vbLQXxmVlTysEuX6RtMoTF46SaEk6eXsm7UeJjjKGjfzjiYhtNXWUJreOIrit0VHwUdbspxzAiPqlZPvmQWBRg7pRyhoMhDnYlb/209vqoKXPjdTk5ZUY17x7oiQQddA74GQ6EmFZVErPttCovR3pVfEaj4qNEUBFSYCRVTCqWT7nHxWA2NpmOw/KZX1cBwK62vqTbtPYMMbXS2iN06ozqY6LfDndbwjK9Oo74VJZE3HbKCCo+irrdlGMYGMeaT5nHmZUTO0MpJhYFaKgrB6C5rT/pNq29Q0yttMTl1BnVAJF1n/CG1XiWz9QqL+19QwQmOJ9cvqPio6jFoxxDeL9OKpZPZYmLoUAo40dpB8bhdqst91BV4qK5PXnLp613xPJZOL0Ct1Mi4hMO245n+UytKiFkSMshdvmMio8SQaPdFLAORnNIaul1quxIr15fZnfyhxOLpuJ2ExEa6irY1Zqc5WOMoa13iDo7NY/X5WTR9KpI0EFzWz8lbgf18dZ8Ksef0qeQUfFR1O2mHEP/cIByj+u4VDGJqCyxXHQ9vsy63sKJRZPNah1mcX0lWw/1JHV0fNeAn+FgKOJ2Azh1ZjWbW7oIBEM0tfYxv64ibpaFsEV0qFvFJxoVH0XdbsoxDA4HjzuRcyyqSrJj+YSj3VLUHk6xgwZGZyqIRTjMOux2A7igsY5eX4D1u4+y5WA3J0+vjNt+1qQyYCRLgmKh4qNESOZboFL49A8HKU8hzBqg0hafnsHMWj7BFM/zCbNsRg0Ab7d0j1FzJLtBXZT4XLiwjhK3g3/71Vba+4Y5e15t3PY1ZW4qvS4Vn1Go+CjqdlOOYWAoEDmdNFmqSsNut8xaPoFxWj7hoIG3D3SNWTccJh1t+ZR6nHzirDlsP9xLqdvJFUunxW0vIsycXMb+JKysYiK1rzdKQaL2jhLNwAlYPhl3u9nhy6laPl6Xk8X1VWzcN7b4tPXZ4jMqoOCfrjyZukov75tVEzOvWzSzJ5eyK4nQ7u5BP1sP9nDO/PiWVKGglo8SQUVIARgYDoxjzce2fDLudrP+TyXaLczZDbVs2tc1Zk661p4hyjzO404pLXE7ufnC+ZzVMLZQzJ5cxv6jA2O6tr/+yy1c94PXeO9I79g3kOeo+CjqdlOOoX84mNIeH7AyHDgk82634DiyWoc5Z34tw8EQb+7tTFivrW/oGJfbeJg1uYyhQGjM7Na/3XIEgE37x7bI8h0VH0UtHuUY+nyB477lj4XDIVR4XfRmONQ6nDQglU2mYc6cOxmXQ/jTrvaE9Vp7fMeEWY+H2ZOtiLfd7Yldb+G1tlSyL+QrKj5KBA12U8CyXsKh06lQVeqmZzDDAQe2+rhTSK8TptzrYvnsSby4oy1hvba+oWMi3cZD4zQrFLspQT65UMjQZZ94WgwbUlV8FHW7KRECwRADw8FIxoJUqCpx051h8fEHQ4iMz+0GcMXSaWw71MOeOBaJMYYj3T6mVp2Y+JxUXUK5x8nOI/HFp3vQj99exFLxsRGRlSKyQ0SaROT2GNe9IvK4fX29iMyNuvYlu3yHiFw5Vp8iMs/uY6fdp2esMezrs0WkT0T+MdWHUOyowaOECbvNwhkLUmFSuZvOgczmLxsOGtxOR0rZGKJZecp0AH797uGY1zsH/PQPByMbRceLiLBgWmXCQIL2vpH1IBUfQEScwL3AVcAS4DoRWTKq2meBTmPMAuAe4G677RJgNbAUWAncJyLOMfq8G7jHGNMIdNp9xx0jinuAXyd748rxGJWhoiccMDAet9vkci9HM5w80x8M4XGO34Ezc1IZp82s5pebD8aMRNtnbwwNr9mcCAunVrCzNb7lEw5GWDS9kiNFcARDMj+1M4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FKxvoqsAtYYY4aMMbuBJru/mH3abS6x+8Du85oxxkBErgGagS3J37qiKKMJh0qPx+1WW+6hoy8L4pNCAtRYXLtiFlsP9bA5RraDsPjMmgDxaZxWQVvvUGRdZzTh/URLT6qmbyjAwHDmj6jIJMn81GYA+6Pet9hlMesYYwJAN1CboG288lqgy+5j9FgxxxCRcuA24OuJbkJEbhKRDSKyoa0t8QJjsaFpdZQwvRHLJ3W3W225h96hAEOBzB0q5w+GxhVsEM01p59EmcfJj1/be9y1fR3WWtCsyaUnNAbAQjvoYNuh2K63sOWzuN6q195b2EcwJCM+sX6yoz+t4tWZqPJEY3wdy02X8HAOY8wDxpgVxpgVdXV1iaoWL6pBRU/Y7VY5DrdbbYW1KJ9J19tQIIT7BNxuYN3rNe+bwdObDx631rL1UA9zastSOlI8HiOH0MXew9PWN4TH6WD+1Ar7fWGv+yTzU2sBZkW9nwkcjFdHRFxANXA0Qdt45e1Ajd3H6LHijXEW8E0R2QP8PfC/ReTWJO5LUZRRjLjdUv+wnVxupZjJpOvNHzQntOYT5uYL5hMKGe57semY8q0He1hSX3XC/YMlzjNqSmO698CydKZUeCIbWsfakJrvJPNTewNotKPQPFgBBGtH1VkL3Gi/vhZ4wVi+nLXAajtSbR7QCLwer0+7zYt2H9h9Pp1oDGPM+caYucaYucB3gW8YY/47hWdQ9KjBo4SJBByMY81nSoUtPhm0fPwTYPkAzK4t46MrZvLY6/vZZe/FOdztY0/HAKfNqjnh/sOcPquGt1viWz5TKr2RPUVFLz72+sqtwHPANuAJY8wWEblDRD5sV3sIa/2lCfgicLvddgvwBLAV+A1wizEmGK9Pu6/bgC/afdXafccdQ5k4VISUnkE/IlAxDjdT2PI52p+5D01/MITbNTE71b5w+UJKPU7+6aebGQ6EeH6rFX594cKJc9Mvm1nN/qODMV2Tbb1D1FV4qS334pDCF5+kfsOMMc8Cz44q+0rUax/w0Tht7wLuSqZPu7wZKxpudHncMaLqfC3RdUVREtNjp9aJdypnIsJrPpl0uw0HJ8byAZhaWcKd15zC3zy2kU/98HV2tfWxbGY1ixIcFJcqy2ZaVtTmli4uPnnqMdfa+4ZYNqMap0OYXO6lLcORg5lGMxwomlZHidDRP0xteeLjAeJRVeLC63JETv7MBP4JFB+AD512Ev92zSls2t9FMGS465pTx72BNRanzarG5RBe3330mPJgyNARlcanrtKrlo9SPKgIKUf7hyLus1QREeqrSzjQlblD0/xBQ6k7tQzcY/HJs+fwsRWzcAi4JlDYAMo8Lk6bVcOruzqOKe/oHyJkiKTxmVLhiez7KVTU8lEUJUJH3/C4xQegvrqUQxkVnxPf5xMLj8sx4cIT5pyGWt450E3f0Mgm0tGnpdZVemkvcMtHxUfRtDpKhKP9JyY+J9WUcqg7c/tThico2i2TnDO/lmDI8EaU6y1s5dTZRzeE3W6FvAE8v35qSlpRESpujDF0DgwzuXz8GZxPqinhSI8vctRBurGi3fLrY+yMOZPwOB280jRyjlDbaMunwstwMJTxk2EzSX791BRFSRs9vgD+oBl3wAFYbreQIWNBB8MnmFg0G5S4nZwzv5bntx6JWDaHun2IcEzAARR2loP8+qkpaaGALXslBcJ7T05ozafGchsdzNC6j88fomSCAw4ywZVLp7Pv6ADbD1t53vZ09HNSdWnkXsLik8nIwUyj4qNEUBEqbsKbQydXjF98Zk2yEnDu7xyYkDmNhc8fpMSdfx9jly+Zhgg8t8XayLq7vZ85tSOZs4shxU7+/dQURUkL4c2hJ+J2mz25HKdDaG6LfTLoROPzByc81DoT1FV6OWveZJ58s4XhQIidR3pZYCcUBairsCxIFR+loFGLRwE40hte9C4Zdx8el4NZk0oj+dHSSSAYSss+n0xxwzlzaekc5Lu/e4/+4SBnzJkUuVZV6sLjdNBewFkOVHyUCKpBxc2Rbh9Oh0TWG8ZLQ11FRiwfX8CKqCv15Kf4XLFkGg115dy3bhcel4MLGkdyyImItdFULR9FUQqdwz0+6iq8OMeR1y2ahinl7G7vJxRK79eZwWHr0Dpvnlo+LqeD+z6xnItPruMbf34qk0a5O+sqvQWd5UDT6yiKAlhHCEyvHr/LLcz8qRUMBUK0dA4yu/bEj5+Oh89viU++ut0AFk2v4oefPi6PMmCJz4EuDbVWioBC3k2tjM3hHh/Tq05cfMKHr205GPvQtIlisADEJxGFnlxUxUdRFMBa85kIy2dRfSVup/D2gfSKT9jyycdQ62Soq/BytH+IYJrdl9miMH9qSkqoxaP0DQXoHQpMiPh4XU5Onl7JO3GOi54owms+hWz5hIyV8boQUfFRIqgEFS/hjAT1EyA+AKfOsI6LTucXm4Gw+ORptNtYhJOMHulW8VEUpUDZ22FlJJhTWz4h/S2fXUOPLxBJH5MOenx+ACpL3GkbI5vMtLNFtGQoW0SmUfFR1OJR2Nth7cuZO0HRaec1TgHgj1GZmyeaXp+V8bmqpDCDdsORgvuOqvgohY6qUNGyt2OAqhIXNWXjT60TTX11KfPryo85NmCiKXTLp6rETXWpO2N58jKNio+iKOw9OjBhLrcw5zfW8VpzBwPD6TmTptcXwO2Ugo12A5g9uYx9RzN3MmwmKdyfmpI0Guym7O04NqvyRHDVKdPx+UM8v/XIhPYbptfnp7LEjcjEH6OdK8yaXEqLut2UQkdPMi1OhgJBWjoHmTvBls/7506mvrqEtZsOTmi/YXoGA1QW6HpPmFmTy2jpHCzIvT4qPopS5Oxq7ScYMpw8vXJC+3U4hFWnz2Dde21pidjqHvRTVaDrPWFmTSpjOBjiSE/hpdlR8VHU4ilydhzpAWDRBIsPwA3nzEGAH/5xz4T33d43dMIZuHOd2ZMLN+JNxUeJoGs/xcn2w714nA7mTplYtxvASTWlfHBZPT9Zv4/D3RP77b2td4i6isIWn3n2zyQT5yNlGhUfRUWnyNlxuJf5UytwO9PzcfAPV5xM0Bj+/dfbJqzPYMjQ0T9c8JbPjJpSSt1OdrVm5mTYTJLUb5uIrBSRHSLSJCK3x7juFZHH7evrRWRu1LUv2eU7ROTKsfoUkXl2HzvtPj2JxhCRM0Vkk/1vs4j8+XgfhqIUG8YYth7sYThRccsAACAASURBVHEaXG5hZk0u4+YL5/P0poP8cvPEBB90DgwTDJmCFx+HQ2ioKy9Oy0dEnMC9wFXAEuA6EVkyqtpngU5jzALgHuBuu+0SYDWwFFgJ3CcizjH6vBu4xxjTCHTafccdA3gXWGGMOd0e4/siUtghMGlCLaDi42C3j9beIU6fXZPWcf7mkgWcMWcStz/1Npv3d51wfwc6rb0v0ybgCIhcZ35dBU2tRSg+wJlAkzGm2RgzDKwBVo2qswp4xH79JHCpWMH3q4A1xpghY8xuoMnuL2afdptL7D6w+7wm0RjGmAFjTHgXWwm6Tz9l9IEVL2/t7QTgfbMmpXUct9PBvR9fzuQKD598aD0b9hw9of722OmAGuomfp0q11gwtYIDXYORLN6FQjLiMwPYH/W+xS6LWccWgm6gNkHbeOW1QFeUmESPFW8MROQsEdkCvAPcHNU+gojcJCIbRGRDW1tbEretKIXPW/s6KXE7WFSfPrdbmOnVJfzkc2czudzDdT94jQdfbiYQDI2rr93t/YiMRIMVMgumVgCFF3SQjPjE2j48+styvDoTVZ5wHsaY9caYpcD7gS+JyHG2uDHmAWPMCmPMirq6uhhdKWoBFR9v7e1k2YyatAUbjGbW5DKevvU8Llw4lTuf2cY19/2RP7zXlvLRCzsO9zJrUhklBXqWTzTz64pXfFqAWVHvZwKjVw0jdez1lmrgaIK28crbgZqoNZvoseKNEcEYsw3oB05J4r4UG13rKU66B/y8c6Cbs+fXZnTc6lI3P7jhDO79+HI6+oa54eHXuea+P/H0pgMMBcZ2LRljeHNvJ+9L8zpVrjB3ShkOgV0Ftu6TjPi8ATTaUWgerACCtaPqrAVutF9fC7xgrK8ya4HVdqTaPKAReD1en3abF+0+sPt8OtEYdh8uABGZA5wM7En6CShKkfJKUzshAxfYxx9kEhHh6mX1rPuni/j3vziVroFh/m7NJj7w7y/wf369PXLEQyx2HOmltXeI98+dnMEZZw+vy8nc2nLeO1JY4jNmVJgxJiAitwLPAU7gYWPMFhG5A9hgjFkLPAQ8KiJNWNbIarvtFhF5AtgKBIBbjDFBgFh92kPeBqwRkTuBjXbfxBsDOA+4XUT8QAj4a2NM+vK4FzB6nHZx8fLONipLXJw+K3sWhNfl5LozZ/OxFbN4pamdH6/fyw9ebub+l3Zx5rzJfGT5DP7s1Ppjjk1Y8/p+nA5h5SnTszbvTLO4vop3D6b3WPJMk1RIsjHmWeDZUWVfiXrtAz4ap+1dwF3J9GmXN2NFw40ujzmGMeZR4NExb0JJgIpOsREKGdbtaOPc+VNwZWi9JxEOh3DBwjouWFjH4W4fT73VwlNvtnDbU+/wlae3cOXS6VywsI4jPT4efW0vf7liJlMKPLtBNIumV/LMO4foHwpQ7i2MnSSFcReKoqTEm/s6Odzj46pTc896mF5dwi0XL+CvL5rPxv1d/OytFn65+RBr7Q2q5zdO4V+uHr3VsLBZVF8FWC7H5bPTGxafKVR8lAhq/xQPv9p8EK/LwaWLp2V7KnEREZbPnsTy2ZP46oeWsv/oAG6ng5mTSgv6DJ9YhJO+bj+k4qMUELrUU1z4gyGeeecwlyyaSkWeuHDcTgcNdshxMTJzUimVXhfbDvVkeyoTRvadvYqiZJTntx6hvW+Ia8+Yme2pKEkiIiyqr2T7YRUfpQBRC6g4+NGre5g5qZSLTp6a7akoKbBoehXbD/UWTFSqio+iaz1FxLZDPbzWfJRPnDUHp6O41k3ynUX1lfQOBWixk6rmOyo+ilJE/Ofvd1LpdXHdmbPGrqzkFIumWxFv2w/3ZnkmE4OKjxKF2kCFzLZDPfz63cN8+ty51JR5sj0dJUVOjkS8Fca6j4qPUjA+ZCU+xhi+8ew2Kr0uPnPevGxPRxkHFV4Xc2rL2FYgQQcqPoraO0XAc1sO8/LOdr54xUK1evKYRdMr2XZI3W5KgaEGUGHS4/Nzxy+3smh6JZ88e062p6OcAIvrq9jT0c/A8HFHluUdKj6KUuB89ektHOkd4ht/cWpO5HFTxs/i+iqMKYygA/1NVNTiKWCe3nSAn288wN9e0lgwaVmKmSV2jrdCyHSg4qNEUA0qLN490M3tT73DijmTuOXi+dmejjIBhNPsbC+AdR8VH0UpQFp7fXz+RxuoKXNz3/XL1d1WIITT7KjloxQEIfW7FRTdg34+/cM36Brw84MbVjC1siTbU1ImkMX1VWw/3EsolN9/tyo+SsTfphqU//QNBfjUD1/nvSO93PeJ5ZwyozrbU1ImmEXTq+grgDQ7Kj6KrvUUCD0+P5/54Ru83dLNf123nIsXaeLQQmRxvZXpYGueu95UfBR1uxUArb0+Pvb913hrXyff/djprDwl904oVSaGk6dXIpL/EW/5cZKUklbC4mPUBspL9rT3c8PDr9PeN8RDn3o/Fy6sy/aUlDRS5nExr7ZcxUfJf9TwyV9e2dnOLT95C4fATz5/NqfPqsn2lJQMsLi+ircPdGV7GieEut0UFZ88xBjDgy83c8PD65leVcIvbjlXhaeIWFxfyf6jg/T6/NmeyrhRy0cZcbupCOUF/UMB/uXn7/CLTQdZuXQ63/7L0yj36p9yMbHYznSw43AvK+ZOzvJsxof+xiq60pNHvHugm795bCN7Ovr54uULufXiBTj0RNKiY3FUmh0VHyVv0fN8ch9jDA+9spu7f7Od2nIvP/nc2Zwzvzbb01KyRH11CVUlLrbmcZodFR+F8EZplaDcpLXHx21Pvc2LO9q4fMk0vvmRZUwq1zN5ihkRYXF9VV5HvKn4KEqOYozhF5sO8LW1W/H5g9yxaimfPHsOIupmUyzX2+Nv7CcYMjjz0PWaVLSbiKwUkR0i0iQit8e47hWRx+3r60VkbtS1L9nlO0TkyrH6FJF5dh877T49icYQkctF5E0Recf+/5LxPgxFyRWO9FiJQb/w+GYap1bw6787nxvOmavCo0RYUl/FoD/I3o7+bE9lXIwpPiLiBO4FrgKWANeJyJJR1T4LdBpjFgD3AHfbbZcAq4GlwErgPhFxjtHn3cA9xphGoNPuO+4YQDvwIWPMqcCNwKOpPQIljK79ZB9jDE+92cLl33mJl3e2869XL+bx/3UODXUV2Z6akmOMBB3k57pPMpbPmUCTMabZGDMMrAFWjaqzCnjEfv0kcKlYX9FWAWuMMUPGmN1Ak91fzD7tNpfYfWD3eU2iMYwxG40xB+3yLUCJiHiTfQCKkiuEMxX8w083s3BaJb/5+wv43PkNeelSUdJP47QKnA7J23WfZNZ8ZgD7o963AGfFq2OMCYhIN1Brl782qu0M+3WsPmuBLmNMIEb9eGO0R/XzEWCjMWZo9E2IyE3ATQCzZ89OfMeKkkGGAkG+/1Iz//1iEx6ng699aAmfPGeuio6SkBK3k4Yp+ZtmJxnxifUXMNo/E69OvPJYFlei+mPOQ0SWYrnirohRD2PMA8ADACtWrFD/kpIT/KmpnX99+l2a2/r54LJ6vvzBJUyr0vN3lORYXF/Fm3s7sz2NcZGM+LQAs6LezwQOxqnTIiIuoBo4OkbbWOXtQI2IuGzrJ7p+vDEQkZnAz4EbjDG7krgnRckqbb1DfOPZbfx84wFmTy7jkc+cqQlBlZRZVF/J2s0H6R7wU13mzvZ0UiKZNZ83gEY7Cs2DFUCwdlSdtViL/QDXAi8Ya/V6LbDajlSbBzQCr8fr027zot0Hdp9PJxpDRGqAZ4AvGWP+mMrNK0qm8QdDPPzKbi799jp+9fZB/vaSBfz2Cxeo8CjjIhJ0cDj/XG9jWj72+sqtwHOAE3jYGLNFRO4ANhhj1gIPAY+KSBOWNbLabrtFRJ4AtgIB4BZjTBAgVp/2kLcBa0TkTmCj3TfxxgBuBRYAXxaRL9tlVxhjWsf3SIoXDXZLL394r407frWVptY+zm+cwlc/tJQFUzWKTRk/S6LS7JzdkF8ZL5LaZGqMeRZ4dlTZV6Je+4CPxml7F3BXMn3a5c1Y0XCjy2OOYYy5E7hzzJtQlCyxp72fO5/Zxu+2HWFObRk/uGEFly2eqnt2lBNmaqWXyeWevAw60AwHipIm+oYC/PcLTTz8ym7cTuG2lYv4zHlz8bqc2Z6aUiBYaXYq2X44//b6qPgoEfQk04khFDL8bOMB7v7Ndtp6h/jI8pnctvJkpmoUm5IGFk6rZM3r+wmFTF5lOFfxUZQJ5I9N7Xzj2W1sOdjDabNqeOCTZ/C+2ZOyPS2lgGmcWsmgP8iBrkFmTS7L9nSSRsVHUSaA7Yd7+Pdnt/PSe23MqCnlux87nQ+fdlJefRNV8pPGaVbQSlNrn4qPkp9otFvqHO728Z3nd/Dkmy1UeF38y58t5pPnzKHEres6SmZotCMm3zvSy8WLpmZ5Nsmj4qMo46DX5+f7LzXz4CvNhELw2fPmccvFC6gp03N2lMxSU+ahrtLLzta+bE8lJVR8FCUF/MEQj72+j//43U46+of58Gkn8U9XnpxX7g6l8Fg4rULFR8k8m/d3sereP/LyP1+c8odg9DEK6naLTzBk+OXmg9zzu/fY2zHAWfMm88OrF7NsZk22p6YoNE6t5Kcb9mOMyZv9Yyo+BcDjG6wE4evea+OTZ89JqW1IBSchxhie33qEb//2PXYc6WVxfRUP3biCSxbpJlEld1gwtYL+4SAHu33MqCnN9nSSQsWnAAjZCuIaR2RVIBSa6OkUDH9sauebz+1g8/4u5k0p57+uex9Xn1qvEWxKzrFwWiUAO4/0qvgomSNgi894zn8JBKPcbhM2o/zmrX2dfOu5HfxpVwcnVZdw90dO5SPLZ+JyJnXqvKJknHDE284jfVx0cn5EvKn4FACBoGW9jM/yUckJs+1QD9/+7Xv8btsRass9fPVDS7juzNkaNq3kPJPKPUyp8LKzNX/S7Kj4FAAnYvkEVXxoau3jP3+/k1++fZAKr4t/vGIhnz53HuVe/fNQ8ofGqfkV8aZ/XQVAMLLmk7pbKGw1wbGRb8VAU2sf//XCTtZuPkiJy8nNF87nf13QoHt1lLykcVoFP3/rQN5EvKn4FAAjls/42xYTo0XnpgsauOn8BmorvNmemqKMm8ZplfQOBTjc46O+OveDDlR8CoCw5TMeHSkmt5uKjlLIRAcdqPgoGSFsvYxHSKItn0KVIRUdpRiIiE9rHxfkwbHsKj4FQNDeqxMax5pN9JpPoaGioxQTk8s9VJe6aW7Lj6ADFZ8CIHgCls9wAYrPOy3d3P/SLp5995CKjlI0iAgNdeU0t/VneypJoeJTAJyI+Pj8UeKTx343Ywyv7urgey/t4uWd7VR6Xdx84Xw+d948FR2laGiYUsHLO9uyPY2kUPEpAAKRgIPU1WPIH5zo6WSUUMjw262H+d66XWxu6WZKhZfbr1rEx8+aTVWJO9vTU5SM0lBXzlNvtdDr81OZ47//Kj4FwIjlk3pbXyA/xWc4EOIXmw5w/0u7aG7rZ/bkMu7681P4yPKZmpFAKVrm15UDsLu9P+czrqv4FADh/GzjsXyi3W4mD/xu/UMBHnt9Hw++vJvDPT6W1FfxX9e9j6tOma6515Sip6HOinhrblPxUTJA8ATcbr48cbsd7R/mkT/t4ZFX99A14OeseZO5+9plXNA4JS92cytKJphTW4ZDyIuINxWfAsBvh1qfcMBBDnKga5AHX25mzev7GfQHuXzJNG6+cD5nzJmU7akpSs7hdTmZOamMXe25H/Gm4lMA+IYt62V84jNi+eRSaredR3q5/6Vmnt50AIBVp8/g5gsbaLTPLVEUJTb5Em6t4lMA+ALj32TaPxSY6OmcEBv3dfK9dbv47dYjlLgdXH/2HD53/jxmTkrteHBFKVYaplTwWnMHoZDJ6YMPk1qhFZGVIrJDRJpE5PYY170i8rh9fb2IzI269iW7fIeIXDlWnyIyz+5jp92nJ9EYIlIrIi+KSJ+I/Pd4H0Q+MxixfFJv2zXon+DZpI4xhpfea2P1A6/y5/f9ifW7j/K3lzbyp9sv5WsfXqrCoygpMH9qOT5/iEM9vmxPJSFjWj4i4gTuBS4HWoA3RGStMWZrVLXPAp3GmAUishq4G/iYiCwBVgNLgZOA34nIQrtNvD7vBu4xxqwRkfvtvr8XbwzAB3wZOMX+V3QM2q6z8Vg+3YN+HGIlJc201y0YMvz63UN8b90uthzsYXpVCf969WKuO3O2nqWjKOOkYUo44q0vp4/UTsbyORNoMsY0G2OGgTXAqlF1VgGP2K+fBC4VKwRpFbDGGDNkjNkNNNn9xezTbnOJ3Qd2n9ckGsMY02+MeQVLhIqO4cCIuRMax5pP96Cf6tLMbkbz+YP8ZP0+Lvn2Om79yUYGh4N88yPLeOmfL+Jz5zeo8CjKCRDe65Pr6z7J/JXPAPZHvW8BzopXxxgTEJFuoNYuf21U2xn261h91gJdxphAjPrxxmhP4h4Klq6B4cjr4HgsnwFLfDoH0u9+6/X5+fH6fTz0ym7aeodYNrOa+69fzuVLpo/rFFZFUY6nrtJLhdeV8+HWyYhPrE+F0Z9y8erEK49lcSWqn+w84iIiNwE3AcyePTvZZjlPa+9Q5PV4LB/r4KkS9nQMpC3arb1viB/+cTc/enUvvb4A5y2Ywn987HTOmV+re3QUZYKJJBjN8XDrZMSnBZgV9X4mcDBOnRYRcQHVwNEx2sYqbwdqRMRlWz/R9eONkRTGmAeABwBWrFiRQ0HFJ0Zb34j4pGr5hEKGQ92DrEjTnpn9Rwd44A/NPLFhP8PBEFedMp2bL5yf8zuvFSXfaZhSzht7OrM9jYQkIz5vAI0iMg84gBVA8PFRddYCNwKvAtcCLxhjjIisBX4iIt/BCjhoBF7HsmKO69Nu86Ldxxq7z6cTjTG+2y4cDnWNLHWlGu3W3jeEP2g4aYIXJXce6eV763bx9OaDOAQ+snwmN13QEEn9oShKemmoq+AXmw4yOByk1JObuQ7HFB97feVW4DnACTxsjNkiIncAG4wxa4GHgEdFpAnLGlltt90iIk8AW4EAcIsxJggQq097yNuANSJyJ7DR7pt4Y9h97QGqAI+IXANcMSoar2DZ2dpLqdtJyJiUo93eO2L5hBvsBcrxRMtFs3l/F/eta+K5LUcodTv59Afm8rnzG5heXXJC/SqKkhoNUQlGl5xUleXZxCapsCJjzLPAs6PKvhL12gd8NE7bu4C7kunTLm/GioYbXZ5ojLkJb6CA2Xaoh8ZpFexq7Us5w8GWg90AETfYeMTHGMNrzUe5b10TL+9sp6rExd9e2sinPjCXyeWelPtTFOXECYdb72rry2/xUXKTgeEAb+3t4sYPzGF3e3/K4rNhbyezJpcypcISiXB27GR5dVcH3/7tDjbs7WRKhZcv2efo5Po5IopS6Mybkvvh1io+ecxzWw4zHAxx8clTeWJDC6ksgfn8QV7Z2c61Z8yMhDkna/m8ta+T7/z2PV5pamdalZc7Vi3lL1fM0nN0FCVHKPU4mVFTSnN77oZbq/jkKYFgiO+/1EzDlHLObqjF6ZCUot3Wbj7IoD/IVaeO7LEJjGE57e3o585ntvH81iPUlnv48geX8ImzZqvoKEoOkusJRlV88pQf/nEP2w/3ct8nluNwCA6RpKPdhgJBvrduF4umV3JOQy1DgcRHMgwMB7j3xSZ+8IfduJ3CP16xkE+fO08zEShKDtMwpZwn37Q8Irm4n04/PfKQ9c0d3P2b7Vy+ZBpXnTIdAKcj+U2m3/3dTna39/Ojz5yJiOAKu91itH9zbyf/8MQm9nQM8Bfvm8FtVy1iWpVGrylKrtNQV0H/cJDW3qGc/JtV8ckzdh7p5a9+/Baza8v49l+eFvlG45Tk3G5PvdnC99btYvX7Z3HBwjqrbQy3mzGGh/+4h7ue2Up9dSk/+fxZfGD+lDTckaIo6SAcbr2rrU/FRzkxmlp7ue4H63E6hIdufD9VUVFlDoeMafk8+tpevvL0u5zTUMsdq0YSgIuIndnaam+M4a5ntvHgK7u5cuk0/u9HTztmLEVRcp/wpu7mtv6c/OKo4pMnvLqrg7/68Zu4HA4e+/zZkVDKMIkCDroH/Xxt7RZ+vvEAly2eyn9/fDkel+O49mHL53sv7eLBV3Zz4zlz+OqHlub0gVSKosSmvqqEErcjZ4MOVHxyHGMMj762lzt+uZU5tWU8/Kn3M6e2/Lh6VsDB8eLzh/fa+Ocn36atb4i/u7SRWy9ZgNt5fF5Xh1iW0+b9XXzruR18cFm9Co+i5DEOh9AwpSJnw61VfHKY1l4f//zk26zb0cbFJ9fxH9e9L677yyEck5W6fyjAv/96G//z2j4WTK3ggRvOSJjQ02VbPv/2q63UVXq5689PVeFRlDynoa6czS1d2Z5GTFR8cpBAMMSP1+/j27/dwVAgxB2rlvLJs+ckDJd0OkYsn7f2dfKFxzex7+gAnztvHv945clj7sVxOITXdx/lnQPd3LFqacYPmFMUZeJpqKvgmXcO4fMHc24/nopPDmGM4ZWmdu56ZhvbD/dy3oIpfO3DS1kwdexs0A472u3pTQf4hyc2M726hDWfP5uzGmqTGtvlEN450I3X5eDaM2ae6K0oipIDzK8rxxjY2zHAydMrsz2dY1DxyQGMMfyxqYPv/u49NuztZOakUu6//gyuXDot6c1hTofw3pFe/vGnrayYO4nvf3JFStZLONz6nPm1lHn010JRCoFwgtHmtj4VH2WEUMjwwvZW7n9pFxv2dlJfXcK/XXMKf7liJl5Xaiay0yHs7Rigwuvi3o8vT9lt1t5nHcd9yaKpKbVTFCV3mWfv9cnFU01VfLJA31CAn27Yz//70x72dgxw0gmITpjwes/Vp9ZTW+Ed99wuPlnFR1EKhQqvi2lVXna15V7Em4pPBtnXMcD/+9MefrphP71DAc6YM4l/uvJkrlw6PWb4cypsOdgDwOVLpo2r/Zc/uIQtB7qZNbnshOahKEpu0TClIif3+qj4pJlQyPDSzjYefXUvL+5oxSnC1cvq+fS58zh9VvzQ51RZdfpJPL3pIOcuGN9O5s+eN2/C5qIoSu7QUFfOLzcfzLkEoyo+aaJrYJifbmjhf9bvZW/HAFMqvPzNxQv4xNlz0pJn6Y5Vp/D3ly3M2fPaFUXJDg11FfT4AnT0DzPlBFzyE42KzwSzu72fB/6wi5+9dYChQIj3z53EP1xxMiuXTj8upc1EUl3q1r05iqIcRzjBaHNbv4pPIdLU2ss9z+/k2XcP4XY6+MjyGdxwzlwW1+fm+emKohQHC+wEo7va+jhz3uQsz2YEFZ8TpG8owLee28Gjr+2l1O3k5gvn85lz51FXmTvfMBRFKV5OqinF43LQnGMRbyo+J8DWgz3c8pO32NvRz3VnzuaLly88oTBnRVGUicbpEObV5t6R2io+42Tjvk5uePh1yjxOHkshjY2iKEqmaagrZ/vh3mxP4xjStwJewLT2+vj8j95kUpmHn/31uSo8iqLkNA115ew7OsBwIJTtqURQ8RkHX1+7lb4hPw/euIIZNaXZno6iKEpCFkytIBgy7M6hNDsqPimy43Avz7xziJvOb2DhtNxK1KcoihKLpSdVA/Duge4sz2QEFZ8U+dlbLbgcwqfP1YwAiqLkB/PrKih1O3lHxSd/eem9Ns6ZX8ukck+2p6IoipIUToew9KSq/LN8RGSliOwQkSYRuT3Gda+IPG5fXy8ic6Oufcku3yEiV47Vp4jMs/vYaffpGe8YE43PH2Rna9+E5mRTFEXJBKfMqGbLwZ5IBvxsM6b4iIgTuBe4ClgCXCciS0ZV+yzQaYxZANwD3G23XQKsBpYCK4H7RMQ5Rp93A/cYYxqBTrvvlMdI9UEkQ1vvEMGQ0czPiqLkHSvmTmLQH+StfZ3ZngqQnOVzJtBkjGk2xgwDa4BVo+qsAh6xXz8JXCpW+tRVwBpjzJAxZjfQZPcXs0+7zSV2H9h9XjPOMSac7kE/gOZQUxQl77hwYR0ep4NHX92b7akAyW0ynQHsj3rfApwVr44xJiAi3UCtXf7aqLYz7Nex+qwFuowxgRj1xzNGBBG5CbjJftsnIh1Ae9y7TsDKu8fTKqeZwjifRQGiz8JCn8MIBfUsdgL/9fFxNZ0CzJmoeSQjPrEOgBjtNIxXJ155LIsrUf3xjHFsgTEPAA+E34vIBmPMihhtiw59FiPos7DQ5zCCPgsL+znMnaj+knG7tQCzot7PBA7GqyMiLqAaOJqgbbzydqDG7mP0WKmOoSiKouQoyYjPG0CjHYXmwVrcXzuqzlrgRvv1tcALxhhjl6+2I9XmAY3A6/H6tNu8aPeB3efT4xxDURRFyVHGdLvZ6yu3As8BTuBhY8wWEbkD2GCMWQs8BDwqIk1Y1shqu+0WEXkC2AoEgFuMMUGAWH3aQ94GrBGRO4GNdt+MZ4wxeGDsKkWDPosR9FlY6HMYQZ+FxYQ+B7GMB0VRFEXJHJrhQFEURck4Kj6KoihKxilK8RkrXVAhICIPi0iriLwbVTZZRJ63Uxc9LyKT7HIRkf+0n8fbIrI8qs2Ndv2dInJjrLFyGRGZJSIvisg2EdkiIn9nlxfVsxCREhF5XUQ228/h63Z5zqazSjd2tpWNIvIr+31RPgsR2SMi74jIJhHZYJel/+/DGFNU/7ACHHYBDYAH2Awsyfa80nCfFwDLgXejyr4J3G6/vh242379Z8CvsfZMnQ2st8snA832/5Ps15OyfW8pPod6YLn9uhJ4DyulU1E9C/t+KuzXbmC9fX9PAKvt8vuBv7Jf/zVwv/16NfC4/XqJ/TfjBebZf0vObN/fOJ/JF4GfAL+y3xflswD2AFNGlaX976MYLZ9k0gXlPcaYP2BFBUYTnaJodOqiHxmL17D2WtUDVwLPG2OOGmM6geex8uflDcaYQ8aYt+zXvcA2rAwYRfUs7Pvps9+67X+GAOvoTgAAAk1JREFUHE5nlU5EZCZwNfCg/T6nU3tlgbT/fRSj+MRKF3RcOp4CZZox5hBYH8rAVLs83jMpqGdlu0veh/Wtv+iehe1m2gS0Yn047CLJdFZAdDqrvH4ONt8F/hkInyuddGovCu9ZGOC3IvKmWGnIIAN/H8mk1yk0kkrHU2ScUOqifEBEKoCngL83xvRYX1xjV41RVhDPwlj7304XkRrg58DiWNXs/wv2OYjIB4FWY8ybInJRuDhG1YJ/FjbnGmMOishU4HkR2Z6g7oQ9i2K0fIo5Hc8R20TG/r/VLk81DVJeISJuLOH5sTHmZ3ZxUT4LAGNMF7AOy2dfjOmszgU+LCJ7sNzul2BZQsX4LDDGHLT/b8X6UnImGfj7KEbxSSZdUKESnaJodOqiG+xIlrOBbtvUfg64QkQm2dEuV9hleYPtm38I2GaM+U7UpaJ6FiJSZ1s8iEgpcBnW+lfRpbMyxnzJGDPTWEkyV2Pd2ycowmchIuUiUhl+jfV7/S6Z+PvIdqRFNv5hRWy8h+Xz/pdszydN9/gYcAjwY30r+SyWn/r3WFnVfw9MtusK1uF+u4B3gBVR/XwGayG1Cfh0tu9rHM/hPCzz/21gk/3vz4rtWQDLsNJVvW1/uHzFLm/A+sBsAn4KeO3yEvt9k329Iaqvf7Gfzw7gqmzf2wk+l4sYiXYrumdh3/Nm+9+W8OdhJv4+NL2OoiiKknGK0e2mKIqiZBkVH0VRFCXjqPgoiqIoGUfFR1EURck4Kj6KoihKxlHxURRFUTKOio+iKIqScf4/2QYEJG6Qy5oAAAAASUVORK5CYII=\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], + "outputs": [], "source": [ "plt.clf()\n", "# plt.plot(x_part, calcs, '.')\n", @@ -743,7 +699,7 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -759,7 +715,7 @@ }, { "cell_type": "code", - "execution_count": 15, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -775,7 +731,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -788,7 +744,7 @@ }, { "cell_type": "code", - "execution_count": 17, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -845,7 +801,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -929,7 +885,7 @@ }, { "cell_type": "code", - "execution_count": 19, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -938,7 +894,7 @@ }, { "cell_type": "code", - "execution_count": 20, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -947,7 +903,7 @@ }, { "cell_type": "code", - "execution_count": 21, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -956,7 +912,7 @@ }, { "cell_type": "code", - "execution_count": 22, + "execution_count": null, "metadata": { "scrolled": false }, @@ -1004,7 +960,7 @@ }, { "cell_type": "code", - "execution_count": 23, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1021,7 +977,7 @@ }, { "cell_type": "code", - "execution_count": 24, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1045,7 +1001,7 @@ }, { "cell_type": "code", - "execution_count": 25, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1072,7 +1028,7 @@ }, { "cell_type": "code", - "execution_count": 26, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1095,7 +1051,7 @@ }, { "cell_type": "code", - "execution_count": 27, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1104,7 +1060,7 @@ }, { "cell_type": "code", - "execution_count": 28, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1120,7 +1076,7 @@ }, { "cell_type": "code", - "execution_count": 29, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1150,7 +1106,7 @@ }, { "cell_type": "code", - "execution_count": 30, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1164,7 +1120,7 @@ }, { "cell_type": "code", - "execution_count": 31, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1182,7 +1138,7 @@ }, { "cell_type": "code", - "execution_count": 32, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1196,7 +1152,7 @@ }, { "cell_type": "code", - "execution_count": 33, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1217,7 +1173,7 @@ }, { "cell_type": "code", - "execution_count": 34, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1227,7 +1183,7 @@ }, { "cell_type": "code", - "execution_count": 35, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1258,7 +1214,7 @@ }, { "cell_type": "code", - "execution_count": 36, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1275,7 +1231,7 @@ }, { "cell_type": "code", - "execution_count": 43, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1289,7 +1245,7 @@ "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(DDstar_m,2))))\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", @@ -1308,7 +1264,64 @@ "\n", "# 3. Constraint - Maximum eta of D contribs\n", "\n", - "constraint3_0 = tf.cond(tf.greater_equal(sum_2, 5.0e-8), lambda: 10000., lambda: 0.)" + "constraint3_0 = tf.cond(tf.greater_equal(tf.abs(Dbar_s), 0.2), lambda: 10000., lambda: 0.)\n", + "\n", + "constraint3_1 = tf.cond(tf.greater_equal(tf.abs(DDstar_s), 0.2), lambda: 10000., 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)\n", + " sigma2 = ztf.constant(0.128)\n", + " sigma3 = ztf.constant(0.548)\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", + "# 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), 0.02), lambda: 10000., lambda: 0.)\n", + "\n", + "#List of all constraints\n", + "\n", + "constraints = [constraint1, constraint2, constraint3_0, constraint3_1, constraint4]" ] }, { @@ -1320,46 +1333,18 @@ }, { "cell_type": "code", - "execution_count": 41, + "execution_count": null, "metadata": { "scrolled": false }, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "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" - ] - }, - { - "ename": "RuntimeError", - "evalue": "exception was raised in user function\nUser function arguments:\n p3770_s = +1.999020\n p4160_s = +2.376152\n p4040_p = -4.437832\n p4040_s = +1.423904\n Dbar_s = +0.342445\n p4160_p = +4.897009\n p4415_p = +4.364139\n p4415_s = +1.177263\n Ctt = -0.212905\n bplus_0 = -0.141110\n psi2s_p = -2.316177\n p3770_p = +0.511010\n bplus_2 = -1.674701\n jpsi_p = +4.878201\n bplus_1 = +1.591894\nOriginal python exception in user function:\nInternalError: Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\nCaused by op 'model_2/_call_unnormalized_pdf/Sqrt_27', 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 \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 \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 365, 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 272, 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 542, 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 \"\", line 71, in \n result = minimizer.minimize(nll)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 207, in minimize\n return self._hook_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 216, in _hook_minimize\n return self._call_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 220, in _call_minimize\n return self._minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 48, in _minimize\n gradients = loss.gradients(params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 158, in gradients\n return self._gradients(params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 251, in _gradients\n gradients = {(p, grad) for p, grad in zip(params_todo, super()._gradients(params_todo))}\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 220, in _gradients\n return tf.gradients(self.value(), params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 198, in value\n return self._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 234, in _value\n loss = super()._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 203, in _value\n constraints=self.constraints)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 266, in _loss_func\n nll = _unbinned_nll_tf(model=model, data=data, fit_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in _unbinned_nll_tf\n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in \n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 43, in _unbinned_nll_tf\n probs = model.pdf(data, norm_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 314, in pdf\n value = self._single_hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 324, in _single_hook_pdf\n return self._hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 327, in _hook_pdf\n return self._norm_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 330, in _norm_pdf\n return self._call_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 338, in _call_pdf\n return self._fallback_pdf(x=x, norm_range=norm_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 341, in _fallback_pdf\n pdf = self._call_unnormalized_pdf(x, name=\"_call_unnormalized_pdf\")\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 289, in _call_unnormalized_pdf\n return self._unnormalized_pdf(x)\n File \"\", line 73, in _unnormalized_pdf\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 File \"\", line 70, in ttau_cusp\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 File \"\", line 15, in h_S\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 File \"\", line 9, in G\n big_bracket = tf.where(tf.math.real(y) > ztf.constant(0.0), inner_rect_bracket(y), inner_right(y))\n File \"\", line 4, in inner_rect_bracket\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 File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_math_ops.py\", line 9251, in sqrt\n \"Sqrt\", x=x, 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\nInternalError (see above for traceback): Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 116, in grad_func\n gradients_values = self.sess.run(gradients)\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 1348, in _do_call\n raise type(e)(node_def, op, message)\n", - "output_type": "error", - "traceback": [ - "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", - "\u001b[1;31mRuntimeError\u001b[0m Traceback (most recent call last)", - "\u001b[1;32m\u001b[0m in \u001b[0;36m\u001b[1;34m\u001b[0m\n\u001b[0;32m 69\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 70\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---> 71\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 72\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 73\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Toy {}: Fitting finished\"\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mformat\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mtoy\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;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 p3770_s = +1.999020\n p4160_s = +2.376152\n p4040_p = -4.437832\n p4040_s = +1.423904\n Dbar_s = +0.342445\n p4160_p = +4.897009\n p4415_p = +4.364139\n p4415_s = +1.177263\n Ctt = -0.212905\n bplus_0 = -0.141110\n psi2s_p = -2.316177\n p3770_p = +0.511010\n bplus_2 = -1.674701\n jpsi_p = +4.878201\n bplus_1 = +1.591894\nOriginal python exception in user function:\nInternalError: Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\nCaused by op 'model_2/_call_unnormalized_pdf/Sqrt_27', 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 \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 \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 365, 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 272, 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 542, 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 \"\", line 71, in \n result = minimizer.minimize(nll)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 207, in minimize\n return self._hook_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 216, in _hook_minimize\n return self._call_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 220, in _call_minimize\n return self._minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 48, in _minimize\n gradients = loss.gradients(params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 158, in gradients\n return self._gradients(params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 251, in _gradients\n gradients = {(p, grad) for p, grad in zip(params_todo, super()._gradients(params_todo))}\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 220, in _gradients\n return tf.gradients(self.value(), params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 198, in value\n return self._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 234, in _value\n loss = super()._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 203, in _value\n constraints=self.constraints)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 266, in _loss_func\n nll = _unbinned_nll_tf(model=model, data=data, fit_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in _unbinned_nll_tf\n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in \n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 43, in _unbinned_nll_tf\n probs = model.pdf(data, norm_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 314, in pdf\n value = self._single_hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 324, in _single_hook_pdf\n return self._hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 327, in _hook_pdf\n return self._norm_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 330, in _norm_pdf\n return self._call_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 338, in _call_pdf\n return self._fallback_pdf(x=x, norm_range=norm_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 341, in _fallback_pdf\n pdf = self._call_unnormalized_pdf(x, name=\"_call_unnormalized_pdf\")\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 289, in _call_unnormalized_pdf\n return self._unnormalized_pdf(x)\n File \"\", line 73, in _unnormalized_pdf\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 File \"\", line 70, in ttau_cusp\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 File \"\", line 15, in h_S\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 File \"\", line 9, in G\n big_bracket = tf.where(tf.math.real(y) > ztf.constant(0.0), inner_rect_bracket(y), inner_right(y))\n File \"\", line 4, in inner_rect_bracket\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 File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_math_ops.py\", line 9251, in sqrt\n \"Sqrt\", x=x, 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\nInternalError (see above for traceback): Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 116, in grad_func\n gradients_values = self.sess.run(gradients)\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 1348, in _do_call\n raise type(e)(node_def, op, message)\n" - ] - } - ], + "outputs": [], "source": [ "# zfit.run.numeric_checks = False \n", "\n", "Ctt_list = []\n", "Ctt_error_list = []\n", "\n", - "nr_of_toys = 1\n", + "nr_of_toys = 2\n", "nevents = int(pdg[\"number_of_decays\"])\n", "nevents = pdg[\"number_of_decays\"]\n", "event_stack = 1000000\n", @@ -1420,7 +1405,7 @@ " 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 = constraint1)\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", diff --git a/__pycache__/pdg_const.cpython-37.pyc b/__pycache__/pdg_const.cpython-37.pyc index 8eb1fd0..e91547c 100644 --- a/__pycache__/pdg_const.cpython-37.pyc +++ b/__pycache__/pdg_const.cpython-37.pyc Binary files differ diff --git a/data/zfit_toys/toy_0/0.pkl b/data/zfit_toys/toy_0/0.pkl index 0c3b19e..b09e1b1 100644 --- a/data/zfit_toys/toy_0/0.pkl +++ b/data/zfit_toys/toy_0/0.pkl Binary files differ diff --git a/data/zfit_toys/toy_0/1.pkl b/data/zfit_toys/toy_0/1.pkl index c21955a..fa9cfc5 100644 --- a/data/zfit_toys/toy_0/1.pkl +++ b/data/zfit_toys/toy_0/1.pkl Binary files differ diff --git a/data/zfit_toys/toy_0/2.pkl b/data/zfit_toys/toy_0/2.pkl index 391e9c8..68e9674 100644 --- a/data/zfit_toys/toy_0/2.pkl +++ b/data/zfit_toys/toy_0/2.pkl Binary files differ diff --git a/raremodel-nb.ipynb b/raremodel-nb.ipynb index 47c3623..529e151 100644 --- a/raremodel-nb.ipynb +++ b/raremodel-nb.ipynb @@ -9,31 +9,9 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "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" - ] - } - ], + "outputs": [], "source": [ "import os\n", "\n", @@ -64,7 +42,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -90,7 +68,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -280,7 +258,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -335,7 +313,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -354,7 +332,7 @@ " '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", - " 'DDstar_mass', 'DDstar_scale', 'DDstar_phase',\n", + " 'Dstar_mass', 'DDstar_scale', 'DDstar_phase', 'D_mass',\n", " 'tau_mass', 'C_tt']\n", "# the name of the parameters\n", "\n", @@ -404,7 +382,7 @@ " \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['DDstar_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", @@ -431,7 +409,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -459,19 +437,9 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "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" - ] - } - ], + "outputs": [], "source": [ "# formfactors\n", "\n", @@ -580,7 +548,7 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -588,16 +556,17 @@ "\n", "Dbar_phase = 0.0\n", "DDstar_phase = 0.0\n", - "DDstar_mass = pdg['D0_M']\n", - " \n", - "Dbar_mass = (pdg['D0_M']+pdg['Dst_M'])/2\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=-1.464, upper_limit=1.464)\n", "Dbar_m = zfit.Parameter(\"Dbar_m\", ztf.constant(Dbar_mass), floating = False)\n", - "Dbar_p = zfit.Parameter(\"Dbar_p\", ztf.constant(Dbar_phase), 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=-2.0, upper_limit=2.0)#, floating = False)\n", - "DDstar_m = zfit.Parameter(\"DDstar_m\", ztf.constant(DDstar_mass), floating = False)\n", - "DDstar_p = zfit.Parameter(\"DDstar_p\", ztf.constant(DDstar_phase), floating = False)" + "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)" ] }, { @@ -609,7 +578,7 @@ }, { "cell_type": "code", - "execution_count": 10, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -626,7 +595,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -639,7 +608,7 @@ " 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", - " DDstar_mass = DDstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p,\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", @@ -663,7 +632,7 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -709,22 +678,9 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ8AAAD8CAYAAACo9anUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO2dd3yc1ZX3v2equmTLsi3cZcu4gCHGoYTezZLE7IZsTEIglZdd2JJsgby7aSxkX/ImIVsghAB5CZtgCCTBCSSEBEwgAYPBNuCGZbnJTcXq0khT7vvH88xoLM+MZmRNP9/Pxx/P3Oe255E0vznnnnuuGGNQFEVRlEziyPYEFEVRlOJDxUdRFEXJOCo+iqIoSsZR8VEURVEyjoqPoiiKknFUfBRFUZSMk5T4iMhKEdkhIk0icnuM614Redy+vl5E5kZd+5JdvkNErhyrTxGZZ/ex0+7Tk8QYy0TkVRHZIiLviEjJeB6GoiiKkhnGFB8RcQL3AlcBS4DrRGTJqGqfBTqNMQuAe4C77bZLgNXAUmAlcJ+IOMfo827gHmNMI9Bp951oDBfwP8DNxpilwEWAP8XnoCiKomSQZCyfM4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FIREbt8jTFmyBizG2iy+4vZp93mErsP7D6vGWOMK4C3jTGbAYwxHcaYYPKPQFEURck0riTqzAD2R71vAc6KV8cYExCRbqDWLn9tVNsZ9utYfdYCXcaYQIz68cZYCBgReQ6owxK7b46+CRG5CbgJoLy8/IxFixYlcevFwY4jvQwHQiycVonXpcuAisXOI314XA7m1JbFrdPj87O3Y4AFdRWUepwZnN0IzW39ADTUlR9Tvv1wL/5giEllHjwuB0d6fJwyoxpJ41yO9Pho7R0CoK7Sy/SqwloBePPNN9uNMXUT0Vcy4hPrZzU6J0+8OvHKY33CJaqfaAwXcB7wfmAA+L2IvGmM+f0xFY15AHgAYMWKFWbDhg0xuitOLv7WOna39/PzL17AgqmV2Z6OkiNc9p2XWDitgvs+cUbcOut2tPKpH77Bj/7qA5wxZ1IGZzfCtd/7E163gx9/7uxjyi/59jqa2/r5+FmzmT25jP/z6+28csfKtIrkXc9s5cfr9+FxOVh12kl8fdUpaRsrG4jI3onqK5mvuS3ArKj3M4GD8erYazDVwNEEbeOVtwM1dh+jx0o0xkvGmHZjzADwLLA8iftSRqFp/pRogiGDy5H4I8LttK4HgqFMTCkmIWNwyPHfTb0uS2Q8TkdknsNpnufAcJBSt5MKr4veocDYDYqYZMTnDaDRjkLzYAUQrB1VZy1wo/36WuAFY2UsXQustiPV5gGNwOvx+rTbvGj3gd3n02OM8RywTETKbFG6ENia/CNQFCUW/mAIlyOxkyoiPqHsfXMJGpAY4uOxXcgelwOP07ruT7P4DPqDlNji0+dT8UnEmG43e33lVqwPeSfwsDFmi4jcAWwwxqwFHgIeFZEmLGtktd12i4g8gSUGAeCWcDBArD7tIW8D1ojIncBGu28SjNEpIt/BEjQDPGuMeeaEnoqiKASCBpczsfiEr6fbokiEMYZYGllii4/bKRGRTLf4+PxByjy2+Kjlk5Bk1nwwxjyL5c6KLvtK1Gsf8NE4be8C7kqmT7u8GSsabnR5ojH+Byvcetz4/X5aWlrw+Xwn0k3OUFJSwsyZM3G73Um3Ua+bEk0gZHA5x3C7OcJut+z99oSMwRnD8qnwWh9vHqdzRHwC6Z3nwHCQUo+Tcq+LroHhtI6V7yQlPsVAS0sLlZWVzJ07N6YJn08YY+jo6KClpYV58+ZlezpKnhIIJeF2c1nXs7rmE4rtdgtbZeVeJ25XZtZ8BodH3G4tnQNpHSvf0bhaG5/PR21tbd4LD1h/iLW1tQVjxSnZIRgcO+AgfD2bbrdQHLdbeG4VXlfG1nzCbrdyr5P+Id1umAgVnygKQXjCjOdeNNpNicYfCo255uN2hi2f7P3yGEPMaDeHrUhetyNjaz7haLdyXfMZExUfRVFi4g+aiLjEYyTaLcuWT4xPstNmVgMwraokY+Iz6B8Jte4fDmD0G11cVHzyjHXr1vHBD34QgKGhIS677DJOP/10Hn/88SzPTCkkgiFDMGTwOBNvyByJdstuwEEsS//GD8zlR585k3Maakf2+aQ54MDnD1JiBxwYY1lCSmw04CCP2bhxI36/n02bNk1If0bj3RSbsIUQDiiIx0i0WzYtn9huN7fTwQULrUwwHldm1nwGhoOU2W43gP6hQOS1cixq+eQQe/bsYdGiRdx4440sW7aMa6+9loGBAX7zm9+waNEizjvvPH72s58B0NrayvXXX8+mTZs4/fTT2bVrV5ZnrxQS4Q9pz1ih1q5cCbVOXGfE8kmf+BhjLLebx0mF17IYdd0nPirJMfj6L7ew9WDPhPa55KQqvvqhpWPW27FjBw899BDnnnsun/nMZ/jOd77D97//fV544QUWLFjAxz72MQCmTp3Kgw8+yLe+9S1+9atfTehcFcVvi4l7DPEJh2JnP9otsfqUuC0xGEqj+AwFQhhjjVXuCVs+6naLh1o+OcasWbM499xzAbj++uvZsGED8+bNo7GxERHh+uuvT9vYujaqhIm43cayfJw5YPnE2ecTTaktPoP+9InBoL2+E85wAGr5JEItnxgkY6Gki9F/RN3d3QUVAq7kB2H31FjRbk6H4JDsRrvFS68TjddtiaQvneJj910yas1HiY1aPjnGvn37ePXVVwF47LHHuOyyy9i9e3dkTeexxx7L5vSUIiHsRvMkcb6Ty+nIststdsBBNGHLJ53iE+67zBMlPsMqPvFQ8ckxFi9ezCOPPMKyZcs4evQoX/jCF3jggQe4+uqrOe+885gzZ07axla3mxIm2YADALdDsh5wMEYihsiaz2AaQ5/Dlo/X5aSyRN1uY6FutxzD4XBw//33H1O2cuVKtm/fflzdiy66iIsuuihDM1OKiXACzrHWfMCKeMv2eT5juabdTgcuh+ALpN/yKfWo2y0Z1PJRFOU4hiP7fJJwuzmy73aLldV6NCVuJ4PD6Zunz2/1XeJyUOYOh1prtFs8VHxyiLlz5/Luu+9mbXzdZKqEGYl2G/tD3etypDWEeSziJRYdTYnbmVbLJ+zSK/U4cTiEco9TLZ8EqPhEUUh5mArpXpTMk8qaj9ftSOvmzbEIhcZ2uwGUuB340rjmExa28PpSudel4pMAFR+bkpISOjo6CuJDO3yeT0lJSbanouQpI6HWSYiPy5lVyydeVuvRZMzyscVHTzNNjAYc2MycOZOWlhba2tqyPZUJIXySaSoUgO4qE0Sym0zBCsfOB7dbqduZ1mg3n/0MwnuK1PJJjIqPjdvt1lM/FcUmnKXaM0ZiUbDXfNK4f2YsQmbk7J5ElLgdkaCAdOAbZfnogXKJUbeboijH4Q+E13wSH6kA2Q84CBpDMklAStzOtKbX8fmPXfNRt1tiVHwURTmOZI9UAGvNJ5sBByaJxKJgr/mkOb2OyyERV2W5faCcEhsVH0VRjiOVNR+v28FQGhfyxyLZfT6laRYfnz8UcbmBrvmMhYqPoijHMZzkkQoAXmd+BByUeZz0pzm9jjdKfNTtlhgVHyWCRrspYYYDqe3zyZb4GGMwZuwjFQAqS1z0+dInBkP+IKWekedV7nHh84eymnool1HxURTlOFLLcJC9NZ/wF6Zk1nwqvG4G/cG0icGgP0iJK9rtZr3WiLfYqPgoinIc/mAIEeu8nrGwot2y8wEbstUnGbdbujNN++wjtMNEDpTToIOYqPgoETS3mxJmOBjC43Qk5c4KbzLNRnaQUNjySUJ9Kmzx6U2T6+14y0czWyciKfERkZUiskNEmkTk9hjXvSLyuH19vYjMjbr2Jbt8h4hcOVafIjLP7mOn3acn0RgiMldEBkVkk/3v2PMIFEVJmeFAKKlgA7AsH2PAn4UzfcKWTzL7fCq96RUfnz9ESSzLR8UnJmP+domIE7gXuApYAlwnIktGVfss0GmMWQDcA9xtt10CrAaWAiuB+0TEOUafdwP3GGMagU6777hj2Owyxpxu/7s5pSegKMpxDAVClLiTFR+n3Sbzrrew+CQTal1Z4gbS7HaLemZq+SQmmd+uM4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FKx7PVVwBpjzJAxZjfQZPcXs0+7zSV2H9h9XjPGGMoEodFuShifPxgRlbEI5zLLRtBBKJWAg8iajz8tc/H5g5HsBhAdcKDiE4tkxGcGsD/qfYtdFrOOMSYAdAO1CdrGK68Fuuw+Ro8VbwyAeSKyUUReEpHzk7gnRVESMBQIRURlLMLh2NkIt07J7ZbmNZ/+4SBlMd1uGu0Wi2QSi8b6sY7+jhyvTrzyWL/VieonGuMQMNsY0yEiZwC/EJGlxpieYyYochNwE8Ds2bNjdKUoSpihUYvniQiLVDbEx9hDJmP5pHvNp38oQLln5CNV3W6JSearTQswK+r9TOBgvDoi4gKqgaMJ2sYrbwdq7D5GjxVzDNul1wFgjHkT2AUsHH0TxpgHjDErjDEr6urqkrjt4kO9bkqYVCyfXFjzSSbUuiKNodahkGFgOBgRHNCAg7FI5rfrDaDRjkLzYAUQrB1VZy1wo/36WuAFY8VdrgVW25Fq84BG4PV4fdptXrT7wO7z6URjiEidHcCAiDTYYzQn/wgUXThTRuNLxfJxZXPNxxafJNSn1O3E7RS6Byd+zWfAzhlXESU+XpcDp0PU8onDmG43Y0xARG4FngOcwMPGmC0icgewwRizFngIeFREmrAsntV22y0i8gSwFQgAtxhjggCx+rSHvA1YIyJ3Ahvtvok3BnABcIeIBIAgcLMx5uj4H0nxoRaPMpqhQIjy8uSO+xqxfDIvPsHIms/Y4iMiTCrzcLRveMLnERaYMu+IYIuI5ndLQFK/XcaYZ4FnR5V9Jeq1D/honLZ3AXcl06dd3owVDTe6POYYxpingKfGvAllTArhCHFlYkjJ8rHdc+nMGB2PkfQ6ydWfXO6ho3/ixScsMNGWD0BNmZuugfRE1+U7muFAUbebchyprPmEjxEYSGPG6HgEQ8nv8wFLfI72D034PAbsiLbogAOAmjIPnQMTL3aFgIqPom435ThSsXzC4cWD2RSfJE0fS3zSZ/lEu90AJpW5VXzioOKjRFARUsKkYvmU2d/2s2n5uJLIvg1Qmybx6Y/jdptc5qGzX91usVDxUdTtphzH6N36iQhnch7IQvbmQCgcap2s5eOlxxeIHBkxUYSPyy4/bs3HQ5daPjFR8VHU4lGOwRhjWT6uZC2f7LndwqHWLkdyc51cbuV365xg6ydewMGkMjf9w8GsHjOeq6j4KBE02E0B6zgFY0ja8nE7HbidEtnrkkkCwfCaT3L1ayu8ALSPI9y6e9DP81uPEAod/4cSPiH1uGi3cg+ARrzFQMVHUbebcgw+v+WSStbyASvibSAL+1lGAg6Sm+u0Kkt8jvT4Uh7rrme28vkfbeCpt1qOu9Y54MfjdByT2w2sNR/rurreRqPio6jbTTmGsIvIm6TlA9ZaR1YCDiJut+S+QtVXlwJwsHsw5bFebe4A4NfvHj7uWvfgMNVl7uM2u04qC7v51PIZjYqPEoXKkAJDtuVTkorl43Fmxe0WDFlzTSa9DsDUSi8OgUNdqVk+3YN+9h+1BGvjvs7jNmR3DfipKXUf165GLZ+4qPgo6nZTjmE8lk+Zx5mlfT7W/8laPi6ng2lVJSlbPvuPDgBwfuMUOgf8tHQe275rwE9N2fHiM9le80lHeHe+o+KjqL2jHINvHJZPmduVpVBra67JbjIFqK8uSdnyae216l+xZBoAm1u6jrneNeiPWDnR1FZYZW29E59VId9R8VEiaLSbAiMbJkfvWUlEadYsn9QyHADU15RyKEXLp7XHEo/zGuvwOB28c6D7mOtdA8Mx3W5up4MpFR5aVXyOQ8VHUbebcgzhwIHRkVuJKPM4s5vbLQXxmVlTysEuX6RtMoTF46SaEk6eXsm7UeJjjKGjfzjiYhtNXWUJreOIrit0VHwUdbspxzAiPqlZPvmQWBRg7pRyhoMhDnYlb/209vqoKXPjdTk5ZUY17x7oiQQddA74GQ6EmFZVErPttCovR3pVfEaj4qNEUBFSYCRVTCqWT7nHxWA2NpmOw/KZX1cBwK62vqTbtPYMMbXS2iN06ozqY6LfDndbwjK9Oo74VJZE3HbKCCo+irrdlGMYGMeaT5nHmZUTO0MpJhYFaKgrB6C5rT/pNq29Q0yttMTl1BnVAJF1n/CG1XiWz9QqL+19QwQmOJ9cvqPio6jFoxxDeL9OKpZPZYmLoUAo40dpB8bhdqst91BV4qK5PXnLp613xPJZOL0Ct1Mi4hMO245n+UytKiFkSMshdvmMio8SQaPdFLAORnNIaul1quxIr15fZnfyhxOLpuJ2ExEa6irY1Zqc5WOMoa13iDo7NY/X5WTR9KpI0EFzWz8lbgf18dZ8Ksef0qeQUfFR1O2mHEP/cIByj+u4VDGJqCyxXHQ9vsy63sKJRZPNah1mcX0lWw/1JHV0fNeAn+FgKOJ2Azh1ZjWbW7oIBEM0tfYxv64ibpaFsEV0qFvFJxoVH0XdbsoxDA4HjzuRcyyqSrJj+YSj3VLUHk6xgwZGZyqIRTjMOux2A7igsY5eX4D1u4+y5WA3J0+vjNt+1qQyYCRLgmKh4qNESOZboFL49A8HKU8hzBqg0hafnsHMWj7BFM/zCbNsRg0Ab7d0j1FzJLtBXZT4XLiwjhK3g3/71Vba+4Y5e15t3PY1ZW4qvS4Vn1Go+CjqdlOOYWAoEDmdNFmqSsNut8xaPoFxWj7hoIG3D3SNWTccJh1t+ZR6nHzirDlsP9xLqdvJFUunxW0vIsycXMb+JKysYiK1rzdKQaL2jhLNwAlYPhl3u9nhy6laPl6Xk8X1VWzcN7b4tPXZ4jMqoOCfrjyZukov75tVEzOvWzSzJ5eyK4nQ7u5BP1sP9nDO/PiWVKGglo8SQUVIARgYDoxjzce2fDLudrP+TyXaLczZDbVs2tc1Zk661p4hyjzO404pLXE7ufnC+ZzVMLZQzJ5cxv6jA2O6tr/+yy1c94PXeO9I79g3kOeo+CjqdlOOoX84mNIeH7AyHDgk82634DiyWoc5Z34tw8EQb+7tTFivrW/oGJfbeJg1uYyhQGjM7Na/3XIEgE37x7bI8h0VH0UtHuUY+nyB477lj4XDIVR4XfRmONQ6nDQglU2mYc6cOxmXQ/jTrvaE9Vp7fMeEWY+H2ZOtiLfd7Yldb+G1tlSyL+QrKj5KBA12U8CyXsKh06lQVeqmZzDDAQe2+rhTSK8TptzrYvnsSby4oy1hvba+oWMi3cZD4zQrFLspQT65UMjQZZ94WgwbUlV8FHW7KRECwRADw8FIxoJUqCpx051h8fEHQ4iMz+0GcMXSaWw71MOeOBaJMYYj3T6mVp2Y+JxUXUK5x8nOI/HFp3vQj99exFLxsRGRlSKyQ0SaROT2GNe9IvK4fX29iMyNuvYlu3yHiFw5Vp8iMs/uY6fdp2esMezrs0WkT0T+MdWHUOyowaOECbvNwhkLUmFSuZvOgczmLxsOGtxOR0rZGKJZecp0AH797uGY1zsH/PQPByMbRceLiLBgWmXCQIL2vpH1IBUfQEScwL3AVcAS4DoRWTKq2meBTmPMAuAe4G677RJgNbAUWAncJyLOMfq8G7jHGNMIdNp9xx0jinuAXyd748rxGJWhoiccMDAet9vkci9HM5w80x8M4XGO34Ezc1IZp82s5pebD8aMRNtnbwwNr9mcCAunVrCzNb7lEw5GWDS9kiNFcARDMj+1M4EmY0yzMWYYWAOsGlVnFfCI/fpJ4FKxvoqsAtYYY4aMMbuBJru/mH3abS6x+8Du85oxxkBErgGagS3J37qiKKMJh0qPx+1WW+6hoy8L4pNCAtRYXLtiFlsP9bA5RraDsPjMmgDxaZxWQVvvUGRdZzTh/URLT6qmbyjAwHDmj6jIJMn81GYA+6Pet9hlMesYYwJAN1CboG288lqgy+5j9FgxxxCRcuA24OuJbkJEbhKRDSKyoa0t8QJjsaFpdZQwvRHLJ3W3W225h96hAEOBzB0q5w+GxhVsEM01p59EmcfJj1/be9y1fR3WWtCsyaUnNAbAQjvoYNuh2K63sOWzuN6q195b2EcwJCM+sX6yoz+t4tWZqPJEY3wdy02X8HAOY8wDxpgVxpgVdXV1iaoWL6pBRU/Y7VY5DrdbbYW1KJ9J19tQIIT7BNxuYN3rNe+bwdObDx631rL1UA9zastSOlI8HiOH0MXew9PWN4TH6WD+1Ar7fWGv+yTzU2sBZkW9nwkcjFdHRFxANXA0Qdt45e1Ajd3H6LHijXEW8E0R2QP8PfC/ReTWJO5LUZRRjLjdUv+wnVxupZjJpOvNHzQntOYT5uYL5hMKGe57semY8q0He1hSX3XC/YMlzjNqSmO698CydKZUeCIbWsfakJrvJPNTewNotKPQPFgBBGtH1VkL3Gi/vhZ4wVi+nLXAajtSbR7QCLwer0+7zYt2H9h9Pp1oDGPM+caYucaYucB3gW8YY/47hWdQ9KjBo4SJBByMY81nSoUtPhm0fPwTYPkAzK4t46MrZvLY6/vZZe/FOdztY0/HAKfNqjnh/sOcPquGt1viWz5TKr2RPUVFLz72+sqtwHPANuAJY8wWEblDRD5sV3sIa/2lCfgicLvddgvwBLAV+A1wizEmGK9Pu6/bgC/afdXafccdQ5k4VISUnkE/IlAxDjdT2PI52p+5D01/MITbNTE71b5w+UJKPU7+6aebGQ6EeH6rFX594cKJc9Mvm1nN/qODMV2Tbb1D1FV4qS334pDCF5+kfsOMMc8Cz44q+0rUax/w0Tht7wLuSqZPu7wZKxpudHncMaLqfC3RdUVREtNjp9aJdypnIsJrPpl0uw0HJ8byAZhaWcKd15zC3zy2kU/98HV2tfWxbGY1ixIcFJcqy2ZaVtTmli4uPnnqMdfa+4ZYNqMap0OYXO6lLcORg5lGMxwomlZHidDRP0xteeLjAeJRVeLC63JETv7MBP4JFB+AD512Ev92zSls2t9FMGS465pTx72BNRanzarG5RBe3330mPJgyNARlcanrtKrlo9SPKgIKUf7hyLus1QREeqrSzjQlblD0/xBQ6k7tQzcY/HJs+fwsRWzcAi4JlDYAMo8Lk6bVcOruzqOKe/oHyJkiKTxmVLhiez7KVTU8lEUJUJH3/C4xQegvrqUQxkVnxPf5xMLj8sx4cIT5pyGWt450E3f0Mgm0tGnpdZVemkvcMtHxUfRtDpKhKP9JyY+J9WUcqg7c/tThico2i2TnDO/lmDI8EaU6y1s5dTZRzeE3W6FvAE8v35qSlpRESpujDF0DgwzuXz8GZxPqinhSI8vctRBurGi3fLrY+yMOZPwOB280jRyjlDbaMunwstwMJTxk2EzSX791BRFSRs9vgD+oBl3wAFYbreQIWNBB8MnmFg0G5S4nZwzv5bntx6JWDaHun2IcEzAARR2loP8+qkpaaGALXslBcJ7T05ozafGchsdzNC6j88fomSCAw4ywZVLp7Pv6ADbD1t53vZ09HNSdWnkXsLik8nIwUyj4qNEUBEqbsKbQydXjF98Zk2yEnDu7xyYkDmNhc8fpMSdfx9jly+Zhgg8t8XayLq7vZ85tSOZs4shxU7+/dQURUkL4c2hJ+J2mz25HKdDaG6LfTLoROPzByc81DoT1FV6OWveZJ58s4XhQIidR3pZYCcUBairsCxIFR+loFGLRwE40hte9C4Zdx8el4NZk0oj+dHSSSAYSss+n0xxwzlzaekc5Lu/e4/+4SBnzJkUuVZV6sLjdNBewFkOVHyUCKpBxc2Rbh9Oh0TWG8ZLQ11FRiwfX8CKqCv15Kf4XLFkGg115dy3bhcel4MLGkdyyImItdFULR9FUQqdwz0+6iq8OMeR1y2ahinl7G7vJxRK79eZwWHr0Dpvnlo+LqeD+z6xnItPruMbf34qk0a5O+sqvQWd5UDT6yiKAlhHCEyvHr/LLcz8qRUMBUK0dA4yu/bEj5+Oh89viU++ut0AFk2v4oefPi6PMmCJz4EuDbVWioBC3k2tjM3hHh/Tq05cfMKHr205GPvQtIlisADEJxGFnlxUxUdRFMBa85kIy2dRfSVup/D2gfSKT9jyycdQ62Soq/BytH+IYJrdl9miMH9qSkqoxaP0DQXoHQpMiPh4XU5Onl7JO3GOi54owms+hWz5hIyV8boQUfFRIqgEFS/hjAT1EyA+AKfOsI6LTucXm4Gw+ORptNtYhJOMHulW8VEUpUDZ22FlJJhTWz4h/S2fXUOPLxBJH5MOenx+ACpL3GkbI5vMtLNFtGQoW0SmUfFR1OJR2Nth7cuZO0HRaec1TgHgj1GZmyeaXp+V8bmqpDCDdsORgvuOqvgohY6qUNGyt2OAqhIXNWXjT60TTX11KfPryo85NmCiKXTLp6rETXWpO2N58jKNio+iKOw9OjBhLrcw5zfW8VpzBwPD6TmTptcXwO2Ugo12A5g9uYx9RzN3MmwmKdyfmpI0Guym7O04NqvyRHDVKdPx+UM8v/XIhPYbptfnp7LEjcjEH6OdK8yaXEqLut2UQkdPMi1OhgJBWjoHmTvBls/7506mvrqEtZsOTmi/YXoGA1QW6HpPmFmTy2jpHCzIvT4qPopS5Oxq7ScYMpw8vXJC+3U4hFWnz2Dde21pidjqHvRTVaDrPWFmTSpjOBjiSE/hpdlR8VHU4ilydhzpAWDRBIsPwA3nzEGAH/5xz4T33d43dMIZuHOd2ZMLN+JNxUeJoGs/xcn2w714nA7mTplYtxvASTWlfHBZPT9Zv4/D3RP77b2td4i6isIWn3n2zyQT5yNlGhUfRUWnyNlxuJf5UytwO9PzcfAPV5xM0Bj+/dfbJqzPYMjQ0T9c8JbPjJpSSt1OdrVm5mTYTJLUb5uIrBSRHSLSJCK3x7juFZHH7evrRWRu1LUv2eU7ROTKsfoUkXl2HzvtPj2JxhCRM0Vkk/1vs4j8+XgfhqIUG8YYth7sYThRccsAACAASURBVHEaXG5hZk0u4+YL5/P0poP8cvPEBB90DgwTDJmCFx+HQ2ioKy9Oy0dEnMC9wFXAEuA6EVkyqtpngU5jzALgHuBuu+0SYDWwFFgJ3CcizjH6vBu4xxjTCHTafccdA3gXWGGMOd0e4/siUtghMGlCLaDi42C3j9beIU6fXZPWcf7mkgWcMWcStz/1Npv3d51wfwc6rb0v0ybgCIhcZ35dBU2tRSg+wJlAkzGm2RgzDKwBVo2qswp4xH79JHCpWMH3q4A1xpghY8xuoMnuL2afdptL7D6w+7wm0RjGmAFjTHgXWwm6Tz9l9IEVL2/t7QTgfbMmpXUct9PBvR9fzuQKD598aD0b9hw9of722OmAGuomfp0q11gwtYIDXYORLN6FQjLiMwPYH/W+xS6LWccWgm6gNkHbeOW1QFeUmESPFW8MROQsEdkCvAPcHNU+gojcJCIbRGRDW1tbEretKIXPW/s6KXE7WFSfPrdbmOnVJfzkc2czudzDdT94jQdfbiYQDI2rr93t/YiMRIMVMgumVgCFF3SQjPjE2j48+styvDoTVZ5wHsaY9caYpcD7gS+JyHG2uDHmAWPMCmPMirq6uhhdKWoBFR9v7e1k2YyatAUbjGbW5DKevvU8Llw4lTuf2cY19/2RP7zXlvLRCzsO9zJrUhklBXqWTzTz64pXfFqAWVHvZwKjVw0jdez1lmrgaIK28crbgZqoNZvoseKNEcEYsw3oB05J4r4UG13rKU66B/y8c6Cbs+fXZnTc6lI3P7jhDO79+HI6+oa54eHXuea+P/H0pgMMBcZ2LRljeHNvJ+9L8zpVrjB3ShkOgV0Ftu6TjPi8ATTaUWgerACCtaPqrAVutF9fC7xgrK8ya4HVdqTaPKAReD1en3abF+0+sPt8OtEYdh8uABGZA5wM7En6CShKkfJKUzshAxfYxx9kEhHh6mX1rPuni/j3vziVroFh/m7NJj7w7y/wf369PXLEQyx2HOmltXeI98+dnMEZZw+vy8nc2nLeO1JY4jNmVJgxJiAitwLPAU7gYWPMFhG5A9hgjFkLPAQ8KiJNWNbIarvtFhF5AtgKBIBbjDFBgFh92kPeBqwRkTuBjXbfxBsDOA+4XUT8QAj4a2NM+vK4FzB6nHZx8fLONipLXJw+K3sWhNfl5LozZ/OxFbN4pamdH6/fyw9ebub+l3Zx5rzJfGT5DP7s1Ppjjk1Y8/p+nA5h5SnTszbvTLO4vop3D6b3WPJMk1RIsjHmWeDZUWVfiXrtAz4ap+1dwF3J9GmXN2NFw40ujzmGMeZR4NExb0JJgIpOsREKGdbtaOPc+VNwZWi9JxEOh3DBwjouWFjH4W4fT73VwlNvtnDbU+/wlae3cOXS6VywsI4jPT4efW0vf7liJlMKPLtBNIumV/LMO4foHwpQ7i2MnSSFcReKoqTEm/s6Odzj46pTc896mF5dwi0XL+CvL5rPxv1d/OytFn65+RBr7Q2q5zdO4V+uHr3VsLBZVF8FWC7H5bPTGxafKVR8lAhq/xQPv9p8EK/LwaWLp2V7KnEREZbPnsTy2ZP46oeWsv/oAG6ng5mTSgv6DJ9YhJO+bj+k4qMUELrUU1z4gyGeeecwlyyaSkWeuHDcTgcNdshxMTJzUimVXhfbDvVkeyoTRvadvYqiZJTntx6hvW+Ia8+Yme2pKEkiIiyqr2T7YRUfpQBRC6g4+NGre5g5qZSLTp6a7akoKbBoehXbD/UWTFSqio+iaz1FxLZDPbzWfJRPnDUHp6O41k3ynUX1lfQOBWixk6rmOyo+ilJE/Ofvd1LpdXHdmbPGrqzkFIumWxFv2w/3ZnkmE4OKjxKF2kCFzLZDPfz63cN8+ty51JR5sj0dJUVOjkS8Fca6j4qPUjA+ZCU+xhi+8ew2Kr0uPnPevGxPRxkHFV4Xc2rL2FYgQQcqPoraO0XAc1sO8/LOdr54xUK1evKYRdMr2XZI3W5KgaEGUGHS4/Nzxy+3smh6JZ88e062p6OcAIvrq9jT0c/A8HFHluUdKj6KUuB89ektHOkd4ht/cWpO5HFTxs/i+iqMKYygA/1NVNTiKWCe3nSAn288wN9e0lgwaVmKmSV2jrdCyHSg4qNEUA0qLN490M3tT73DijmTuOXi+dmejjIBhNPsbC+AdR8VH0UpQFp7fXz+RxuoKXNz3/XL1d1WIITT7KjloxQEIfW7FRTdg34+/cM36Brw84MbVjC1siTbU1ImkMX1VWw/3EsolN9/tyo+SsTfphqU//QNBfjUD1/nvSO93PeJ5ZwyozrbU1ImmEXTq+grgDQ7Kj6KrvUUCD0+P5/54Ru83dLNf123nIsXaeLQQmRxvZXpYGueu95UfBR1uxUArb0+Pvb913hrXyff/djprDwl904oVSaGk6dXIpL/EW/5cZKUklbC4mPUBspL9rT3c8PDr9PeN8RDn3o/Fy6sy/aUlDRS5nExr7ZcxUfJf9TwyV9e2dnOLT95C4fATz5/NqfPqsn2lJQMsLi+ircPdGV7GieEut0UFZ88xBjDgy83c8PD65leVcIvbjlXhaeIWFxfyf6jg/T6/NmeyrhRy0cZcbupCOUF/UMB/uXn7/CLTQdZuXQ63/7L0yj36p9yMbHYznSw43AvK+ZOzvJsxof+xiq60pNHvHugm795bCN7Ovr54uULufXiBTj0RNKiY3FUmh0VHyVv0fN8ch9jDA+9spu7f7Od2nIvP/nc2Zwzvzbb01KyRH11CVUlLrbmcZodFR+F8EZplaDcpLXHx21Pvc2LO9q4fMk0vvmRZUwq1zN5ihkRYXF9VV5HvKn4KEqOYozhF5sO8LW1W/H5g9yxaimfPHsOIupmUyzX2+Nv7CcYMjjz0PWaVLSbiKwUkR0i0iQit8e47hWRx+3r60VkbtS1L9nlO0TkyrH6FJF5dh877T49icYQkctF5E0Recf+/5LxPgxFyRWO9FiJQb/w+GYap1bw6787nxvOmavCo0RYUl/FoD/I3o7+bE9lXIwpPiLiBO4FrgKWANeJyJJR1T4LdBpjFgD3AHfbbZcAq4GlwErgPhFxjtHn3cA9xphGoNPuO+4YQDvwIWPMqcCNwKOpPQIljK79ZB9jDE+92cLl33mJl3e2869XL+bx/3UODXUV2Z6akmOMBB3k57pPMpbPmUCTMabZGDMMrAFWjaqzCnjEfv0kcKlYX9FWAWuMMUPGmN1Ak91fzD7tNpfYfWD3eU2iMYwxG40xB+3yLUCJiHiTfQCKkiuEMxX8w083s3BaJb/5+wv43PkNeelSUdJP47QKnA7J23WfZNZ8ZgD7o963AGfFq2OMCYhIN1Brl782qu0M+3WsPmuBLmNMIEb9eGO0R/XzEWCjMWZo9E2IyE3ATQCzZ89OfMeKkkGGAkG+/1Iz//1iEx6ng699aAmfPGeuio6SkBK3k4Yp+ZtmJxnxifUXMNo/E69OvPJYFlei+mPOQ0SWYrnirohRD2PMA8ADACtWrFD/kpIT/KmpnX99+l2a2/r54LJ6vvzBJUyr0vN3lORYXF/Fm3s7sz2NcZGM+LQAs6LezwQOxqnTIiIuoBo4OkbbWOXtQI2IuGzrJ7p+vDEQkZnAz4EbjDG7krgnRckqbb1DfOPZbfx84wFmTy7jkc+cqQlBlZRZVF/J2s0H6R7wU13mzvZ0UiKZNZ83gEY7Cs2DFUCwdlSdtViL/QDXAi8Ya/V6LbDajlSbBzQCr8fr027zot0Hdp9PJxpDRGqAZ4AvGWP+mMrNK0qm8QdDPPzKbi799jp+9fZB/vaSBfz2Cxeo8CjjIhJ0cDj/XG9jWj72+sqtwHOAE3jYGLNFRO4ANhhj1gIPAY+KSBOWNbLabrtFRJ4AtgIB4BZjTBAgVp/2kLcBa0TkTmCj3TfxxgBuBRYAXxaRL9tlVxhjWsf3SIoXDXZLL394r407frWVptY+zm+cwlc/tJQFUzWKTRk/S6LS7JzdkF8ZL5LaZGqMeRZ4dlTZV6Je+4CPxml7F3BXMn3a5c1Y0XCjy2OOYYy5E7hzzJtQlCyxp72fO5/Zxu+2HWFObRk/uGEFly2eqnt2lBNmaqWXyeWevAw60AwHipIm+oYC/PcLTTz8ym7cTuG2lYv4zHlz8bqc2Z6aUiBYaXYq2X44//b6qPgoEfQk04khFDL8bOMB7v7Ndtp6h/jI8pnctvJkpmoUm5IGFk6rZM3r+wmFTF5lOFfxUZQJ5I9N7Xzj2W1sOdjDabNqeOCTZ/C+2ZOyPS2lgGmcWsmgP8iBrkFmTS7L9nSSRsVHUSaA7Yd7+Pdnt/PSe23MqCnlux87nQ+fdlJefRNV8pPGaVbQSlNrn4qPkp9otFvqHO728Z3nd/Dkmy1UeF38y58t5pPnzKHEres6SmZotCMm3zvSy8WLpmZ5Nsmj4qMo46DX5+f7LzXz4CvNhELw2fPmccvFC6gp03N2lMxSU+ahrtLLzta+bE8lJVR8FCUF/MEQj72+j//43U46+of58Gkn8U9XnpxX7g6l8Fg4rULFR8k8m/d3sereP/LyP1+c8odg9DEK6naLTzBk+OXmg9zzu/fY2zHAWfMm88OrF7NsZk22p6YoNE6t5Kcb9mOMyZv9Yyo+BcDjG6wE4evea+OTZ89JqW1IBSchxhie33qEb//2PXYc6WVxfRUP3biCSxbpJlEld1gwtYL+4SAHu33MqCnN9nSSQsWnAAjZCuIaR2RVIBSa6OkUDH9sauebz+1g8/4u5k0p57+uex9Xn1qvEWxKzrFwWiUAO4/0qvgomSNgi894zn8JBKPcbhM2o/zmrX2dfOu5HfxpVwcnVZdw90dO5SPLZ+JyJnXqvKJknHDE284jfVx0cn5EvKn4FACBoGW9jM/yUckJs+1QD9/+7Xv8btsRass9fPVDS7juzNkaNq3kPJPKPUyp8LKzNX/S7Kj4FAAnYvkEVXxoau3jP3+/k1++fZAKr4t/vGIhnz53HuVe/fNQ8ofGqfkV8aZ/XQVAMLLmk7pbKGw1wbGRb8VAU2sf//XCTtZuPkiJy8nNF87nf13QoHt1lLykcVoFP3/rQN5EvKn4FAAjls/42xYTo0XnpgsauOn8BmorvNmemqKMm8ZplfQOBTjc46O+OveDDlR8CoCw5TMeHSkmt5uKjlLIRAcdqPgoGSFsvYxHSKItn0KVIRUdpRiIiE9rHxfkwbHsKj4FQNDeqxMax5pN9JpPoaGioxQTk8s9VJe6aW7Lj6ADFZ8CIHgCls9wAYrPOy3d3P/SLp5995CKjlI0iAgNdeU0t/VneypJoeJTAJyI+Pj8UeKTx343Ywyv7urgey/t4uWd7VR6Xdx84Xw+d948FR2laGiYUsHLO9uyPY2kUPEpAAKRgIPU1WPIH5zo6WSUUMjw262H+d66XWxu6WZKhZfbr1rEx8+aTVWJO9vTU5SM0lBXzlNvtdDr81OZ47//Kj4FwIjlk3pbXyA/xWc4EOIXmw5w/0u7aG7rZ/bkMu7681P4yPKZmpFAKVrm15UDsLu9P+czrqv4FADh/GzjsXyi3W4mD/xu/UMBHnt9Hw++vJvDPT6W1FfxX9e9j6tOma6515Sip6HOinhrblPxUTJA8ATcbr48cbsd7R/mkT/t4ZFX99A14OeseZO5+9plXNA4JS92cytKJphTW4ZDyIuINxWfAsBvh1qfcMBBDnKga5AHX25mzev7GfQHuXzJNG6+cD5nzJmU7akpSs7hdTmZOamMXe25H/Gm4lMA+IYt62V84jNi+eRSaredR3q5/6Vmnt50AIBVp8/g5gsbaLTPLVEUJTb5Em6t4lMA+ALj32TaPxSY6OmcEBv3dfK9dbv47dYjlLgdXH/2HD53/jxmTkrteHBFKVYaplTwWnMHoZDJ6YMPk1qhFZGVIrJDRJpE5PYY170i8rh9fb2IzI269iW7fIeIXDlWnyIyz+5jp92nJ9EYIlIrIi+KSJ+I/Pd4H0Q+MxixfFJv2zXon+DZpI4xhpfea2P1A6/y5/f9ifW7j/K3lzbyp9sv5WsfXqrCoygpMH9qOT5/iEM9vmxPJSFjWj4i4gTuBS4HWoA3RGStMWZrVLXPAp3GmAUishq4G/iYiCwBVgNLgZOA34nIQrtNvD7vBu4xxqwRkfvtvr8XbwzAB3wZOMX+V3QM2q6z8Vg+3YN+HGIlJc201y0YMvz63UN8b90uthzsYXpVCf969WKuO3O2nqWjKOOkYUo44q0vp4/UTsbyORNoMsY0G2OGgTXAqlF1VgGP2K+fBC4VKwRpFbDGGDNkjNkNNNn9xezTbnOJ3Qd2n9ckGsMY02+MeQVLhIqO4cCIuRMax5pP96Cf6tLMbkbz+YP8ZP0+Lvn2Om79yUYGh4N88yPLeOmfL+Jz5zeo8CjKCRDe65Pr6z7J/JXPAPZHvW8BzopXxxgTEJFuoNYuf21U2xn261h91gJdxphAjPrxxmhP4h4Klq6B4cjr4HgsnwFLfDoH0u9+6/X5+fH6fTz0ym7aeodYNrOa+69fzuVLpo/rFFZFUY6nrtJLhdeV8+HWyYhPrE+F0Z9y8erEK49lcSWqn+w84iIiNwE3AcyePTvZZjlPa+9Q5PV4LB/r4KkS9nQMpC3arb1viB/+cTc/enUvvb4A5y2Ywn987HTOmV+re3QUZYKJJBjN8XDrZMSnBZgV9X4mcDBOnRYRcQHVwNEx2sYqbwdqRMRlWz/R9eONkRTGmAeABwBWrFiRQ0HFJ0Zb34j4pGr5hEKGQ92DrEjTnpn9Rwd44A/NPLFhP8PBEFedMp2bL5yf8zuvFSXfaZhSzht7OrM9jYQkIz5vAI0iMg84gBVA8PFRddYCNwKvAtcCLxhjjIisBX4iIt/BCjhoBF7HsmKO69Nu86Ldxxq7z6cTjTG+2y4cDnWNLHWlGu3W3jeEP2g4aYIXJXce6eV763bx9OaDOAQ+snwmN13QEEn9oShKemmoq+AXmw4yOByk1JObuQ7HFB97feVW4DnACTxsjNkiIncAG4wxa4GHgEdFpAnLGlltt90iIk8AW4EAcIsxJggQq097yNuANSJyJ7DR7pt4Y9h97QGqAI+IXANcMSoar2DZ2dpLqdtJyJiUo93eO2L5hBvsBcrxRMtFs3l/F/eta+K5LUcodTv59Afm8rnzG5heXXJC/SqKkhoNUQlGl5xUleXZxCapsCJjzLPAs6PKvhL12gd8NE7bu4C7kunTLm/GioYbXZ5ojLkJb6CA2Xaoh8ZpFexq7Us5w8GWg90AETfYeMTHGMNrzUe5b10TL+9sp6rExd9e2sinPjCXyeWelPtTFOXECYdb72rry2/xUXKTgeEAb+3t4sYPzGF3e3/K4rNhbyezJpcypcISiXB27GR5dVcH3/7tDjbs7WRKhZcv2efo5Po5IopS6Mybkvvh1io+ecxzWw4zHAxx8clTeWJDC6ksgfn8QV7Z2c61Z8yMhDkna/m8ta+T7/z2PV5pamdalZc7Vi3lL1fM0nN0FCVHKPU4mVFTSnN77oZbq/jkKYFgiO+/1EzDlHLObqjF6ZCUot3Wbj7IoD/IVaeO7LEJjGE57e3o585ntvH81iPUlnv48geX8ImzZqvoKEoOkusJRlV88pQf/nEP2w/3ct8nluNwCA6RpKPdhgJBvrduF4umV3JOQy1DgcRHMgwMB7j3xSZ+8IfduJ3CP16xkE+fO08zEShKDtMwpZwn37Q8Irm4n04/PfKQ9c0d3P2b7Vy+ZBpXnTIdAKcj+U2m3/3dTna39/Ojz5yJiOAKu91itH9zbyf/8MQm9nQM8Bfvm8FtVy1iWpVGrylKrtNQV0H/cJDW3qGc/JtV8ckzdh7p5a9+/Baza8v49l+eFvlG45Tk3G5PvdnC99btYvX7Z3HBwjqrbQy3mzGGh/+4h7ue2Up9dSk/+fxZfGD+lDTckaIo6SAcbr2rrU/FRzkxmlp7ue4H63E6hIdufD9VUVFlDoeMafk8+tpevvL0u5zTUMsdq0YSgIuIndnaam+M4a5ntvHgK7u5cuk0/u9HTztmLEVRcp/wpu7mtv6c/OKo4pMnvLqrg7/68Zu4HA4e+/zZkVDKMIkCDroH/Xxt7RZ+vvEAly2eyn9/fDkel+O49mHL53sv7eLBV3Zz4zlz+OqHlub0gVSKosSmvqqEErcjZ4MOVHxyHGMMj762lzt+uZU5tWU8/Kn3M6e2/Lh6VsDB8eLzh/fa+Ocn36atb4i/u7SRWy9ZgNt5fF5Xh1iW0+b9XXzruR18cFm9Co+i5DEOh9AwpSJnw61VfHKY1l4f//zk26zb0cbFJ9fxH9e9L677yyEck5W6fyjAv/96G//z2j4WTK3ggRvOSJjQ02VbPv/2q63UVXq5689PVeFRlDynoa6czS1d2Z5GTFR8cpBAMMSP1+/j27/dwVAgxB2rlvLJs+ckDJd0OkYsn7f2dfKFxzex7+gAnztvHv945clj7sVxOITXdx/lnQPd3LFqacYPmFMUZeJpqKvgmXcO4fMHc24/nopPDmGM4ZWmdu56ZhvbD/dy3oIpfO3DS1kwdexs0A472u3pTQf4hyc2M726hDWfP5uzGmqTGtvlEN450I3X5eDaM2ae6K0oipIDzK8rxxjY2zHAydMrsz2dY1DxyQGMMfyxqYPv/u49NuztZOakUu6//gyuXDot6c1hTofw3pFe/vGnrayYO4nvf3JFStZLONz6nPm1lHn010JRCoFwgtHmtj4VH2WEUMjwwvZW7n9pFxv2dlJfXcK/XXMKf7liJl5Xaiay0yHs7Rigwuvi3o8vT9lt1t5nHcd9yaKpKbVTFCV3mWfv9cnFU01VfLJA31CAn27Yz//70x72dgxw0gmITpjwes/Vp9ZTW+Ed99wuPlnFR1EKhQqvi2lVXna15V7Em4pPBtnXMcD/+9MefrphP71DAc6YM4l/uvJkrlw6PWb4cypsOdgDwOVLpo2r/Zc/uIQtB7qZNbnshOahKEpu0TClIif3+qj4pJlQyPDSzjYefXUvL+5oxSnC1cvq+fS58zh9VvzQ51RZdfpJPL3pIOcuGN9O5s+eN2/C5qIoSu7QUFfOLzcfzLkEoyo+aaJrYJifbmjhf9bvZW/HAFMqvPzNxQv4xNlz0pJn6Y5Vp/D3ly3M2fPaFUXJDg11FfT4AnT0DzPlBFzyE42KzwSzu72fB/6wi5+9dYChQIj3z53EP1xxMiuXTj8upc1EUl3q1r05iqIcRzjBaHNbv4pPIdLU2ss9z+/k2XcP4XY6+MjyGdxwzlwW1+fm+emKohQHC+wEo7va+jhz3uQsz2YEFZ8TpG8owLee28Gjr+2l1O3k5gvn85lz51FXmTvfMBRFKV5OqinF43LQnGMRbyo+J8DWgz3c8pO32NvRz3VnzuaLly88oTBnRVGUicbpEObV5t6R2io+42Tjvk5uePh1yjxOHkshjY2iKEqmaagrZ/vh3mxP4xjStwJewLT2+vj8j95kUpmHn/31uSo8iqLkNA115ew7OsBwIJTtqURQ8RkHX1+7lb4hPw/euIIZNaXZno6iKEpCFkytIBgy7M6hNDsqPimy43Avz7xziJvOb2DhtNxK1KcoihKLpSdVA/Duge4sz2QEFZ8U+dlbLbgcwqfP1YwAiqLkB/PrKih1O3lHxSd/eem9Ns6ZX8ukck+2p6IoipIUToew9KSq/LN8RGSliOwQkSYRuT3Gda+IPG5fXy8ic6Oufcku3yEiV47Vp4jMs/vYaffpGe8YE43PH2Rna9+E5mRTFEXJBKfMqGbLwZ5IBvxsM6b4iIgTuBe4ClgCXCciS0ZV+yzQaYxZANwD3G23XQKsBpYCK4H7RMQ5Rp93A/cYYxqBTrvvlMdI9UEkQ1vvEMGQ0czPiqLkHSvmTmLQH+StfZ3ZngqQnOVzJtBkjGk2xgwDa4BVo+qsAh6xXz8JXCpW+tRVwBpjzJAxZjfQZPcXs0+7zSV2H9h9XjPOMSac7kE/gOZQUxQl77hwYR0ep4NHX92b7akAyW0ynQHsj3rfApwVr44xJiAi3UCtXf7aqLYz7Nex+qwFuowxgRj1xzNGBBG5CbjJftsnIh1Ae9y7TsDKu8fTKqeZwjifRQGiz8JCn8MIBfUsdgL/9fFxNZ0CzJmoeSQjPrEOgBjtNIxXJ155LIsrUf3xjHFsgTEPAA+E34vIBmPMihhtiw59FiPos7DQ5zCCPgsL+znMnaj+knG7tQCzot7PBA7GqyMiLqAaOJqgbbzydqDG7mP0WKmOoSiKouQoyYjPG0CjHYXmwVrcXzuqzlrgRvv1tcALxhhjl6+2I9XmAY3A6/H6tNu8aPeB3efT4xxDURRFyVHGdLvZ6yu3As8BTuBhY8wWEbkD2GCMWQs8BDwqIk1Y1shqu+0WEXkC2AoEgFuMMUGAWH3aQ94GrBGRO4GNdt+MZ4wxeGDsKkWDPosR9FlY6HMYQZ+FxYQ+B7GMB0VRFEXJHJrhQFEURck4Kj6KoihKxilK8RkrXVAhICIPi0iriLwbVTZZRJ63Uxc9LyKT7HIRkf+0n8fbIrI8qs2Ndv2dInJjrLFyGRGZJSIvisg2EdkiIn9nlxfVsxCREhF5XUQ228/h63Z5zqazSjd2tpWNIvIr+31RPgsR2SMi74jIJhHZYJel/+/DGFNU/7ACHHYBDYAH2Awsyfa80nCfFwDLgXejyr4J3G6/vh242379Z8CvsfZMnQ2st8snA832/5Ps15OyfW8pPod6YLn9uhJ4DyulU1E9C/t+KuzXbmC9fX9PAKvt8vuBv7Jf/zVwv/16NfC4/XqJ/TfjBebZf0vObN/fOJ/JF4GfAL+y3xflswD2AFNGlaX976MYLZ9k0gXlPcaYP2BFBUYTnaJodOqiHxmL17D2WtUDVwLPG2OOGmM6geex8uflDcaYQ8aYt+zXvcA2rAwYRfUs7Pvps9+67X+GAOvoTgAAAk1JREFUHE5nlU5EZCZwNfCg/T6nU3tlgbT/fRSj+MRKF3RcOp4CZZox5hBYH8rAVLs83jMpqGdlu0veh/Wtv+iehe1m2gS0Yn047CLJdFZAdDqrvH4ONt8F/hkInyuddGovCu9ZGOC3IvKmWGnIIAN/H8mk1yk0kkrHU2ScUOqifEBEKoCngL83xvRYX1xjV41RVhDPwlj7304XkRrg58DiWNXs/wv2OYjIB4FWY8ybInJRuDhG1YJ/FjbnGmMOishU4HkR2Z6g7oQ9i2K0fIo5Hc8R20TG/r/VLk81DVJeISJuLOH5sTHmZ3ZxUT4LAGNMF7AOy2dfjOmszgU+LCJ7sNzul2BZQsX4LDDGHLT/b8X6UnImGfj7KEbxSSZdUKESnaJodOqiG+xIlrOBbtvUfg64QkQm2dEuV9hleYPtm38I2GaM+U7UpaJ6FiJSZ1s8iEgpcBnW+lfRpbMyxnzJGDPTWEkyV2Pd2ycowmchIuUiUhl+jfV7/S6Z+PvIdqRFNv5hRWy8h+Xz/pdszydN9/gYcAjwY30r+SyWn/r3WFnVfw9MtusK1uF+u4B3gBVR/XwGayG1Cfh0tu9rHM/hPCzz/21gk/3vz4rtWQDLsNJVvW1/uHzFLm/A+sBsAn4KeO3yEvt9k329Iaqvf7Gfzw7gqmzf2wk+l4sYiXYrumdh3/Nm+9+W8OdhJv4+NL2OoiiKknGK0e2mKIqiZBkVH0VRFCXjqPgoiqIoGUfFR1EURck4Kj6KoihKxlHxURRFUTKOio+iKIqScf4/2QYEJG6Qy5oAAAAASUVORK5CYII=\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], + "outputs": [], "source": [ "plt.clf()\n", "# plt.plot(x_part, calcs, '.')\n", @@ -743,7 +699,7 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -759,7 +715,7 @@ }, { "cell_type": "code", - "execution_count": 15, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -775,7 +731,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -788,7 +744,7 @@ }, { "cell_type": "code", - "execution_count": 17, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -845,7 +801,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -929,7 +885,7 @@ }, { "cell_type": "code", - "execution_count": 19, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -938,7 +894,7 @@ }, { "cell_type": "code", - "execution_count": 20, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -947,7 +903,7 @@ }, { "cell_type": "code", - "execution_count": 21, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -956,7 +912,7 @@ }, { "cell_type": "code", - "execution_count": 22, + "execution_count": null, "metadata": { "scrolled": false }, @@ -1004,7 +960,7 @@ }, { "cell_type": "code", - "execution_count": 23, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1021,7 +977,7 @@ }, { "cell_type": "code", - "execution_count": 24, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1045,7 +1001,7 @@ }, { "cell_type": "code", - "execution_count": 25, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1072,7 +1028,7 @@ }, { "cell_type": "code", - "execution_count": 26, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1095,7 +1051,7 @@ }, { "cell_type": "code", - "execution_count": 27, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1104,7 +1060,7 @@ }, { "cell_type": "code", - "execution_count": 28, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1120,7 +1076,7 @@ }, { "cell_type": "code", - "execution_count": 29, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1150,7 +1106,7 @@ }, { "cell_type": "code", - "execution_count": 30, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1164,7 +1120,7 @@ }, { "cell_type": "code", - "execution_count": 31, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1182,7 +1138,7 @@ }, { "cell_type": "code", - "execution_count": 32, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1196,7 +1152,7 @@ }, { "cell_type": "code", - "execution_count": 33, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1217,7 +1173,7 @@ }, { "cell_type": "code", - "execution_count": 34, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1227,7 +1183,7 @@ }, { "cell_type": "code", - "execution_count": 35, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1258,7 +1214,7 @@ }, { "cell_type": "code", - "execution_count": 36, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1275,7 +1231,7 @@ }, { "cell_type": "code", - "execution_count": 43, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -1289,7 +1245,7 @@ "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(DDstar_m,2))))\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", @@ -1308,7 +1264,64 @@ "\n", "# 3. Constraint - Maximum eta of D contribs\n", "\n", - "constraint3_0 = tf.cond(tf.greater_equal(sum_2, 5.0e-8), lambda: 10000., lambda: 0.)" + "constraint3_0 = tf.cond(tf.greater_equal(tf.abs(Dbar_s), 0.2), lambda: 10000., lambda: 0.)\n", + "\n", + "constraint3_1 = tf.cond(tf.greater_equal(tf.abs(DDstar_s), 0.2), lambda: 10000., 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)\n", + " sigma2 = ztf.constant(0.128)\n", + " sigma3 = ztf.constant(0.548)\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", + "# 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), 0.02), lambda: 10000., lambda: 0.)\n", + "\n", + "#List of all constraints\n", + "\n", + "constraints = [constraint1, constraint2, constraint3_0, constraint3_1, constraint4]" ] }, { @@ -1320,46 +1333,18 @@ }, { "cell_type": "code", - "execution_count": 41, + "execution_count": null, "metadata": { "scrolled": false }, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "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" - ] - }, - { - "ename": "RuntimeError", - "evalue": "exception was raised in user function\nUser function arguments:\n p3770_s = +1.999020\n p4160_s = +2.376152\n p4040_p = -4.437832\n p4040_s = +1.423904\n Dbar_s = +0.342445\n p4160_p = +4.897009\n p4415_p = +4.364139\n p4415_s = +1.177263\n Ctt = -0.212905\n bplus_0 = -0.141110\n psi2s_p = -2.316177\n p3770_p = +0.511010\n bplus_2 = -1.674701\n jpsi_p = +4.878201\n bplus_1 = +1.591894\nOriginal python exception in user function:\nInternalError: Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\nCaused by op 'model_2/_call_unnormalized_pdf/Sqrt_27', 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 \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 \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 365, 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 272, 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 542, 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 \"\", line 71, in \n result = minimizer.minimize(nll)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 207, in minimize\n return self._hook_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 216, in _hook_minimize\n return self._call_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 220, in _call_minimize\n return self._minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 48, in _minimize\n gradients = loss.gradients(params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 158, in gradients\n return self._gradients(params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 251, in _gradients\n gradients = {(p, grad) for p, grad in zip(params_todo, super()._gradients(params_todo))}\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 220, in _gradients\n return tf.gradients(self.value(), params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 198, in value\n return self._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 234, in _value\n loss = super()._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 203, in _value\n constraints=self.constraints)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 266, in _loss_func\n nll = _unbinned_nll_tf(model=model, data=data, fit_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in _unbinned_nll_tf\n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in \n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 43, in _unbinned_nll_tf\n probs = model.pdf(data, norm_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 314, in pdf\n value = self._single_hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 324, in _single_hook_pdf\n return self._hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 327, in _hook_pdf\n return self._norm_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 330, in _norm_pdf\n return self._call_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 338, in _call_pdf\n return self._fallback_pdf(x=x, norm_range=norm_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 341, in _fallback_pdf\n pdf = self._call_unnormalized_pdf(x, name=\"_call_unnormalized_pdf\")\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 289, in _call_unnormalized_pdf\n return self._unnormalized_pdf(x)\n File \"\", line 73, in _unnormalized_pdf\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 File \"\", line 70, in ttau_cusp\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 File \"\", line 15, in h_S\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 File \"\", line 9, in G\n big_bracket = tf.where(tf.math.real(y) > ztf.constant(0.0), inner_rect_bracket(y), inner_right(y))\n File \"\", line 4, in inner_rect_bracket\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 File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_math_ops.py\", line 9251, in sqrt\n \"Sqrt\", x=x, 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\nInternalError (see above for traceback): Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 116, in grad_func\n gradients_values = self.sess.run(gradients)\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 1348, in _do_call\n raise type(e)(node_def, op, message)\n", - "output_type": "error", - "traceback": [ - "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", - "\u001b[1;31mRuntimeError\u001b[0m Traceback (most recent call last)", - "\u001b[1;32m\u001b[0m in \u001b[0;36m\u001b[1;34m\u001b[0m\n\u001b[0;32m 69\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 70\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---> 71\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 72\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 73\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Toy {}: Fitting finished\"\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mformat\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mtoy\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;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 p3770_s = +1.999020\n p4160_s = +2.376152\n p4040_p = -4.437832\n p4040_s = +1.423904\n Dbar_s = +0.342445\n p4160_p = +4.897009\n p4415_p = +4.364139\n p4415_s = +1.177263\n Ctt = -0.212905\n bplus_0 = -0.141110\n psi2s_p = -2.316177\n p3770_p = +0.511010\n bplus_2 = -1.674701\n jpsi_p = +4.878201\n bplus_1 = +1.591894\nOriginal python exception in user function:\nInternalError: Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\nCaused by op 'model_2/_call_unnormalized_pdf/Sqrt_27', 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 \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 \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 365, 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 272, 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 542, 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 \"\", line 71, in \n result = minimizer.minimize(nll)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 207, in minimize\n return self._hook_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 216, in _hook_minimize\n return self._call_minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\baseminimizer.py\", line 220, in _call_minimize\n return self._minimize(loss=loss, params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 48, in _minimize\n gradients = loss.gradients(params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 158, in gradients\n return self._gradients(params=params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 251, in _gradients\n gradients = {(p, grad) for p, grad in zip(params_todo, super()._gradients(params_todo))}\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 220, in _gradients\n return tf.gradients(self.value(), params)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 198, in value\n return self._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 234, in _value\n loss = super()._value()\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 203, in _value\n constraints=self.constraints)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 266, in _loss_func\n nll = _unbinned_nll_tf(model=model, data=data, fit_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in _unbinned_nll_tf\n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 39, in \n for p, d, r in zip(model, data, fit_range)]\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py\", line 43, in _unbinned_nll_tf\n probs = model.pdf(data, norm_range=fit_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 314, in pdf\n value = self._single_hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 324, in _single_hook_pdf\n return self._hook_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 327, in _hook_pdf\n return self._norm_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 330, in _norm_pdf\n return self._call_pdf(x=x, norm_range=norm_range, name=name)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 338, in _call_pdf\n return self._fallback_pdf(x=x, norm_range=norm_range)\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 341, in _fallback_pdf\n pdf = self._call_unnormalized_pdf(x, name=\"_call_unnormalized_pdf\")\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\basepdf.py\", line 289, in _call_unnormalized_pdf\n return self._unnormalized_pdf(x)\n File \"\", line 73, in _unnormalized_pdf\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 File \"\", line 70, in ttau_cusp\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 File \"\", line 15, in h_S\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 File \"\", line 9, in G\n big_bracket = tf.where(tf.math.real(y) > ztf.constant(0.0), inner_rect_bracket(y), inner_right(y))\n File \"\", line 4, in inner_rect_bracket\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 File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\tensorflow\\python\\ops\\gen_math_ops.py\", line 9251, in sqrt\n \"Sqrt\", x=x, 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\nInternalError (see above for traceback): Dst tensor is not initialized.\n\t [[node model_2/_call_unnormalized_pdf/Sqrt_27 (defined at :4) ]]\n\t [[node gradients_1/AddN_11 (defined at C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\core\\loss.py:220) ]]\n\n File \"C:\\Users\\sa_li\\.conda\\envs\\rmd\\lib\\site-packages\\zfit\\minimizers\\minimizer_minuit.py\", line 116, in grad_func\n gradients_values = self.sess.run(gradients)\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 1348, in _do_call\n raise type(e)(node_def, op, message)\n" - ] - } - ], + "outputs": [], "source": [ "# zfit.run.numeric_checks = False \n", "\n", "Ctt_list = []\n", "Ctt_error_list = []\n", "\n", - "nr_of_toys = 1\n", + "nr_of_toys = 2\n", "nevents = int(pdg[\"number_of_decays\"])\n", "nevents = pdg[\"number_of_decays\"]\n", "event_stack = 1000000\n", @@ -1420,7 +1405,7 @@ " 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 = constraint1)\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", diff --git a/raremodel-nb.py b/raremodel-nb.py index acbf931..fab7c0c 100644 --- a/raremodel-nb.py +++ b/raremodel-nb.py @@ -284,7 +284,7 @@ # ## Build pdf -# In[5]: +# In[11]: class total_pdf(zfit.pdf.ZPDF): @@ -302,7 +302,7 @@ 'omega_mass', 'omega_scale', 'omega_phase', 'omega_width', 'phi_mass', 'phi_scale', 'phi_phase', 'phi_width', 'Dbar_mass', 'Dbar_scale', 'Dbar_phase', - 'DDstar_mass', 'DDstar_scale', 'DDstar_phase', + 'Dstar_mass', 'DDstar_scale', 'DDstar_phase', 'D_mass', 'tau_mass', 'C_tt'] # the name of the parameters @@ -352,7 +352,7 @@ def P2_D(q): 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)) - 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['DDstar_mass'], q)) + 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))) return Dbar_contrib + DDstar_contrib def ttau_cusp(q): @@ -372,7 +372,7 @@ # ## Load data -# In[6]: +# In[12]: x_min = 2*pdg['muon_M'] @@ -392,7 +392,7 @@ # ## Setup parameters -# In[7]: +# In[13]: # formfactors @@ -495,38 +495,29 @@ # ## Dynamic generation of 2 particle contribution -# In[8]: +# In[14]: m_c = 1300 -DDstar_eta = 0 -Dbar_phase = np.random.uniform(0.0, 2*np.pi) -DDstar_phase = np.random.uniform(0.0, 2*np.pi) -DDstar_mass = pdg['D0_M'] - -if Dbar_phase < np.pi: - Dbar_phase = 0.0 -else: - Dbar_phase = np.pi - -Dbar_mass = (pdg['D0_M']+pdg['Dst_M'])/2 - +Dbar_phase = 0.0 +DDstar_phase = 0.0 +Dstar_mass = pdg['Dst_M'] +Dbar_mass = pdg['D0_M'] +D_mass = pdg['D0_M'] Dbar_s = zfit.Parameter("Dbar_s", ztf.constant(0.0), lower_limit=-1.464, upper_limit=1.464) Dbar_m = zfit.Parameter("Dbar_m", ztf.constant(Dbar_mass), floating = False) Dbar_p = zfit.Parameter("Dbar_p", ztf.constant(Dbar_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)#, floating = False) DDstar_s = zfit.Parameter("DDstar_s", ztf.constant(0.0), lower_limit=-2.0, upper_limit=2.0)#, floating = False) -DDstar_m = zfit.Parameter("DDstar_m", ztf.constant(DDstar_mass), floating = False) +Dstar_m = zfit.Parameter("Dstar_m", ztf.constant(Dstar_mass), floating = False) +D_m = zfit.Parameter("D_m", ztf.constant(D_mass), floating = False) DDstar_p = zfit.Parameter("DDstar_p", ztf.constant(DDstar_phase), lower_limit=-2*np.pi, upper_limit=2*np.pi)#, floating = False) -# Dbar_s.set_value(0.0) -# DDstar_s.set_value(0.0) - # ## Tau parameters -# In[9]: +# In[15]: tau_m = zfit.Parameter("tau_m", ztf.constant(pdg['tau_M']), floating = False) @@ -535,7 +526,7 @@ # ## Setup pdf -# In[10]: +# In[16]: total_f = total_pdf(obs=obs, jpsi_mass = jpsi_m, jpsi_scale = jpsi_s, jpsi_phase = jpsi_p, jpsi_width = jpsi_w, @@ -547,7 +538,7 @@ rho_mass = rho_m, rho_scale = rho_s, rho_phase = rho_p, rho_width = rho_w, omega_mass = omega_m, omega_scale = omega_s, omega_phase = omega_p, omega_width = omega_w, phi_mass = phi_m, phi_scale = phi_s, phi_phase = phi_p, phi_width = phi_w, - DDstar_mass = DDstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p, + Dstar_mass = Dstar_m, DDstar_scale = DDstar_s, DDstar_phase = DDstar_p, D_mass = D_m, Dbar_mass = Dbar_m, Dbar_scale = Dbar_s, Dbar_phase = Dbar_p, tau_mass = tau_m, C_tt = Ctt, b0_0 = b0_0, b0_1 = b0_1, b0_2 = b0_2, bplus_0 = bplus_0, bplus_1 = bplus_1, bplus_2 = bplus_2, @@ -564,7 +555,7 @@ # ## Test if graphs actually work and compute values -# In[11]: +# In[17]: # def total_test_tf(xq): @@ -607,7 +598,7 @@ # fT_y = zfit.run(tf.math.real(formfactor(test_q,"T", b0, bplus, bT))) -# In[12]: +# In[18]: plt.clf() @@ -625,7 +616,7 @@ # print(jpsi_width) -# In[13]: +# In[19]: @@ -638,7 +629,7 @@ # plt.semilogy(test_q, calcs_test, label = 'pdf') -# In[14]: +# In[20]: # 0.213/(0.00133+0.213+0.015) @@ -646,7 +637,7 @@ # ## Adjust scaling of different parts -# In[15]: +# In[21]: total_f.update_integration_options(draws_per_dim=200000, mc_sampler=None) @@ -656,7 +647,7 @@ # print(pdg["jpsi_BR"]/pdg["NR_BR"], inte_fl*pdg["psi2s_auc"]/pdg["NR_auc"]) -# In[16]: +# In[22]: # # print("jpsi:", inte_fl) @@ -705,7 +696,7 @@ # # Sampling # ## Mixture distribution for sampling -# In[17]: +# In[23]: @@ -786,25 +777,25 @@ return sample, thresholds, weights, weights_max, n_to_produce -# In[18]: +# In[24]: total_f._sample_and_weights = UniformSampleAndWeights -# In[19]: +# In[25]: # 0.00133/(0.00133+0.213+0.015)*(x_max-3750)/(x_max-x_min) -# In[20]: +# In[26]: # zfit.settings.set_verbosity(10) -# In[21]: +# In[27]: # # zfit.run.numeric_checks = False @@ -847,7 +838,7 @@ # pkl.dump(sam, f, pkl.HIGHEST_PROTOCOL) -# In[22]: +# In[28]: # with open(r"data/zfit_toys/toy_0/0.pkl", "rb") as input_file: @@ -861,7 +852,7 @@ # print(np.sum(sam-sam2)) -# In[23]: +# In[29]: # print("Time to generate full toy: {} s".format(int(time.time()-start))) @@ -882,7 +873,7 @@ # print(total_samp[:nevents].shape) -# In[24]: +# In[30]: # plt.clf() @@ -906,7 +897,7 @@ # plt.savefig('test2.png') -# In[25]: +# In[31]: # sampler = total_f.create_sampler(n=nevents) @@ -926,13 +917,13 @@ # minimum = minimizer.minimize(nll) -# In[26]: +# In[32]: # jpsi_width -# In[27]: +# In[33]: # plt.hist(sample, weights=1 / prob(sample)) @@ -940,7 +931,7 @@ # # Fitting -# In[28]: +# In[34]: # start = time.time() @@ -967,7 +958,7 @@ # print("Hesse errors:", result.hesse()) -# In[29]: +# In[35]: # print("Time taken for fitting: {}".format(display_time(int(time.time()-start)))) @@ -978,7 +969,7 @@ # res_y = zfit.run(jpsi_res(test_q)) -# In[30]: +# In[36]: # plt.clf() @@ -993,7 +984,7 @@ # # print(jpsi_width) -# In[31]: +# In[37]: # _tot = 4.37e-7+6.02e-5+4.97e-6 @@ -1004,7 +995,7 @@ # print(_probs) -# In[32]: +# In[38]: # dtype = 'float64' @@ -1022,14 +1013,14 @@ # # print(zfit.run(mixture.prob(mixture.sample((10, 1))))) -# In[33]: +# In[39]: # print((zfit.run(jpsi_p)%(2*np.pi))/np.pi) # print((zfit.run(psi2s_p)%(2*np.pi))/np.pi) -# In[34]: +# In[40]: # def jpsi_res(q): @@ -1057,7 +1048,7 @@ # phase = p4415_phase, width = p4415_width) -# In[35]: +# In[41]: # 0.15**2*4.2/1000 @@ -1066,12 +1057,11 @@ # ## Constraints -# In[36]: +# In[62]: # 1. Constraint - Real part of sum of Psi contrib and D contribs - sum_list = [] 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)))) @@ -1080,7 +1070,7 @@ 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)))) 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)))) 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)))) -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(DDstar_m,2)))) +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))))) 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)))) sum_ru_1 = ztf.to_complex(ztf.constant(0.0)) @@ -1089,8 +1079,74 @@ sum_ru_1 += part sum_1 = tf.math.real(sum_ru_1) -constraint1 = zfit.constraint.GaussianConstraint(params = [sum_1], mu = [ztf.constant(1.7*10**-8)], - sigma = [ztf.constant(2.2*10**-8)]) +constraint1 = zfit.constraint.GaussianConstraint(params = sum_1, mu = ztf.constant(1.7*10**-8), + sigma = ztf.constant(2.2*10**-8)) + +# 2. Constraint - Abs. of sum of Psi contribs and D contribs + +sum_2 = tf.abs(sum_ru_1) +constraint2 = tf.cond(tf.greater_equal(sum_2, 5.0e-8), lambda: 10000., lambda: 0.) + +# 3. Constraint - Maximum eta of D contribs + +constraint3_0 = tf.cond(tf.greater_equal(tf.abs(Dbar_s), 0.2), lambda: 10000., lambda: 0.) + +constraint3_1 = tf.cond(tf.greater_equal(tf.abs(DDstar_s), 0.2), lambda: 10000., lambda: 0.) + +# 4. Constraint - Formfactor multivariant gaussian covariance fplus + +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)], + [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)], + [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)], + [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)], + [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)], + [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)], + [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)], + [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)], + [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.)]] + +def triGauss(val1,val2,val3,m = Cov_matrix): + + mean1 = ztf.constant(0.466) + mean2 = ztf.constant(-0.885) + mean3 = ztf.constant(-0.213) + sigma1 = ztf.constant(0.014) + sigma2 = ztf.constant(0.128) + sigma3 = ztf.constant(0.548) + x1 = (val1-mean1)/sigma1 + x2 = (val2-mean2)/sigma2 + x3 = (val3-mean3)/sigma3 + rho12 = m[0][1] + rho13 = m[0][2] + rho23 = m[1][2] + 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)) + d = 2*(rho12*rho12+rho13*rho13+rho23*rho23-2*rho12*rho13*rho23-1) + + fcn = -w/d + chisq = -2*fcn + return chisq + +constraint4 = triGauss(bplus_0, bplus_1, bplus_2) + +# 5. Constraint - Abs. of sum of light contribs + +sum_list_5 = [] + +sum_list_5.append(rho_s*rho_w/rho_m) +sum_list_5.append(omega_s*omega_w/omega_m) +sum_list_5.append(phi_s*phi_w/phi_m) + + +sum_ru_5 = ztf.constant(0.0) + +for part in sum_list_5: + sum_ru_5 += part + +constraint5 = tf.cond(tf.greater_equal(tf.abs(sum_ru_5), 0.02), lambda: 10000., lambda: 0.) + +#List of all constraints + +constraints = [constraint1, constraint2, constraint3_0, constraint3_1, constraint4] # # Analysis @@ -1103,7 +1159,7 @@ Ctt_list = [] Ctt_error_list = [] -nr_of_toys = 1 +nr_of_toys = 2 nevents = int(pdg["number_of_decays"]) nevents = pdg["number_of_decays"] event_stack = 1000000 @@ -1164,7 +1220,7 @@ for param in total_f.get_dependents(): param.randomize() - nll = zfit.loss.UnbinnedNLL(model=total_f, data=data, fit_range = (x_min, x_max), constraints = constraint1) + nll = zfit.loss.UnbinnedNLL(model=total_f, data=data, fit_range = (x_min, x_max), constraints = constraints) minimizer = zfit.minimize.MinuitMinimizer(verbosity = 5) # minimizer._use_tfgrad = False