{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "

11: Metodo Monte Carlo

\n", "\n", "Introduzione\n", "------------\n", "Il Metodo Monte Carlo si riferisce a un insieme di tecniche di calcolo e simulazione che si basano sull'uso di numeri casuali." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 11.1 Integrazione numerica con metodo MC\n", "\n", "L'idea di usare numeri casuali per calcolare un integrale numerico si basa sul teorema della media in Analisi che afferma che l'integrale $\\int_a^b f(x)\\,dx$ è uguale al prodotto della lunghezza (con segno) dell'intervallo di integrazione $(b-a)$ per il valor medio della funzione $f$, $$, nell'intervallo. Il valor medio può essere calcolato\n", "facendo la media dei valori di $f$ in un set di punti scelti a caso in $[a,b]$.
\n", "Esempio: $\\int_1^2 x^2 dx = \\frac{7}{3}$
" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Stima dell'integrale: 2.520720064255817\n", "Risultato esatto: 2.3333333333333335\n" ] } ], "source": [ "import numpy as np\n", "from numpy.random import default_rng\n", "\n", "rng = default_rng()\n", "\n", "npoints=10\n", "xval=rng.uniform(size=npoints)+1. #punti casuali in [1,2]\n", "#print(xval)\n", "\n", "y_ave=np.mean(xval*xval)\n", "my_int=y_ave*(2.-1.)\n", "print(\"Stima dell'integrale:\",my_int)\n", "print(\"Risultato esatto:\",7/3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Aumentiamo il numero di punti:" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Stima dell'integrale: 2.2300646377312843\n", "Risultato esatto: 2.3333333333333335\n" ] } ], "source": [ "npoints=100\n", "xval=rng.uniform(size=npoints)+1. #punti casuali in [1,2]\n", "#print(xval)\n", "\n", "y_ave=np.mean(xval*xval)\n", "my_int=y_ave*(2.-1.)\n", "print(\"Stima dell'integrale:\",my_int)\n", "print(\"Risultato esatto:\",7/3)" ] }, { "cell_type": "code", "execution_count": 30, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Stima dell'integrale: 2.3142078477545596\n", "Risultato esatto: 2.3333333333333335\n" ] } ], "source": [ "npoints=1000\n", "xval=rng.uniform(size=npoints)+1. #punti casuali in [1,2]\n", "#print(xval)\n", "\n", "y_ave=np.mean(xval*xval)\n", "my_int=y_ave*(2.-1.)\n", "print(\"Stima dell'integrale:\",my_int)\n", "print(\"Risultato esatto:\",7/3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Confrontate il risultato con quello fornito dal metodo dei trapezi e da `quad` nel notebook 10_scipy." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Esempio in più dimensioni: il volume di una semisfera di raggio 1: \n", "$$\\int_{-1}^1 dx \\int_{-\\sqrt{1-x^2}}^{\\sqrt{1-x^2}}dy\\,\\sqrt{1-x^2-y^2} = \\frac{2}{3}\\pi$$" ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [], "source": [ "#npoints = 100\n", "npoints = 10000\n", "xval=2*rng.uniform(size=npoints)-1.\n", "#print(xval)\n", "yval=2*rng.uniform(size=npoints)-1.\n", "#print(yval)\n", "points=[pair for pair in zip(xval,yval) if pair[0]**2 + pair[1]**2 < 1.] #zip returns tuples\n", "points=np.array(points)\n", "#type(points)\n", "#len(points)\n", "#print(points)" ] }, { "cell_type": "code", "execution_count": 34, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Stima dell'integrale: 2.099909570041918\n", "Risultato esatto: 2.0943951023931953\n" ] } ], "source": [ "def my_f(point):\n", " return np.sqrt(1.- (point*point).sum(axis=1)).mean()\n", "\n", "base_area=np.pi\n", "print(\"Stima dell'integrale:\",base_area*my_f(points))\n", "print(\"Risultato esatto:\",2/3*np.pi)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Il grande vantaggio del metodo Monte Carlo è che si può estendere immediatamente a qualunque numero di dimensioni." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 11.2 Generazione di eventi distribuiti secondo una densità di probabilità fissata\n", "\n", "### Metodo hit and miss\n", "\n", "\n", "Supponiamo di avere la seguente distribuzione di probabilità:\n", "$$P(x) = \\begin{cases} 0.6 \\quad x<0.5\\\\ 1.7 \\quad x>0.5 \\end{cases}$$\n", "\n", "
\n", " \n", "
\n", "\n", "e di voler generare valori di $x$ distribuiti secondo $P(x)$\n", "\n", "Procediamo nel modo seguente:\n", "\n", "- Troviamo un numero $Y$ tale che $P(x) < Y, \\,\\forall x$. Nel nostro caso potremmo prendere $Y = 2$, per esempio. \n", "- Generiamo un numero casuale $x_0$ in [0,1], in modo uniforme.\n", "- Generiamo un secondo numero casuale $y_0$ in [0,$Y$]\n", "- Teniamo il numero $x_0$ nel sample se $y_0 < P(x_0)$ altrimenti lo scartiamo e generiamo un nuovo valore $x_0$.\n", "\n", "La probabilità di tenere un punto nel bin di sinistra è:\n", "$$P_{sx}=P(x < 0.5)\\cdot P(y < 0.6) = 0.5\\cdot\\frac{0.6}{2}$$\n", "La probabilità di tenere un punto nel bin di destra è:\n", "$$P_{dx}=P(x > 0.5)\\cdot P(y < 1.7) = 0.5\\cdot\\frac{1.4}{2}$$\n", "Quindi il rapporto $R$ fra le due probabilità, ovvero il rapporto fra il numero di punti accettati in $x < 0.5$ e il numero di punti accettati in $x > 0.5$ è, come deve:\n", "$$ R = \\frac{P_{sx}}{P_{dx}} = \\frac{0.6}{1.7}$$" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Numero di punti in [0,0.5]: 136\n", "Numero di punti in [0.5,1]: 362\n", "rapporto fra i numeri di punti: 2.661764705882353 rapporto atteso=7/3: 2.3333333333333335\n" ] } ], "source": [ "npoints=1000\n", "# npoints = 1000000\n", "xval=rng.uniform(size=npoints)\n", "# Separo i valori minori/maggiori di 0.5\n", "x1val=[x for x in xval if x < 0.5]\n", "x2val=[x for x in xval if x > 0.5]\n", "\n", "y1=0.6\n", "y2=1.4\n", "\n", "x1val=[x for x in x1val if y1 > rng.uniform()*2]\n", "x2val=[x for x in x2val if y2 > rng.uniform()*2]\n", "\n", "print(\"Numero di punti in [0,0.5]:\",len(x1val))\n", "print(\"Numero di punti in [0.5,1]:\",len(x2val))\n", "print(\"rapporto fra i numeri di punti:\",len(x2val)/len(x1val),\"rapporto atteso=7/3:\",7/3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Test della generazione Monte Carlo di punti secondo una distribuzione data\n", "\n", "#### Una funzione generica per generare eventi secondo una distribuzione qualsiasi" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "def generate_f(f,a,b,max):\n", " '''\n", " Generate random points distributed according to an input distribution\n", " Input:\n", " f: the probability distribution\n", " a,b: extrema of the range for the desired random number\n", " max: a real number such that f(x) < max for each x\n", " Returns:\n", " ranval: random number\n", " '''\n", " not_found=True\n", " while not_found:\n", " x=rng.uniform(a,b)\n", " y=rng.uniform(0,max)\n", " if f(x)>y:\n", " not_found=False\n", " \n", " return x\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Esempio 1" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [ "def f1(x):\n", " if x < 0.5:\n", " val=0.6\n", " else:\n", " val=1.4\n", " \n", " return val\n", " " ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Numero di punti in [0,0.5]: 300661\n", "Numero di punti in [0.5,1]: 699339\n", "rapporto fra i numeri di punti: 2.3260050355716237 rapporto atteso=7/3: 2.3333333333333335\n" ] } ], "source": [ "#npoints=1000\n", "npoints = 1000000\n", "xval = [generate_f(f1,0,1,2) for n in range(npoints)]\n", "x1val=[x for x in xval if x < 0.5]\n", "x2val=[x for x in xval if x > 0.5]\n", "\n", "print(\"Numero di punti in [0,0.5]:\",len(x1val))\n", "print(\"Numero di punti in [0.5,1]:\",len(x2val))\n", "print(\"rapporto fra i numeri di punti:\",len(x2val)/len(x1val),\"rapporto atteso=7/3:\",7/3)" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD4CAYAAADy46FuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAV2UlEQVR4nO3df6zd9X3f8edrdkPJMogBw5gNs1u8toAaLXiO125VVm/YSauYSSA5a4uVWbLKWJZNmxpopTE1shS0abSog8oKHoZFgOVmw1tHUwuWZVPBxKRJjKGUu9CBixs7NaOsFTQm7/1xPlc7vrn+3Ot77g8u9/mQjs73vL+fz/d+PjI6r/P9fr7nkKpCkqSz+QsLPQBJ0jubQSFJ6jIoJEldBoUkqcugkCR1LV/oAcy2Sy65pNasWbPQw5CkReWZZ575dlWtnGzfuy4o1qxZw+HDhxd6GJK0qCT532fb56UnSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUNWVQJNmT5ESSZyfZ9y+SVJJLhmq3JxlL8kKSzUP165IcafvuTpJWPy/JI61+KMmaoT7bk7zYHttHnq0k6ZxN54zifmDLxGKSK4C/B7w8VLsa2AZc0/rck2RZ230vsBNY1x7jx9wBvFZVVwF3AXe2Y10E3AF8CNgA3JFkxblNT5I0qim/mV1VXx7+lD/kLuAXgEeHaluBh6vqLeClJGPAhiR/AFxQVU8CJHkAuAF4rPX5V63/fuDX2tnGZuBgVZ1qfQ4yCJeHzm2K0jvHmtt+c0H+7h989qcW5O/q3WFGP+GR5GPAH1bV19sVpHGrgKeGXh9rte+07Yn18T6vAFTV6SSvAxcP1yfpM3E8OxmcrXDllVfOZErSu9pCBdRCWshwfLd9IDjnoEjyXuCXgOsn2z1JrTr1mfY5s1i1G9gNsH79ev/frpKWZDjOlZnc9fSDwFrg6+2S0mrgq0n+MoNP/VcMtV0NvNrqqyepM9wnyXLgQuBU51iSpHl0zkFRVUeq6tKqWlNVaxi8oX+wqv4IOABsa3cyrWWwaP10VR0H3kiysa0/3Mz/X9s4AIzf0XQj8ERVFfBF4PokK9oi9vWtJkmaR1NeekryEPBh4JIkx4A7quq+ydpW1dEk+4DngNPArVX1dtt9C4M7qM5nsIj9WKvfBzzYFr5PMbhriqo6leQzwFdau18eX9iWJM2f6dz19PEp9q+Z8HoXsGuSdoeBayepvwncdJZj7wH2TDVGSdLc8ZvZkqQug0KS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklS15RBkWRPkhNJnh2q/eskv5fkG0n+Y5L3D+27PclYkheSbB6qX5fkSNt3d5K0+nlJHmn1Q0nWDPXZnuTF9tg+W5OWJE3fdM4o7ge2TKgdBK6tqh8Ffh+4HSDJ1cA24JrW554ky1qfe4GdwLr2GD/mDuC1qroKuAu4sx3rIuAO4EPABuCOJCvOfYqSpFFMGRRV9WXg1ITab1fV6fbyKWB1294KPFxVb1XVS8AYsCHJ5cAFVfVkVRXwAHDDUJ+9bXs/sKmdbWwGDlbVqap6jUE4TQwsSdIcm401in8IPNa2VwGvDO071mqr2vbE+hl9Wvi8DlzcOdb3SLIzyeEkh0+ePDnSZCRJZxopKJL8EnAa+Px4aZJm1anPtM+ZxardVbW+qtavXLmyP2hJ0jmZcVC0xeWfBn6mXU6Cwaf+K4aarQZebfXVk9TP6JNkOXAhg0tdZzuWJGkezSgokmwBPg18rKr+bGjXAWBbu5NpLYNF66er6jjwRpKNbf3hZuDRoT7jdzTdCDzRgueLwPVJVrRF7OtbTZI0j5ZP1SDJQ8CHgUuSHGNwJ9LtwHnAwXaX61NV9fNVdTTJPuA5Bpekbq2qt9uhbmFwB9X5DNY0xtc17gMeTDLG4ExiG0BVnUryGeArrd0vV9UZi+qSpLk3ZVBU1ccnKd/Xab8L2DVJ/TBw7ST1N4GbznKsPcCeqcYoSZo7fjNbktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnLoJAkdRkUkqQug0KS1GVQSJK6DApJUpdBIUnqMigkSV1TBkWSPUlOJHl2qHZRkoNJXmzPK4b23Z5kLMkLSTYP1a9LcqTtuztJWv28JI+0+qEka4b6bG9/48Uk22dt1pKkaZvOGcX9wJYJtduAx6tqHfB4e02Sq4FtwDWtzz1JlrU+9wI7gXXtMX7MHcBrVXUVcBdwZzvWRcAdwIeADcAdw4EkSZofUwZFVX0ZODWhvBXY27b3AjcM1R+uqreq6iVgDNiQ5HLggqp6sqoKeGBCn/Fj7Qc2tbONzcDBqjpVVa8BB/newJIkzbGZrlFcVlXHAdrzpa2+CnhlqN2xVlvVtifWz+hTVaeB14GLO8f6Hkl2Jjmc5PDJkydnOCVJ0mRmezE7k9SqU59pnzOLVburan1VrV+5cuW0BipJmp6ZBsW32uUk2vOJVj8GXDHUbjXwaquvnqR+Rp8ky4ELGVzqOtuxJEnzaKZBcQAYvwtpO/DoUH1bu5NpLYNF66fb5ak3kmxs6w83T+gzfqwbgSfaOsYXgeuTrGiL2Ne3miRpHi2fqkGSh4APA5ckOcbgTqTPAvuS7ABeBm4CqKqjSfYBzwGngVur6u12qFsY3EF1PvBYewDcBzyYZIzBmcS2dqxTST4DfKW1++WqmrioLkmaY1MGRVV9/Cy7Np2l/S5g1yT1w8C1k9TfpAXNJPv2AHumGqMkae74zWxJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKlrpKBI8s+SHE3ybJKHknx/kouSHEzyYnteMdT+9iRjSV5Isnmofl2SI23f3UnS6ucleaTVDyVZM8p4JUnnbsZBkWQV8E+A9VV1LbAM2AbcBjxeVeuAx9trklzd9l8DbAHuSbKsHe5eYCewrj22tPoO4LWqugq4C7hzpuOVJM3MqJeelgPnJ1kOvBd4FdgK7G379wI3tO2twMNV9VZVvQSMARuSXA5cUFVPVlUBD0zoM36s/cCm8bMNSdL8mHFQVNUfAv8GeBk4DrxeVb8NXFZVx1ub48Clrcsq4JWhQxxrtVVte2L9jD5VdRp4Hbh44liS7ExyOMnhkydPznRKkqRJjHLpaQWDT/xrgb8C/MUkP9vrMkmtOvVenzMLVburan1VrV+5cmV/4JKkczLKpae/C7xUVSer6jvAF4AfA77VLifRnk+09seAK4b6r2ZwqepY255YP6NPu7x1IXBqhDFLks7RKEHxMrAxyXvbusEm4HngALC9tdkOPNq2DwDb2p1MaxksWj/dLk+9kWRjO87NE/qMH+tG4Im2jiFJmifLZ9qxqg4l2Q98FTgN/C6wG3gfsC/JDgZhclNrfzTJPuC51v7Wqnq7He4W4H7gfOCx9gC4D3gwyRiDM4ltMx2vJGlmZhwUAFV1B3DHhPJbDM4uJmu/C9g1Sf0wcO0k9TdpQSNJWhh+M1uS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoaKSiSvD/J/iS/l+T5JH8zyUVJDiZ5sT2vGGp/e5KxJC8k2TxUvy7Jkbbv7iRp9fOSPNLqh5KsGWW8kqRzN+oZxa8Cv1VVPwx8AHgeuA14vKrWAY+31yS5GtgGXANsAe5Jsqwd515gJ7CuPba0+g7gtaq6CrgLuHPE8UqSztGMgyLJBcBPAPcBVNWfV9X/AbYCe1uzvcANbXsr8HBVvVVVLwFjwIYklwMXVNWTVVXAAxP6jB9rP7Bp/GxDkjQ/lo/Q9weAk8C/T/IB4BngU8BlVXUcoKqOJ7m0tV8FPDXU/1irfadtT6yP93mlHet0kteBi4FvjzDud6Q1t/3mQg9BkiY1SlAsBz4IfLKqDiX5VdplprOY7EygOvVenzMPnOxkcOmKK6+8sjfmKfmGLUlnGmWN4hhwrKoOtdf7GQTHt9rlJNrziaH2Vwz1Xw282uqrJ6mf0SfJcuBC4NTEgVTV7qpaX1XrV65cOcKUJEkTzTgoquqPgFeS/FArbQKeAw4A21ttO/Bo2z4AbGt3Mq1lsGj9dLtM9UaSjW394eYJfcaPdSPwRFvHkCTNk1EuPQF8Evh8kvcA3wQ+wSB89iXZAbwM3ARQVUeT7GMQJqeBW6vq7XacW4D7gfOBx9oDBgvlDyYZY3AmsW3E8UqSztFIQVFVXwPWT7Jr01na7wJ2TVI/DFw7Sf1NWtBIkhaG38yWJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnLoJAkdRkUkqQug0KS1GVQSJK6DApJUtfIQZFkWZLfTfJf2uuLkhxM8mJ7XjHU9vYkY0leSLJ5qH5dkiNt391J0urnJXmk1Q8lWTPqeCVJ52Y2zig+BTw/9Po24PGqWgc83l6T5GpgG3ANsAW4J8my1udeYCewrj22tPoO4LWqugq4C7hzFsYrSToHIwVFktXATwGfGypvBfa27b3ADUP1h6vqrap6CRgDNiS5HLigqp6sqgIemNBn/Fj7gU3jZxuSpPkx6hnFrwC/AHx3qHZZVR0HaM+Xtvoq4JWhdsdabVXbnlg/o09VnQZeBy6eOIgkO5McTnL45MmTI05JkjRsxkGR5KeBE1X1zHS7TFKrTr3X58xC1e6qWl9V61euXDnN4UiSpmP5CH1/HPhYko8C3w9ckOQ/AN9KcnlVHW+XlU609seAK4b6rwZebfXVk9SH+xxLshy4EDg1wpglSedoxmcUVXV7Va2uqjUMFqmfqKqfBQ4A21uz7cCjbfsAsK3dybSWwaL10+3y1BtJNrb1h5sn9Bk/1o3tb3zPGYUkae6MckZxNp8F9iXZAbwM3ARQVUeT7AOeA04Dt1bV263PLcD9wPnAY+0BcB/wYJIxBmcS2+ZgvJKkjlkJiqr6EvCltv3HwKaztNsF7Jqkfhi4dpL6m7SgkSQtDL+ZLUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnLoJAkdc04KJJckeS/JXk+ydEkn2r1i5IcTPJie14x1Of2JGNJXkiyeah+XZIjbd/dSdLq5yV5pNUPJVkzwlwlSTMwyhnFaeCfV9WPABuBW5NcDdwGPF5V64DH22vavm3ANcAW4J4ky9qx7gV2AuvaY0ur7wBeq6qrgLuAO0cYryRpBmYcFFV1vKq+2rbfAJ4HVgFbgb2t2V7ghra9FXi4qt6qqpeAMWBDksuBC6rqyaoq4IEJfcaPtR/YNH62IUmaH7OyRtEuCf114BBwWVUdh0GYAJe2ZquAV4a6HWu1VW17Yv2MPlV1GngduHiSv78zyeEkh0+ePDkbU5IkNSMHRZL3Ab8B/NOq+pNe00lq1an3+pxZqNpdVeurav3KlSunGrIk6RyMFBRJvo9BSHy+qr7Qyt9ql5Nozyda/RhwxVD31cCrrb56kvoZfZIsBy4ETo0yZknSuRnlrqcA9wHPV9W/Hdp1ANjetrcDjw7Vt7U7mdYyWLR+ul2eeiPJxnbMmyf0GT/WjcATbR1DkjRPlo/Q98eBnwOOJPlaq/0i8FlgX5IdwMvATQBVdTTJPuA5BndM3VpVb7d+twD3A+cDj7UHDILowSRjDM4kto0wXknSDMw4KKrqfzL5GgLAprP02QXsmqR+GLh2kvqbtKCRJC0Mv5ktSeoyKCRJXQaFJKnLoJAkdRkUkqQug0KS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpK5FERRJtiR5IclYktsWejyStJS844MiyTLg3wEfAa4GPp7k6oUdlSQtHe/4oAA2AGNV9c2q+nPgYWDrAo9JkpaM5Qs9gGlYBbwy9PoY8KHhBkl2Ajvby/+b5IUR/t4lwLdH6L8YLbU5L7X5gnNeEnLnSHP+q2fbsRiCIpPU6owXVbuB3bPyx5LDVbV+No61WCy1OS+1+YJzXirmas6L4dLTMeCKodergVcXaCyStOQshqD4CrAuydok7wG2AQcWeEyStGS84y89VdXpJP8Y+CKwDNhTVUfn8E/OyiWsRWapzXmpzRec81IxJ3NOVU3dSpK0ZC2GS0+SpAVkUEiSupZkUEz1kyAZuLvt/0aSDy7EOGfTNOb8M22u30jyO0k+sBDjnE3T/emXJH8jydtJbpzP8c2F6cw5yYeTfC3J0ST/fb7HONum8d/2hUn+c5Kvtzl/YiHGOVuS7ElyIsmzZ9k/++9fVbWkHgwWxP8X8APAe4CvA1dPaPNR4DEG3+HYCBxa6HHPw5x/DFjRtj+yFOY81O4J4L8CNy70uOfh3/n9wHPAle31pQs97nmY8y8Cd7btlcAp4D0LPfYR5vwTwAeBZ8+yf9bfv5biGcV0fhJkK/BADTwFvD/J5fM90Fk05Zyr6neq6rX28ikG31dZzKb70y+fBH4DODGfg5sj05nzPwC+UFUvA1TVYp/3dOZcwF9KEuB9DILi9PwOc/ZU1ZcZzOFsZv39aykGxWQ/CbJqBm0Wk3Odzw4Gn0gWsynnnGQV8PeBX5/Hcc2l6fw7/zVgRZIvJXkmyc3zNrq5MZ05/xrwIwy+qHsE+FRVfXd+hrcgZv396x3/PYo5MOVPgkyzzWIy7fkk+TsMguJvzemI5t505vwrwKer6u3Bh81FbzpzXg5cB2wCzgeeTPJUVf3+XA9ujkxnzpuBrwE/CfwgcDDJ/6iqP5njsS2UWX//WopBMZ2fBHm3/WzItOaT5EeBzwEfqao/nqexzZXpzHk98HALiUuAjyY5XVX/aV5GOPum+9/2t6vqT4E/TfJl4APAYg2K6cz5E8Bna3ABfyzJS8APA0/PzxDn3ay/fy3FS0/T+UmQA8DN7e6BjcDrVXV8vgc6i6acc5IrgS8AP7eIP10Om3LOVbW2qtZU1RpgP/CPFnFIwPT+234U+NtJlid5L4NfYn5+nsc5m6Yz55cZnEGR5DLgh4Bvzuso59esv38tuTOKOstPgiT5+bb/1xncAfNRYAz4MwafSBatac75XwIXA/e0T9inaxH/8uY05/yuMp05V9XzSX4L+AbwXeBzVTXpbZaLwTT/nT8D3J/kCIPLMp+uqkX78+NJHgI+DFyS5BhwB/B9MHfvX/6EhySpayleepIknQODQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnr/wEmo+c/ft+QMAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "nbins = 10\n", "xrange = (0,1)\n", "\n", "fig, ax = plt.subplots()\n", "nevent, bins, patches = ax.hist(xval, nbins, range=xrange)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Esempio 2" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "def f2(x):\n", " val=3./8.*(1+x**2)\n", " \n", " return val " ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlwAAAGbCAYAAAARGU4hAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAe7klEQVR4nO3df4zc9Z3f8ef77MA5djDmSPYcTGtXddMD3OTiFeUu5bou9NgELqZSqRylwbRU7iFyTSRSYVqpvepk1W3FSQcEVDdENoXL1soltZXEl3LuWVF7EM6OSBxDOJzgEIPP1vHDwSniavruH/PlMmfWu7Prfc/O7vf5kEYz857vd76f93736335+/3OdyIzkSRJUp2fme0BSJIkzXcGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGILZ3sAk7n44otz5cqVpcv4yU9+wuLFi0uXMaja3Du0u/829w7t7r/NvUO7+7f3+t4PHDjwZ5n57jPrAx+4Vq5cyf79+0uXsW/fPkZGRkqXMaja3Du0u/829w7t7r/NvUO7+7f3kfLlRMQPx6t7SFGSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKmYgUuSJKnYwtkegCRJareVm79avozto4vLlzER93BJkiQVM3BJkiQVM3BJkiQVM3BJkiQVM3BJkiQVM3BJkiQVM3BJkiQVmzRwRcT7IuLJrtuPI+LTEXFRRDwaEc8298u65rkrIg5HxDMRcV1XfW1EHGxeuycioqoxSZKkQTHphU8z8xngAwARsQB4AfgysBnYm5lbI2Jz8/zOiLgM2ABcDrwX+IOI+BuZ+SbwALAJeBz4GjAK7Jnppqbq4AsnuaX4omtHtl5f+v6SJGlwTfWQ4jXA9zPzh8B6YEdT3wHc2DxeD4xl5huZ+RxwGLgyIpYDF2TmY5mZwENd80iSJM1b0ck+PU4c8XngW5l5X0S8mpkXdr32SmYui4j7gMcz8+Gm/iCdvVhHgK2ZeW1Tvxq4MzNvGGc5m+jsCWNoaGjt2NjYdPvryYmXT3L89dJFsOaSpbULmKZTp06xZMmS2R7GrGlz/23uHdrdf5t7h3b3P6i9H3zhZPkyVi1d0Jfe161bdyAzh8+s9/xdihFxHvBR4K7JJh2nlhPU317M3AZsAxgeHs6RkZFehzkt9z6yi7sP1n6t5JGPj5S+/3Tt27eP6p/vIGtz/23uHdrdf5t7h3b3P6i9V5/WA53vUpzN3qdySPHDdPZuHW+eH28OE9Lcn2jqR4FLu+ZbAbzY1FeMU5ckSZrXphK4PgZ8oev5bmBj83gjsKurviEizo+IVcBq4InMPAa8FhFXNZ9OvLlrHkmSpHmrp+NoEfFO4O8D/7yrvBXYGRG3As8DNwFk5qGI2Ak8BZwGbm8+oQhwG7AdWETnvK5Z/4SiJElStZ4CV2b+H+Dnzqi9ROdTi+NNvwXYMk59P3DF1IcpSZI0d3mleUmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGIGLkmSpGILZ3sAbbFy81f7spwjW6/vy3IkSVLvDFySJOms+rXDYL7zkKIkSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVKxngJXRFwYEV+MiO9FxNMR8UsRcVFEPBoRzzb3y7qmvysiDkfEMxFxXVd9bUQcbF67JyKioilJkqRB0usert8Bfj8z/ybwfuBpYDOwNzNXA3ub50TEZcAG4HJgFLg/IhY07/MAsAlY3dxGZ6gPSZKkgTVp4IqIC4BfAR4EyMw/z8xXgfXAjmayHcCNzeP1wFhmvpGZzwGHgSsjYjlwQWY+lpkJPNQ1jyRJ0rwVnewzwQQRHwC2AU/R2bt1APgU8EJmXtg13SuZuSwi7gMez8yHm/qDwB7gCLA1M69t6lcDd2bmDeMscxOdPWEMDQ2tHRsbO7cuJ3Hi5ZMcf710EX2z5pKlU5r+1KlTLFmypGg0g6/N/be5d2h3/23uHdrd/3R6P/jCyaLR9NeqpQv6st7XrVt3IDOHz6wv7GHehcAHgd/IzG9GxO/QHD48i/HOy8oJ6m8vZm6jE/IYHh7OkZGRHoY5ffc+sou7D/byoxh8Rz4+MqXp9+3bR/XPd5C1uf829w7t7r/NvUO7+59O77ds/mrNYPps++jiWV3vvZzDdRQ4mpnfbJ5/kU4AO94cJqS5P9E1/aVd868AXmzqK8apS5IkzWuTBq7M/FPgRxHxvqZ0DZ3Di7uBjU1tI7Crebwb2BAR50fEKjonxz+RmceA1yLiqubTiTd3zSNJkjRv9Xoc7TeARyLiPOAHwD+hE9Z2RsStwPPATQCZeSgidtIJZaeB2zPzzeZ9bgO2A4vonNe1Z4b6kCSpVVZO41DfHWtOz5tDhHNNT4ErM58E3nYCGJ29XeNNvwXYMk59P3DFFManKZrqBjidje/I1uunNL0kSW3nleYlSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKKLZztAWjuWbn5q31ZzpGt1/dlOZIkVTNwSZI0w/r1H1PNHR5SlCRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKtZT4IqIIxFxMCKejIj9Te2iiHg0Ip5t7pd1TX9XRByOiGci4rqu+trmfQ5HxD0RETPfkiRJ0mCZyh6udZn5gcwcbp5vBvZm5mpgb/OciLgM2ABcDowC90fEgmaeB4BNwOrmNnruLUiSJA22czmkuB7Y0TzeAdzYVR/LzDcy8zngMHBlRCwHLsjMxzIzgYe65pEkSZq3opN9Jpko4jngFSCB/5yZ2yLi1cy8sGuaVzJzWUTcBzyemQ839QeBPcARYGtmXtvUrwbuzMwbxlneJjp7whgaGlo7NjZ2bl1O4sTLJzn+eukiBtbQIga29zWXLC1fxqlTp1iyZEn5cgZRm3uHdvff5t6hP/0ffOFk6ftP1yD/m19t1dIFffm9X7du3YGuo4F/YWGP838oM1+MiPcAj0bE9yaYdrzzsnKC+tuLmduAbQDDw8M5MjLS4zCn595HdnH3wV5/FPPLHWtOD2zvRz4+Ur6Mffv2Uf37Naja3Du0u/829w796f+WzV8tff/pGuR/86ttH108q7/3PR1SzMwXm/sTwJeBK4HjzWFCmvsTzeRHgUu7Zl8BvNjUV4xTlyRJmtcmDVwRsTgi3vXWY+BXge8Cu4GNzWQbgV3N493Ahog4PyJW0Tk5/onMPAa8FhFXNZ9OvLlrHkmSpHmrl/2KQ8CXmys4LAR+NzN/PyL+GNgZEbcCzwM3AWTmoYjYCTwFnAZuz8w3m/e6DdgOLKJzXteeGexFkiRpIE0auDLzB8D7x6m/BFxzlnm2AFvGqe8Hrpj6MCVJkuaudp45pzlhZR9OOt0+urh8GZIk+dU+kiRJxQxckiRJxQxckiRJxQxckiRJxQxckiRJxQxckiRJxQxckiRJxQxckiRJxbzwqSSpVQ6+cJJb+nBhZambe7gkSZKKGbgkSZKKGbgkSZKKGbgkSZKKGbgkSZKK+SlFSdJAWNmnTw7esaYvi5H+EvdwSZIkFTNwSZIkFTNwSZIkFTNwSZIkFfOkebVaP77i48jW60vfX5I0+NzDJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMwLn0qSJrWy+ALB0nznHi5JkqRiBi5JkqRiHlKUivXrUIzf2ShJg8s9XJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScV6vtJ8RCwA9gMvZOYNEXER8N+AlcAR4B9l5ivNtHcBtwJvAv8iM7/e1NcC24FFwNeAT2VmzlQzUptN9Yr2d6w5zS1TnMer2UvS9ExlD9engKe7nm8G9mbmamBv85yIuAzYAFwOjAL3N2EN4AFgE7C6uY2e0+glSZLmgJ4CV0SsAK4HPtdVXg/saB7vAG7sqo9l5huZ+RxwGLgyIpYDF2TmY81erYe65pEkSZq3opcjehHxReDfA+8CPtMcUnw1My/smuaVzFwWEfcBj2fmw039QWAPncOOWzPz2qZ+NXBnZt4wzvI20dkTxtDQ0NqxsbFz63ISJ14+yfHXSxcxsIYW0dreod39T6f3NZcsrRnMLDh16hRLliyZ7WHMiun0fvCFk0Wj6T+3+9kexexYtXRBX7b5devWHcjM4TPrk57DFRE3ACcy80BEjPSwrBinlhPU317M3AZsAxgeHs6RkV4WO333PrKLuw/2fDrbvHLHmtOt7R3a3f90ej/y8ZGawcyCffv2Uf1vy6CaTu9TPd9vkLndt7P37aOLZ3Wb7+Wn/iHgoxHxEeBngQsi4mHgeEQsz8xjzeHCE830R4FLu+ZfAbzY1FeMU5ckTdNUPywB0/vAhKRzM2ngysy7gLsAmj1cn8nMfxwR/wnYCGxt7nc1s+wGfjcifht4L52T45/IzDcj4rWIuAr4JnAzcO/MtiOp0nT+uE+Vn4SUNB+dy37FrcDOiLgVeB64CSAzD0XETuAp4DRwe2a+2cxzGz+9LMSe5iZJkjSvTSlwZeY+YF/z+CXgmrNMtwXYMk59P3DFVAcpSXNRP/YISpobvNK8JElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSsXZe319S6x184aRXW5fUN+7hkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKmbgkiRJKuZ1uCQNlJV9ujbWHWv6shhJAtzDJUmSVM7AJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVMzAJUmSVGzSwBURPxsRT0TEtyPiUET8u6Z+UUQ8GhHPNvfLuua5KyIOR8QzEXFdV31tRBxsXrsnIqKmLUmSpMHRyx6uN4C/l5nvBz4AjEbEVcBmYG9mrgb2Ns+JiMuADcDlwChwf0QsaN7rAWATsLq5jc5cK5IkSYNp0sCVHaeap+9obgmsB3Y09R3Ajc3j9cBYZr6Rmc8Bh4ErI2I5cEFmPpaZCTzUNY8kSdK8FZ3sM8lEnT1UB4C/Dnw2M++MiFcz88KuaV7JzGURcR/weGY+3NQfBPYAR4CtmXltU78auDMzbxhneZvo7AljaGho7djY2Ll1OYkTL5/k+OulixhYQ4tobe/Q7v7b3Du0u/829w7t7r/Nva9auoAlS5aUL2fdunUHMnP4zPrCXmbOzDeBD0TEhcCXI+KKCSYf77ysnKA+3vK2AdsAhoeHc2RkpJdhTtu9j+zi7oM9/SjmnTvWnG5t79Du/tvcO7S7/zb3Du3uv829bx9dTHWemMiUPqWYma8C++ice3W8OUxIc3+imewocGnXbCuAF5v6inHqkiRJ81ovn1J8d7Nni4hYBFwLfA/YDWxsJtsI7Goe7wY2RMT5EbGKzsnxT2TmMeC1iLiq+XTizV3zSJIkzVu97FdcDuxozuP6GWBnZn4lIh4DdkbErcDzwE0AmXkoInYCTwGngdubQ5IAtwHbgUV0zuvaM5PNSJIkDaJJA1dmfgf4xXHqLwHXnGWeLcCWcer7gYnO/5IkSZp3vNK8JElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSMQOXJElSsUkDV0RcGhF/GBFPR8ShiPhUU78oIh6NiGeb+2Vd89wVEYcj4pmIuK6rvjYiDjav3RMRUdOWJEnS4OhlD9dp4I7M/AXgKuD2iLgM2AzszczVwN7mOc1rG4DLgVHg/ohY0LzXA8AmYHVzG53BXiRJkgbSpIErM49l5reax68BTwOXAOuBHc1kO4Abm8frgbHMfCMznwMOA1dGxHLggsx8LDMTeKhrHkmSpHkrOtmnx4kjVgLfAK4Ans/MC7teeyUzl0XEfcDjmflwU38Q2AMcAbZm5rVN/Wrgzsy8YZzlbKKzJ4yhoaG1Y2Nj02quVydePsnx10sXMbCGFtHa3qHd/be5d2h3/23uHdrdf5t7X7V0AUuWLClfzrp16w5k5vCZ9YW9vkFELAF+D/h0Zv54gtOvxnshJ6i/vZi5DdgGMDw8nCMjI70Oc1rufWQXdx/s+Ucxr9yx5nRre4d299/m3qHd/be5d2h3/23uffvoYqrzxER6+pRiRLyDTth6JDO/1JSPN4cJae5PNPWjwKVds68AXmzqK8apS5IkzWu9fEoxgAeBpzPzt7te2g1sbB5vBHZ11TdExPkRsYrOyfFPZOYx4LWIuKp5z5u75pEkSZq3etmv+CHgE8DBiHiyqf0rYCuwMyJuBZ4HbgLIzEMRsRN4is4nHG/PzDeb+W4DtgOL6JzXtWdm2pAkSRpckwauzPxfjH/+FcA1Z5lnC7BlnPp+OifcS5IktYZXmpckSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSpm4JIkSSo2aeCKiM9HxImI+G5X7aKIeDQinm3ul3W9dldEHI6IZyLiuq762og42Lx2T0TEzLcjSZI0eHrZw7UdGD2jthnYm5mrgb3NcyLiMmADcHkzz/0RsaCZ5wFgE7C6uZ35npIkSfPSpIErM78BvHxGeT2wo3m8A7ixqz6WmW9k5nPAYeDKiFgOXJCZj2VmAg91zSNJkjSvTfccrqHMPAbQ3L+nqV8C/KhruqNN7ZLm8Zl1SZKkeW/hDL/feOdl5QT18d8kYhOdw48MDQ2xb9++GRnc2QwtgjvWnC5dxqBqc+/Q7v7b3Du0u/829w7t7r/NvZ86dao8T0xkuoHreEQsz8xjzeHCE039KHBp13QrgBeb+opx6uPKzG3ANoDh4eEcGRmZ5jB7c+8ju7j74Exnz7nhjjWnW9s7tLv/NvcO7e6/zb1Du/tvc+/bRxdTnScmMt1DiruBjc3jjcCurvqGiDg/IlbROTn+ieaw42sRcVXz6cSbu+aRJEma1yaNuRHxBWAEuDgijgL/FtgK7IyIW4HngZsAMvNQROwEngJOA7dn5pvNW91G5xOPi4A9zU2SJGnemzRwZebHzvLSNWeZfguwZZz6fuCKKY1OkiRpHvBK85IkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScUMXJIkScX6HrgiYjQinomIwxGxud/LlyRJ6re+Bq6IWAB8FvgwcBnwsYi4rJ9jkCRJ6rd+7+G6EjicmT/IzD8HxoD1fR6DJElSX0Vm9m9hEf8QGM3Mf9Y8/wTwtzPzk2dMtwnY1Dx9H/BM8dAuBv6seBmDqs29Q7v7b3Pv0O7+29w7tLt/e6/3VzPz3WcWF/Zhwd1inNrbEl9mbgO21Q+nIyL2Z+Zwv5Y3SNrcO7S7/zb3Du3uv829Q7v7t/fZ673fhxSPApd2PV8BvNjnMUiSJPVVvwPXHwOrI2JVRJwHbAB293kMkiRJfdXXQ4qZeToiPgl8HVgAfD4zD/VzDGfRt8OXA6jNvUO7+29z79Du/tvcO7S7f3ufJX09aV6SJKmNvNK8JElSMQOXJElSsdYEroi4KSIORcT/i4izfiz0bF89FBEXRcSjEfFsc7+sPyM/d72MPSLeFxFPdt1+HBGfbl77zYh4oeu1j/S9iXPQ67qLiCMRcbDpcf9U5x9EPa77SyPiDyPi6WYb+VTXa3Nu3U/29WHRcU/z+nci4oO9zjsX9ND/x5u+vxMRfxQR7+96bdxtYK7oofeRiDjZ9fv8b3qdd9D10Pu/7Or7uxHxZkRc1Lw219f75yPiRER89yyvD8Y2n5mtuAG/QOciqvuA4bNMswD4PvDXgPOAbwOXNa/9R2Bz83gz8B9mu6cp9D6lsTc/hz+lc/E2gN8EPjPbfVT3DxwBLj7Xn98g3XoZO7Ac+GDz+F3An3T93s+pdT/RNtw1zUeAPXSuC3gV8M1e5x30W4/9/zKwrHn84bf6b56Puw3MhVuPvY8AX5nOvIN8m+r4gV8D/ud8WO/N+H8F+CDw3bO8PhDbfGv2cGXm05k52RXrJ/rqofXAjubxDuDGkoHWmOrYrwG+n5k/rBxUH53rupvX6z4zj2Xmt5rHrwFPA5f0a4AzrJevD1sPPJQdjwMXRsTyHucddJP2kJl/lJmvNE8fp3M9xPngXNbfXF/3Ux3/x4Av9GVkfZCZ3wBenmCSgdjmWxO4enQJ8KOu50f56R+eocw8Bp0/UMB7+jy2czHVsW/g7RvjJ5tdsZ+fS4fUGr32n8D/iIgD0fl6qanOP4imNPaIWAn8IvDNrvJcWvcTbcOTTdPLvINuqj3cSud//m852zYwF/Ta+y9FxLcjYk9EXD7FeQdVz+OPiHcCo8DvdZXn8nrvxUBs8/3+ap9SEfEHwM+P89K/zsxdvbzFOLU5cd2MiXqf4vucB3wUuKur/ADwW3R+Fr8F3A380+mNtMYM9f+hzHwxIt4DPBoR32v+5zTQZnDdL6Hzj/CnM/PHTXng1/0ZetmGzzbNnN3+u/TcQ0SsoxO4/k5XeU5uA41eev8WnVMlTjXnI/53YHWP8w6yqYz/14D/nZnde4Tm8nrvxUBs8/MqcGXmtef4FhN99dDxiFiemceaXZEnznFZM2qi3iNiKmP/MPCtzDze9d5/8Tgi/gvwlZkY80yaif4z88Xm/kREfJnO7uZv0IJ1HxHvoBO2HsnML3W998Cv+zP08vVhZ5vmvB7mHXQ9fX1aRPwt4HPAhzPzpbfqE2wDc8GkvXf9R4LM/FpE3B8RF/cy74CbyvjfdgRjjq/3XgzENu8hxb9soq8e2g1sbB5vBHrZYzYopjL2tx3bb/5Qv+UfAON+EmSATdp/RCyOiHe99Rj4VX7a57xe9xERwIPA05n522e8NtfWfS9fH7YbuLn55NJVwMnmcOt8+OqxSXuIiL8CfAn4RGb+SVd9om1gLuil959vft+JiCvp/A18qZd5B1xP44+IpcDfpevfgXmw3nsxGNt81dn4g3aj88fiKPAGcBz4elN/L/C1ruk+QudTWt+ncyjyrfrPAXuBZ5v7i2a7pyn0Pu7Yx+n9nXT+8Vl6xvz/FTgIfKf5ZVw+2z3NdP90PqXy7eZ2qE3rns4hpWzW75PN7SNzdd2Ptw0Dvw78evM4gM82rx+k61PLZ9v+59Kth/4/B7zSta73N/WzbgNz5dZD759sevs2nQ8M/PJ8WfeT9d48vwUYO2O++bDevwAcA/4vnb/ztw7iNu9X+0iSJBXzkKIkSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVIxA5ckSVKx/w8KW4jjIz0VRAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "#npoints=1000\n", "npoints = 100000\n", "xval = [generate_f(f2,-1,1,1) for n in range(npoints)]\n", "\n", "nbins = 20\n", "xrange = (-1,1)\n", "\n", "fig, ax = plt.subplots(figsize=(10,7))\n", "nevent, bins, patches = ax.hist(xval, nbins, range=xrange)\n", "ax.grid(True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 11.3 Propagazione dell'errore con metodo MC" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Supponiamo di aver fatto un esperimento per misurare l'accelerazione di gravità $g$ misurando il periodo $T$ di oscillazione di un pendolo di lunghezza variabile $L$. Sappiamo che:\n", "$$ T = 2 \\pi \\sqrt{\\frac{L}{g}} \\quad \\rightarrow \\quad g = \\frac{ 4 \\pi^2 L}{T^2}$$\n", "Supponiamo di aver misurato $L = 1.000 \\pm 0.004\\, m$, $T = 2.000 \\pm 0.005\\, s$ e che le misure siano distribuite gaussianamente.\n", "Per la propagazine degli errori ci servono:\n", "$$ \\frac{\\partial g}{\\partial T} = \\frac{ - 8 \\pi^2 L}{T^3},\\qquad \\frac{\\partial g}{\\partial L} = \\frac{ 4 \\pi^2}{T^2}$$\n", "\n", "$$ \\sigma_g = \\sqrt{\\left. \\sigma_T \\, \\frac{\\partial g}{\\partial T}\\right|^2_{\\overline{T},\\overline{L}} +\n", "\\left. \\sigma_L \\, \\frac{\\partial g}{\\partial L}\\right|^2_{\\overline{T},\\overline{L}}\n", "}$$" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "g_exp = 9.869604401089358 +\\- 0.06319630315448918\n" ] } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "def g(L,T):\n", " return 4*np.pi**2*L/T**2\n", "\n", "def g_err(L,L_err,T,T_err):\n", " return 4*np.pi**2/T**2*np.sqrt(4*(L/T*T_err)**2 + L_err**2)\n", "\n", "L_mean =1.0\n", "L_err = 0.004\n", "T_mean = 2.0\n", "T_err = 0.005\n", "\n", "g_m = g(L_mean,T_mean)\n", "g_e = g_err(L_mean,L_err,T_mean,T_err)\n", "print(f\"g_exp = {g_m} +\\- {g_e}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Per valutare l'errore su $g$ con il metodo Monte Carlo generiamo 1000 valori per L e 1000 per T da distribuzioni gaussiane con i rispettivi valori medi e deviazioni standard. Calcoliamo i corrispondenti valori di $g$ e valutiamo il valor medio e la deviazione standard della distribuzione ottenuta." ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "g_MC = 9.867338107189237 +\\- 0.06329906777231029\n" ] } ], "source": [ "npoints = 1000\n", "T_test = rng.normal(T_mean,T_err,npoints)\n", "L_test = rng.normal(L_mean,L_err,npoints)\n", "g_test = g(L_test,T_test)\n", "g_m1 = g_test.mean()\n", "g_e1 = g_test.std()\n", "print(f\"g_MC = {g_m1} +\\- {g_e1}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plot di controllo\n", "\n", "Nota sulle normalizzazioni:
\n", " - Dato un istogramma con $N_{tot}$ valori, i valori in ciascun bin sommano a $N_{tot}$. \n", " - Se si vuole sovrapporre una gaussiana $g(x)$ all'istogramma bisogna ricordare che: $\\int g(x)\\,dx = \\sum_i g(x_i) w =1$, dove $w$ è la larghezza di ciascun bin. Di conseguenza i valori di $g(x)$ vanno moltiplicati per $w \\cdot N_{tot}$ per avere un confronto significativo." ] }, { "cell_type": "code", "execution_count": 68, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 68, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAOxElEQVR4nO3df4xlZ13H8ffHri0WQrpLZzdrt7rFbISlUWkmtUhCSFZCRcIWkyZLgk6wZkNSsBiNbuWP8g+miBo1EZKVVlYlJU0t6YYEZF3F6h8tTmkp227LFhrbpevuYCOoJMDi1z/uabhMZ3Zm7rl3fjzzfiU355znnDP3+9yz+5lnzr3n3FQVkqS2/MhaFyBJGj/DXZIaZLhLUoMMd0lqkOEuSQ3astYFAFx++eW1e/futS5DkjaUhx566BtVNbXQunUR7rt372Z2dnaty5CkDSXJvy+2ztMyktQgw12SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1aMlwT3JnknNJTgy1fTjJE0keTfKpJJcNrbs1yVNJnkzy5kkVLkla3HJG7h8Hrp/Xdgy4uqp+BvgKcCtAkr3AAeA13T4fSXLR2KqVJC3LkuFeVfcDz89r+1xVne8WHwB2dfP7gU9W1Xeq6mngKeDaMdYrjSQZ/SFtROM45/7rwGe6+SuAZ4fWne7aXiTJwSSzSWbn5ubGUIYk6QW9wj3J+4HzwCdeaFpgswW/x6+qDlfVdFVNT00teN8bSdKIRr5xWJIZ4K3AvvrBF7GeBq4c2mwX8Nzo5UmSRjHSyD3J9cDvAW+rqm8PrToKHEhySZKrgD3AF/qXKUlaiSVH7knuAt4IXJ7kNHAbg0/HXAIcy+Adpweq6t1V9ViSu4HHGZyuubmqvj+p4iVJC8sPzqisnenp6fJ+7pqkPp96WQf/RaQFJXmoqqYXWrcuvqxDWs/8xaCNyNsPSFKDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSg7xCVavGKz2l1ePIXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktQgw12SGuQVqtI65RW96sORuyQ1yHCXpAYZ7pLUIMNdkhq0ZLgnuTPJuSQnhtq2JTmW5FQ33Tq07tYkTyV5MsmbJ1W4JGlxyxm5fxy4fl7bIeB4Ve0BjnfLJNkLHABe0+3zkSQXja1aSdKyLBnuVXU/8Py85v3AkW7+CHDDUPsnq+o7VfU08BRw7Zhq1SaW9HtIm82o59x3VNUZgG66vWu/Anh2aLvTXZskaRWN+w3VhcZIC15OkeRgktkks3Nzc2MuQ5I2t1HD/WySnQDd9FzXfhq4cmi7XcBzC/2AqjpcVdNVNT01NTViGZKkhYwa7keBmW5+BrhvqP1AkkuSXAXsAb7Qr0RJ0koteW+ZJHcBbwQuT3IauA24Hbg7yU3AM8CNAFX1WJK7gceB88DNVfX9CdUuSVrEkuFeVe9YZNW+Rbb/IPDBPkVJkvrxClVJapDhLkkN8n7u0gR5AZXWiiN3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktQgw12SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhrUK9yT/FaSx5KcSHJXkpck2ZbkWJJT3XTruIqVJC3PyOGe5ArgN4HpqroauAg4ABwCjlfVHuB4tyxJWkV9T8tsAX4syRbgUuA5YD9wpFt/BLih53NIklZo5HCvqq8DfwQ8A5wBvllVnwN2VNWZbpszwPaF9k9yMMlsktm5ublRy5AkLaDPaZmtDEbpVwE/Drw0yTuXu39VHa6q6aqanpqaGrUMSdIC+pyW+UXg6aqaq6rvAfcCvwCcTbIToJue61+mJGkl+oT7M8B1SS5NEmAfcBI4Csx028wA9/UrUZK0UltG3bGqHkxyD/BF4DzwMHAYeBlwd5KbGPwCuHEchUqSlm/kcAeoqtuA2+Y1f4fBKF6StEa8QlWSGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUoF4XMWnzSda6AknL4chdkhpkuEtSgwx3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktSgXuGe5LIk9yR5IsnJJK9Lsi3JsSSnuunWcRUrSVqeviP3PwM+W1WvAn4WOAkcAo5X1R7geLcsSVpFI4d7kpcDbwDuAKiq71bVfwH7gSPdZkeAG/oWKUlamT4j91cCc8BfJXk4yceSvBTYUVVnALrp9oV2TnIwyWyS2bm5uR5lSJLm6xPuW4BrgI9W1WuB/2UFp2Cq6nBVTVfV9NTUVI8yJEnz9Qn308DpqnqwW76HQdifTbIToJue61eiJGmlRg73qvoP4NkkP9017QMeB44CM13bDHBfrwolSSu2pef+7wU+keRi4GvAuxj8wrg7yU3AM8CNPZ9D0gol/favGk8dWju9wr2qHgGmF1i1r8/PlST14xWqktQgw12SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSgwx3SWqQ4S5JDTLcJalBhrskNchwl6QG9f0mJkkN6vNNTn6L0/rgyF2SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUoN7hnuSiJA8n+XS3vC3JsSSnuunW/mVKklZiHCP3W4CTQ8uHgONVtQc43i1LklZRr3BPsgv4ZeBjQ837gSPd/BHghj7PIUlaub4j9z8Ffhf4v6G2HVV1BqCbbl9oxyQHk8wmmZ2bm+tZhiRp2MjhnuStwLmqemiU/avqcFVNV9X01NTUqGVIkhbQ566QrwfeluQtwEuAlyf5W+Bskp1VdSbJTuDcOArV+PS545+kjWHkkXtV3VpVu6pqN3AA+MeqeidwFJjpNpsB7utdpSRpRSbxOffbgTclOQW8qVuWJK2isXxZR1V9Hvh8N/+fwL5x/FxJ0mi8QlWSGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1aCwfhdTq8gpTrWd9/n1Wja+Ozc6RuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSgwx3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBo0c7kmuTPJPSU4meSzJLV37tiTHkpzqplvHV64kaTn6jNzPA79dVa8GrgNuTrIXOAQcr6o9wPFuWZK0ikYO96o6U1Vf7Ob/GzgJXAHsB450mx0BbuhbpCRpZcZyzj3JbuC1wIPAjqo6A4NfAMD2cTyHJGn5eod7kpcBfwe8r6q+tYL9DiaZTTI7NzfXtwxJ0pBe4Z7kRxkE+yeq6t6u+WySnd36ncC5hfatqsNVNV1V01NTU33KkCTN0+fTMgHuAE5W1Z8MrToKzHTzM8B9o5cnSRrFlh77vh74VeDLSR7p2n4fuB24O8lNwDPAjf1KlCSt1MjhXlX/CmSR1ftG/bmSpP68QlWSGtTntMyml8X+blmGqvHVIUnzOXKXpAY5cl8jfUb9krQUR+6S1CBH7pLWDd/HGh9H7pLUIMNdkhpkuEtSgwx3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ3yrpCSNr2+36+wHu9IuenD3S/NkNQiT8tIUoOaGLk7+pZkDvwwR+6S1CDDXZIa1MRpGUlaS+vxu18nNnJPcn2SJ5M8leTQpJ5HkvRiEwn3JBcBfwH8ErAXeEeSvZN4LknSi01q5H4t8FRVfa2qvgt8Etg/oeeSJM0zqXPuVwDPDi2fBn5+eIMkB4GD3eL/JHlyQrWMy+XAN9a6iDWymfsOm7v/m7nvsAr97/kRzp9cbMWkwn2hcn/obYOqOgwcntDzj12S2aqaXus61sJm7jts7v5v5r7Dxu7/pE7LnAauHFreBTw3oeeSJM0zqXD/N2BPkquSXAwcAI5O6LkkSfNM5LRMVZ1P8h7g74GLgDur6rFJPNcq2jCnkCZgM/cdNnf/N3PfYQP3P7Ue71UpSerF2w9IUoMMd0lq0KYM96VujZBka5JPJXk0yReSXD207pYkJ5I8luR9Q+3bkhxLcqqbbl2t/qzUhPr/gSRfT/JI93jLavVnuZLcmeRckhOLrE+SP+9el0eTXDO0bsHXbIMd90n0f90fd+jd9wX3XffHvqo21YPBG7xfBV4JXAx8Cdg7b5sPA7d1868CjnfzVwMngEsZvBn9D8Cebt0fAoe6+UPAh9a6r6vc/w8Av7PW/Vui728ArgFOLLL+LcBnGFyncR3w4FKv2UY57hPs/7o/7n36fqF91/ux34wj9+XcGmEvcBygqp4AdifZAbwaeKCqvl1V54F/Bt7e7bMfONLNHwFumGw3Rjap/q97VXU/8PwFNtkP/HUNPABclmQnF37NNspxn1T/N4Qefb/Qvuv62G/GcF/o1ghXzNvmS8CvACS5lsElvrsYjFrfkOQVSS5l8Nv+hYu1dlTVGYBuun1iPehnUv0HeE/3J+2d6+5P1OVZ7LW50Gu2UY77cozSf9j4xx2W9/9ivnV97DdjuC95awTgdmBrkkeA9wIPA+er6iTwIeAY8FkGIXh+grVOwqT6/1Hgp4CfA84Afzz+0idusddmOa9ZC0bpfwvHHRo8xpvxyzqWvDVCVX0LeBcM3mgBnu4eVNUdwB3duj/ofh7A2SQ7q+pM9+fcuUl2ooeJ9L+qzr6wf5K/BD49sR5MzmKvzcWLtMPGOe7LseL+N3LcYbRbpqzrY78ZR+5L3hohyWXdOoDfAO7vAo8k27vpTzA4dXFXt91RYKabnwHum2gvRjeR/r9wfrLzdgancDaao8CvdZ+cuA74Zvfn9oVes41y3Jdjxf1v5LjD4n1fap/1e+zX+h3dtXgwOFf8FQafAHh/1/Zu4N3d/OuAU8ATwL3A1qF9/wV4nMEpiX1D7a9g8CbkqW66ba37ucr9/xvgy8CjDP7R71zrfi7Q77sYnDr4HoOR2k3z+h0GXzLz1a4v0xd6zTbgcZ9E/9f9cR9D31+070Y49t5+QJIatBlPy0hS8wx3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1KD/By9l0OMdVhrYAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAQp0lEQVR4nO3df4xlZX3H8fenoChawyKzZGVpF81GRdJWMqGojTHZEgg1LjYhWRPbraXZmKBFU2OhJsX+QYPV2h9JNdkKdasEslEMGxOt262WNing8NOFBXeVFBZWdiyp2pqga7/94x7S2/Hu7sw9986Px/crmdxznnPO3O+zz+xnnjn33nNSVUiS2vJzK12AJGnyDHdJapDhLkkNMtwlqUGGuyQ16NSVLgDgrLPOqk2bNq10GZK0ptx7773fraqZUdtWRbhv2rSJubm5lS5DktaUJP9+vG2elpGkBp003JPcnORokv0jtn0gSSU5a6jtuiSHkjyW5NJJFyxJOrnFzNw/DVy2sDHJucAlwBNDbecD24DXdcd8IskpE6lUkrRoJw33qroTeHbEpr8APggMX79gK3BbVT1XVY8Dh4CLJlGoJGnxxjrnnuRtwFNV9eCCTecATw6tH+7aRn2PHUnmkszNz8+PU4Yk6TiWHO5JTgc+BPzxqM0j2kZemayqdlbVbFXNzsyMfCePJGlM47wV8lXAecCDSQA2AvcluYjBTP3coX03Ak/3LVKStDRLnrlX1Teqan1VbaqqTQwC/cKq+g6wB9iW5LQk5wGbgXsmWrEk6aQW81bIW4F/A16d5HCSq463b1U9DOwGHgG+DFxdVT+ZVLGSpMU56WmZqnrHSbZvWrB+A3BDv7KkVSSjXkpaJG+GoxXiJ1QlqUGGuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSgwx3SWqQ4S5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQePcIFtae/rcTUlag5y5S1KDDHdJapDhLkkNOmm4J7k5ydEk+4faPprk0SQPJflCkjOGtl2X5FCSx5JcOq3CJUnHt5iZ+6eByxa07QUuqKpfAr4JXAeQ5HxgG/C67phPJDllYtVKkhblpOFeVXcCzy5o+0pVHetW7wI2dstbgduq6rmqehw4BFw0wXolSYswiXPuvwt8qVs+B3hyaNvhru2nJNmRZC7J3Pz8/ATKkCQ9r1e4J/kQcAy45fmmEbvVqGOramdVzVbV7MzMTJ8yJEkLjP0hpiTbgbcCW6rq+QA/DJw7tNtG4Onxy5MkjWOsmXuSy4A/BN5WVT8c2rQH2JbktCTnAZuBe/qXKUlaipPO3JPcCrwFOCvJYeB6Bu+OOQ3Ym8HHuu+qqndX1cNJdgOPMDhdc3VV/WRaxUuSRsv/nVFZObOzszU3N7fSZahlK3VtmVXw/0vtSnJvVc2O2uYnVCWpQYa7JDXIcJekBhnuktQgw12SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSg8a+E5O0rFbqkr3SGuXMXZIaZLhLUoMMd0lqkOEuSQ3yBVVpmvq8EOz9V9WDM3dJatBJwz3JzUmOJtk/1HZmkr1JDnaP64a2XZfkUJLHklw6rcIlSce3mJn7p4HLFrRdC+yrqs3Avm6dJOcD24DXdcd8IskpE6tWkrQoJw33qroTeHZB81ZgV7e8C7hiqP22qnquqh4HDgEXTahWSdIijXvO/eyqOgLQPa7v2s8Bnhza73DXJklaRpN+QXXUWwNGvuSfZEeSuSRz8/PzEy5Dkn62jRvuzyTZANA9Hu3aDwPnDu23EXh61Deoqp1VNVtVszMzM2OWIUkaZdxw3wNs75a3A3cMtW9LclqS84DNwD39SpQkLdVJP8SU5FbgLcBZSQ4D1wM3AruTXAU8AVwJUFUPJ9kNPAIcA66uqp9MqXZJ0nGcNNyr6h3H2bTlOPvfANzQpyhJUj9+QlWSGmS4S1KDDHdJapDhLkkN8pK/Wj7eB1VaNs7cJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoN8K6TUor5vO62Rt2HQGuLMXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktQgw12SGmS4S1KDeoV7kvcneTjJ/iS3JnlRkjOT7E1ysHtcN6liJUmLM3a4JzkH+H1gtqouAE4BtgHXAvuqajOwr1uXJC2jvqdlTgVenORU4HTgaWArsKvbvgu4oudzSJKWaOxwr6qngI8BTwBHgO9V1VeAs6vqSLfPEWD9qOOT7Egyl2Rufn5+3DIkSSP0OS2zjsEs/TzgFcBLkrxzscdX1c6qmq2q2ZmZmXHLkCSN0Oe0zK8Dj1fVfFX9GLgdeCPwTJINAN3j0f5lSpKWok+4PwFcnOT0JAG2AAeAPcD2bp/twB39SpQkLdXYd2KqqruTfA64DzgG3A/sBF4K7E5yFYNfAFdOolBJ0uL1us1eVV0PXL+g+TkGs3hJ0grxE6qS1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDWo1/vcJU1RstIVaA1z5i5JDTLcJalBhrskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktSgXuGe5Iwkn0vyaJIDSd6Q5Mwke5Mc7B7XTapYSdLi9J25/xXw5ap6DfDLwAHgWmBfVW0G9nXrkqRlNHa4J3kZ8GbgJoCq+lFV/SewFdjV7bYLuKJvkZKkpekzc38lMA/8XZL7k3wqyUuAs6vqCED3uH7UwUl2JJlLMjc/P9+jDEnSQn3C/VTgQuCTVfV64L9ZwimYqtpZVbNVNTszM9OjDEnSQn3C/TBwuKru7tY/xyDsn0myAaB7PNqvREnSUo0d7lX1HeDJJK/umrYAjwB7gO1d23bgjl4VSpKWrO8Nst8L3JLkhcC3gXcx+IWxO8lVwBPAlT2fQ5K0RL3CvaoeAGZHbNrS5/tKkvrpO3PXz5pkpSuQtAhefkCSGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhpkuEtSgwx3SWqQ15b5WeT1YaTmOXOXpAYZ7pLUIMNdkhpkuEtSg3xBVdJP6/Oie9Xk6tDYnLlLUoMMd0lqkOEuSQ3qHe5JTklyf5IvdutnJtmb5GD3uK5/mZKkpZjEzP0a4MDQ+rXAvqraDOzr1iVJy6hXuCfZCPwG8Kmh5q3Arm55F3BFn+eQJC1d35n7XwIfBP5nqO3sqjoC0D2uH3Vgkh1J5pLMzc/P9yxDkjRs7HBP8lbgaFXdO87xVbWzqmaranZmZmbcMiRJI/T5ENObgLcluRx4EfCyJJ8FnkmyoaqOJNkAHJ1EoZKkxRt75l5V11XVxqraBGwD/qmq3gnsAbZ3u20H7uhdpSRpSabxPvcbgUuSHAQu6dYlSctoIteWqaqvAV/rlv8D2DKJ7ytJGo+fUJWkBhnuktQgL/m7FnkPVEkn4cxdkhpkuEtSgwx3SWqQ4S5JDTLcJalBvltG0mR5c+1VwZm7JDXIcJekBhnuktQgw12SGmS4S1KDDHdJapDhLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkho0drgnOTfJV5McSPJwkmu69jOT7E1ysHtcN7lyJUmL0Wfmfgz4g6p6LXAxcHWS84FrgX1VtRnY161LkpbR2OFeVUeq6r5u+QfAAeAcYCuwq9ttF3BF3yIlSUszkXPuSTYBrwfuBs6uqiMw+AUArJ/Ec0iSFq93uCd5KfB54H1V9f0lHLcjyVySufn5+b5lSJKG9Ar3JC9gEOy3VNXtXfMzSTZ02zcAR0cdW1U7q2q2qmZnZmb6lCFJWqDPu2UC3AQcqKqPD23aA2zvlrcDd4xfniRpHH3uofom4LeAbyR5oGv7I+BGYHeSq4AngCv7ldioPveZlKSTGDvcq+pfgeMl1JZxv68kqT8/oSpJDTLcJalBfc65S9Jk9XktqmpydTTAmbskNchwl6QGGe6S1CDDXZIaZLhLUoMMd0lqkOEuSQ0y3CWpQYa7JDXIcJekBhnuktQgry3Th9dkl7RKOXOXpAY5c5ekvn+Fr8IrUjpzl6QGGe6S1CDDXZIa5Dl33/Eiqa9VeAcpZ+6S1KCphXuSy5I8luRQkmun9TySBAxmz+N+NWgqp2WSnAL8DXAJcBj4epI9VfXINJ6v1cGRpHFNa+Z+EXCoqr5dVT8CbgO2Tum5JEkLTOsF1XOAJ4fWDwO/OrxDkh3Ajm71v5I81uP5zgK+2+P41aKVfoB9WY1a6Qe01JekT19+8XgbphXuo86T/L+XhKtqJ7BzIk+WzFXV7CS+10pqpR9gX1ajVvoB9mUxpnVa5jBw7tD6RuDpKT2XJGmBaYX714HNSc5L8kJgG7BnSs8lSVpgKqdlqupYkvcA/wCcAtxcVQ9P47k6Ezm9swq00g+wL6tRK/0A+3JSqVV4NTNJUj9+QlWSGmS4S1KDVlW4J7k5ydEk+4+zfV2SLyR5KMk9SS4Y2vb+JA8n2Z/k1iQv6to/nOSpJA90X5evgb5c0/Xj4STvG2o/M8neJAe7x3VruC/LPi5Jzk3y1SQHunquGbFPkvx1d9mMh5JcOLRt5CU1VmJcptiXtTguI38+l3tcptiP8cakqlbNF/Bm4EJg/3G2fxS4vlt+DbCvWz4HeBx4cbe+G/idbvnDwAfWUF8uAPYDpzN4wfsfgc3dtj8Dru2WrwU+sob7suzjAmwALuyWfx74JnD+gn0uB77E4LMaFwN3d+2nAN8CXgm8EHjw+WNXYlym2Jc1NS4n+vlc7nGZYj/GGpNVNXOvqjuBZ0+wy/nAvm7fR4FNSc7utp0KvDjJqQzCZEXfV9+jL68F7qqqH1bVMeCfgbd3x2wFdnXLu4ArplH7QlPqy7KrqiNVdV+3/APgAIOJwbCtwN/XwF3AGUk2cOJLaiz7uEyxL8uuZ19O9PO5rOMyxX6MZVWF+yI8CPwmQJKLGHz0dmNVPQV8DHgCOAJ8r6q+MnTce7o/gW5erlMZizCyLwxmum9O8vIkpzP4Tf/8B8LOrqojMPhBAtYve9WjjdMXWMFxSbIJeD1w94JNoy6dcc4J2mGFx2XCfYG1NS4nsmLjMuF+wBhjstbC/UZgXZIHgPcC9wPHus5uBc4DXgG8JMk7u2M+CbwK+BUGwf/ny171aCP7UlUHgI8Ae4EvMwjOYytW5eKM05cVG5ckLwU+D7yvqr6/cPOIQ+oE7StqCn1Za+Oy6kyhH2ONyZq6E1P3D/UuGLwwweA8++PApcDjVTXfbbsdeCPw2ap65vnjk/wt8MXlrnuUE/SFqroJuKnb9qcMfrsDPJNkQ1Ud6f6UO7rshY8wTl9WalySvIDBf7xbqur2Ebsc79IZLzxOO6zQuEyjL2twXE5k2cdlGv0Yd0zW1Mw9yRkZXM4A4PeAO7tgeQK4OMnpXbhsYXC+i+fPZ3XezuBUwYo7QV9Isr57/AUGpztu7fbbA2zvlrcDdyxfxcc3Tl9WYly6n42bgANV9fHj7LYH+O3uXQ0XMzjFd4QTX1Jj2cdlWn1Zg+NyIss6LtPqx9hj0vcV4kl+MfiPfwT4MYPfcFcB7wbe3W1/A3AQeBS4HVg3dOyfdO37gc8Ap3XtnwG+ATzU/cNuWAN9+RfgEQanMbYMtb+cwQuXB7vHM9dwX5Z9XIBfY/An8EPAA93X5Qv6EgY3mvlWV9/s0PGXM3gHxLeAD63kuEyxL2txXH7q53MlxmWK/RhrTLz8gCQ1aE2dlpEkLY7hLkkNMtwlqUGGuyQ1yHCXpAYZ7pLUIMNdkhr0v9GsQJg7KNA7AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3de5zN5d7/8ddnxiGkwiCRQ+1BhNgjpzZFEnfb4W53Zxdmy2nKodztbMVOtg5C3XRyymF+NZGUjGoXpo1EGOfzoSKHwUQOoWmM6/fHtdSUwcw6XevweT4e67HW+q61rLfvzPrMta7v9b0uMcaglFIqssS4DqCUUsr/tLgrpVQE0uKulFIRSIu7UkpFIC3uSikVgQq5DgAQFxdnqlat6jqGUkqFldWrV39vjCmb12MhUdyrVq1Kenq66xhKKRVWRGTPxR7TbhmllIpAWtyVUioCaXFXSqkIpMVdKaUikBZ3pZSKQJct7iIyVUQOi8imXNtGi8g2EdkgInNE5Jpcjz0pIrtEZLuItAlUcKWUUheXn5b7dODu321bANxsjKkL7ACeBBCRWkBnoLbnNW+ISKzf0iqllMqXyxZ3Y8wS4Ojvts03xpz13P0KqOS53QGYaYzJMsZ8C+wCbvVjXqWCwhjDmTNnfrmcPXv28i9SKoT44ySmh4B3PbcrYov9efs82y4gIr2B3gCVK1f2QwylfHPw4EFSU1NZvHgxixYt4sCBA788VqRIERo1asTtt99O69atue222xARh2mVujSfiruIDAHOAinnN+XxtDxXAzHGTAImASQkJOiKIcqZffv28eKLLzJ58mSysrK49tpradGiBXXr1iU21vYqHj58mCVLlvDcc88xYsQIEhISePrpp7nnnnu0yKuQ5HVxF5FE4B6glfl1Oad9wPW5nlYJOPD71yoVbDI8jwJ8FvgcWIFtgtQDmsDBsgd5V97l3ex3IRvMsF/bHsePH2fWrFm88MILtG/fngYNGjB16lTq1asXnP+IUvnk1VBIEbkb+AfQ3hhzOtdDqUBnESkqItWAeGCl7zGV8rMjwBRgGVAH6I89YlSOvL9/elx99dX06tWL7du3M336dDIyMmjUqBGvvfYaumSlCiX5GQo5A1gO1BCRfSLSA3gNKAksEJF1IjIBwBizGZgFbAE+BfoaY3ICll4pb2wEJgLHsGO7OgKlCvZPFC5cmMTERNavX0+rVq3o378/nTp14sSJE36Pq5Q3JBRaGwkJCUZnhVSB9Eu3zErgE6AycC9w9eVfm7tbJs/HjWHcuHE88cQTNGjQgE8//ZRSpQr410IpL4jIamNMQl6P6RmqKnoswxb26kBX8lXY80NEeOyxx3j//fdZt24dLVu2JDMz0z//uFJe0uKuosMSYD5QC7gfKOz/t2jfvj3z5s1j+/bttGjRgkOHDvn/TZTKJy3uKuIlJyfbUTF1sV0xATxn+q677uLf//43u3fvpkOHDpw5cyZwb6bUJWhxVxHtiy++oFevXlANOxomCJNhtGjRgpSUFFasWEH37t11FI1yQou7ili7du2iU6dO3HDDDfA/BKWwn9epUydGjhzJu+++yzPPPBO8N1bKQ4u7ikinTp2iffv2GGP46KOPoFjwMwwaNIju3bvzr3/9i9mzZwc/gIpqWtxVRHrsscfYtm0bs2bN4g9/+IOTDCLChAkTaNSoET179mTPnouuZayU32lxVxFn9uzZvPnmmwwePJhWrVo5zVKkSBHeeecdzp07R5cuXcjJ0XP6VHDoSUwqLOQ5N0xejgPjgdJAD4Laz56X8ydAvf3223Tt2pV//etf/POf/3QbSkUMPYlJRYdzwAee6wAPeSyoLl268OCDDzJ8+HCWL1/uOo6KAlrcVeRYDewB2gJlHGfJwxtvvEHFihXp0aMHWVlZruOoCKfFXUWGE8BC7Hj2WxxnuYirrrqK8ePHs3XrVkaNGuU6jopwWtxVZPg3kINdYSCE185o164d999/P88++yzbt293HUdFMC3uKvxtA7YCLQjJ7pjfGzt2LMWLF6dPnz569qoKGC3uKrxlYWd6LAc0dZwln6699lpGjx7N4sWLmT59uus4KkJpcVfh7Utsf/s9hNTomMt56KGHaNq0KU8++SQnT550HUdFIC3uKnwdw87RfjN28Y0wEhMTw9ixYzl06BAvvPCC6zgqAmlxV+Froef6TqcpvNawYUO6dOnCyy+/zO7du13HURFGi7sKT3uBTdh+9mscZ/HBCy+8QExMDP/4xz9cR1ERRou7Cj8Gu/z6lUAzx1l8VKlSJQYNGsSsWbP48ssvXcdREUSLuwo/W4D9QCugqOMsfvDEE09QsWJFnnjiCR0aqfxGi7sKL+eA/wBxQD3HWfykRIkSPP300yxfvpxPPvnEdRwVIbS4q/CyAfgeaElE/fZ2796dG2+8kSFDhnDu3DnXcVQEiKCPh4p4Z4FFQAXgJrdR/K1w4cIMHz6c9evX66pNyi+0uKvwsQY7tr0lIT1/jLc6d+5M7dq1+ec//8nZs2ddx1FhTou7Cg8/A0uwJyu5WTUv4GJjYxkxYgQ7duzgrbfech1HhbnLFncRmSoih0VkU65tpUVkgYjs9FyXyvXYkyKyS0S2i0ibQAVXUWYN8CMR22o/r2PHjiQkJDBixAhtvSufFMrHc6YDrwH/L9e2wUCaMWakiAz23P+HiNQCOgO1geuAhSJS3RijC0cqr2VlZdk5ZKoAVR2HKaB8Lw+Y2x+AdJgxYwZdu3b1eyYVHS7bcjfGLAGO/m5zByDZczsZ6Jhr+0xjTJYx5ltgF3Crn7KqKDVt2jQ4CTR3nSRIagDl4bnnntMFtZXXvO1zL2+MyQDwXJfzbK+IPTH8vH2ebRcQkd4iki4i6ZmZmV7GUJEuOzubkSNHQiXgBtdpgkSA5rB9+3bef/9912lUmPL3AdW8voPmecqdMWaSMSbBGJNQtmxZP8dQkSIlJYU9e/bYVnsE97Vf4CaoWbMmzz77rI57V17xtrgfEpEKAJ7rw57t+4Drcz2vEnDA+3gqmuXk5PD8889Tv359iHedJshi4KmnnmLjxo3MmzfPdRoVhrwt7qlAoud2IjA31/bOIlJURKphP5IrfYuootXs2bPZuXMnQ4YMia5Wu8df//pXbrzxRp577jmdc0YVWH6GQs4AlgM1RGSfiPQARgKtRWQn0NpzH2PMZmAWdmqnT4G+OlJGecMYw6hRo6hRowadOnVyHceJQoUK8fe//51Vq1axZMkS13FUmMnPaJm/GmMqGGMKG2MqGWOmGGOOGGNaGWPiPddHcz3/OWPMjcaYGsaYfwc2vopU//nPf1izZg2PP/44MTHRe65dYmIiZcuWZfTo0a6jqDATvZ8aFdJGjRpF+fLlo36cd7Fixejfvz8ff/wxmzdvdh1HhREt7irkbNiwgc8++4wBAwZwxRVXuI7j3COPPELx4sUZM2aM6ygqjGhxVyFn9OjRlChRgocffth1lJBQpkwZevToQUpKCvv373cdR4UJLe4qpOzdu5eZM2fSq1cvSpUqdfkXRImBAweSk5PD2LFjXUdRYUKLuwopr7zyCsYYHnvsMddRQkq1atW47777mDRpEidPnnQdR4UBLe4qZPz4449MnjyZe++9lypVqriOE3IGDhzIiRMnmD59uusoKgzkZ1ZIpYJi+vTpHD9+nIEDB7qOEhLynFGyEgwYNoAB3w+4ZNPMDNOTnqKdttxVSDh37hzjxo2jUaNGNG7c2HWc0NUY+AHY4TqICnVa3FVI+OSTT9i1a5f2tV/OTcBVwFeug6hQp8VdhYT/+7//o1KlStx7772uo4S2WKARsBvIcBtFhTYt7sq5DRs28Pnnn9OvXz8KFy7sOk7oawAUBla4DqJCmRZ35dwrr7xCsWLF6NWrl+so4aEYcAuwEbuurFJ50OKunDp69CgpKSl06dKF0qVLu44TPm4FcrALhyuVBx0KqYImz6F9XwI/weTYyUwePjnomcJWWeyyg+lAM2xfvFK5aMtduXMOWAVUAa51nCUc3QqcALa7DqJCkRZ35c4O4Bi2SKmCqw5cgx5YVXnS4q7cWQmUBGq6DhKmYoCGwB7goOMsKuRocVduZALfYIuT9hd7rz72yNkq10FUqNHirtxYhS3qDVwHCXPFgTrAeuCM4ywqpGhxV8GXBawDagFXOs4SCRoCZ7EFXikPLe4q+DYCP2OLkvLddUAl7LchnQxSeWhxV8FlsEWoPHC94yyRpCFwBPjWdRAVKrS4q+DaCxzCFqM8zmlSXqqFnZZAD6wqDy3uKrhWAUWxBwGV/xTGjpzZhj2xSUU9Le4qeH4EtgD1sAVe+VcCtttrtesgKhT4VNxFZKCIbBaRTSIyQ0SuEJHSIrJARHZ6rnUJe2WtxU52pQdSA6M08AdgNWRnZ7tOoxzzuriLSEVgAJBgjLkZO2q5MzAYSDPGxANpnvsqyp07d862KKtiJ71SgdEQ+BHmzZvnOolyzNdumUJAMREphD2d4gDQAUj2PJ4MdPTxPVQEmD9/vp1HJsF1kggXD1wF48ePd51EOeZ1cTfG7AfGAN9hF/w6boyZD5Q3xmR4npMBlMvr9SLSW0TSRSQ9MzPT2xgqTIwfPx5KoPPIBFoM8EdYuHAhO3fudJ1GOeRLt0wpbCu9GvY0ihIi0iW/rzfGTDLGJBhjEsqW1e/pkWzv3r189NFHv86DogKrARQqVIhJkya5TqIc8qVb5k7gW2NMpjEmG/gAaAocEpEKAJ7rw77HVOHszTffxBgDf3SdJEqUhI4dOzJt2jR++ukn12mUI74U9++AxiJSXEQEaAVsBVKBRM9zEoG5vkVU4Sw7O5vJkyfTtm1b0HFTQZOUlMSRI0eYPXu26yjKEV/63FcAs7GrOG70/FuTgJFAaxHZCbT23FdR6qOPPiIjI4OkpCTXUaLKHXfcQXx8PBMmTHAdRTni02gZY8wwY0xNY8zNxpiuxpgsY8wRY0wrY0y85/qov8Kq8DNhwgSuv/562rVr5zpKVImJiSEpKYkvv/ySjRs3uo6jHNAzVFXAfP3118yfP59evXoRG6srcgRbYmIiRYsWZeLEia6jKAe0uKuAmTRpErGxsfTo0cN1lKhUpkwZ7rvvPt566y1OnTrlOo4KMi3uKiCysrKYNm0a7du357rrrnMdJ2r16dOHEydOMHPmTNdRVJBpcVcBMWfOHDIzM/VAqmPNmjWjdu3a2jUThbS4q4CYOHEiN9xwA3feeafrKFFNROjTpw+rVq1izZo1ruOoINLirvxu27ZtLFq0iN69exMTo79irnXt2pVixYpp6z3K6CdP+d3EiRMpXLgw3bt3dx1FAddccw2dO3cmJSWFEyd0JY9oocVd+dWZM2dITk6mU6dOlCuX55xxyoGkpCROnTrFO++84zqKChIt7sqv3nvvPX744Qc9kBpiGjZsyC233MKECRPsPD8q4mlxV341ceJEqlevzu233+46ispFREhKSmL9+vWsXLnSdRwVBFrcld9s3LiRZcuW0adPH+xcciqUPPDAA1x55ZV6YDVKaHFXfjNx4kSKFi1KYmLi5Z+sgq5kyZI88MADzJw5k2PHjrmOowJMl05QBSLDL9Ii/xl4E6gBca/FBTOSKoCkpCQmTZrEW2+9Rf/+/V3HUQGkLXflH5uALHRBjhBXv359GjZsqAdWo4AWd+Uf6UBZoLLrIOpykpKS2LJlC0uXLnUdRQWQFnfluwOeSwKgx1FD3v3338/VV1+tC3lEOC3uynfpQGGgnusgKj9KlChBt27dmD17NpmZma7jqADR4q588xN2kcWbgSscZ1H51qdPH37++WemT5/uOooKEC3uyjcbgGxsl4wKG7Vr1+ZPf/oTEydO5Ny5c67jqADQ4q68Z7BdMhWAio6zqAJLSkri66+/ZuHCha6jqADQ4q68txc4DDR0HUR549577yUuLk4PrEYoLe7Ke+lAUWx/uwo7RYsW5aGHHiI1NZX9+/e7jqP8TIu78s4pYDNQFyjiOIvyWu/evcnJyWHy5Mmuoyg/0+KuvLMWyEG7ZMLcjTfeSJs2bZg8eTLZ2dmu4yg/0rllVMGdw3bJVAF0PY6QdNE5gPJSCjgARboUgVp2kxmmUxOEO225q4L7GjiGttojRXXgamCV6yDKn3wq7iJyjYjMFpFtIrJVRJqISGkRWSAiOz3XpfwVVoWIVUAJoKbrIMovYrATvn0LfO84i/IbX1vu44BPjTE1sSefbwUGA2nGmHggzXNfRYofgB3YYqCdepGjAbYapLsOovzF6+IuIlcBzYEpAMaYn40xx4AOQLLnaclAR19DqhCyGjs5mE7tG1muxPa3r8POza/Cni8t9xuATGCaiKwVkTdFpARQ3hiTAeC5zvOQm4j0FpF0EUnXyYvCQ1ZWFqzh1z5aFVkaYucK2uQ6iPIHX4p7IeyXufHGmPrYkc/57oIxxkwyxiQYYxLKli3rQwwVLO+99x6cRg+kRqrK2KbYSnQhjwjgS3HfB+wzxqzw3J+NLfaHRKQCgOf6sG8RVah47bXXoAz2O5uKPALcChyE5cuXu06jfOR1cTfGHAT2ikgNz6ZWwBYgFTi/QnIiMNenhCokpKens2LFCttq1wG0kasOUNTzh1yFNV8/pv2BFBHZANwCPA+MBFqLyE6gtee+CnOvv/46JUqUsD9lFbmKAvVtF1xGRobrNMoHPhV3Y8w6T795XWNMR2PMD8aYI8aYVsaYeM/1UX+FVW58//33zJgxg27duumCHNGgIZw9e1bnmwlz+gVbXdaUKVPIysqib9++rqOoYCgDbdq0YcKECTrfTBjT4q4uKScnh/Hjx3PHHXdQu3Zt13FUkPTr14+MjAzmzJnjOorykhZ3dUmpqans2bOHfv36uY6igqht27bccMMNvPLKK66jKC9pcVeXNG7cOKpUqUL79u1dR1FBFBsbS79+/fjyyy9ZvXq16zjKCzo7iLqodevWsXjxYkaNGkWhQvqrEk1kuMAZoDAk9EiATvl/rU4XHBq05a4u6tVXX6V48eL07NnTdRTlQjHs0NdNwI+Os6gC0+Ku8pSZmUlKSgrdunWjVCmdtTlqNcKuuKWzRYYdLe4qT5MmTSIrK4v+/fu7jqJcigP+gJ3D/6zjLKpAtLirC2RnZ/PGG2/QunVratWq5TqOcq0Rvy6IrsKGFnd1gVmzZnHgwAEeffRR11FUKLgR24L/CtBjpWFDi7v6DWMML7/8MjVr1qRt27au46hQEAM0BjKAPY6zqHzT4q5+Y8mSJaxZs4aBAwcSE6O/HsqjHnb0jM4EHDb006t+46WXXiIuLo6uXbu6jqJCSWHsdM/b0UW0w4QWd/WLHTt2MG/ePB555BGKFSvmOo4KNbcCsdi+dxXytLirX4wdO5YiRYrwyCOPuI6iQtGVQF3sItqnHWdRl6XFXQFw5MgRpk+fTpcuXShfvrzrOCpUNcaOd9eTmkKeThgShWS4XLhxEXAGpl4xlanDpwY7kgoX5bEnNa0AmmD74lVI0pa7gp+xH9bqQDnHWVToa4Y9qWmd6yDqUrS4K1iLnQGwmesgKixUBa4DlgHn3EZRF6fFPdrlYMcuXw9UcZxFhQcBbgN+ALY6zqIuSvvco90W4Bhwt+sgKqzUBEoDS4Fa2ILvkecxnXzSueD9R1vu0cxgP5xx2P52pfIrBtuNlwF86ziLypMW92i2CziE/ZDqb4IqqHrYse9LXQdRedGPdLQywBLgKqCO4ywqPBXCDof8BtjnOIu6gBb3aLUb2Is9MKZHXpS3ErATii1xHUT9ns/FXURiRWStiHzkuV9aRBaIyE7Pta7RFoqWYL9S13cdRIW1otizVndg+99VyPBHy/1RfjsgajCQZoyJB9I891Uo2Ys9CNYUPcNQ+e5WbJH/wnUQlZtPxV1EKgH/BbyZa3MHINlzOxno6Mt7qABYgv0qneA6iIoIxbAFfgtw2HEW9QtfW+5jgUH89jy18saYDADPtZ7QHkLWrFkDO7EHwoq4TqMiRmPst0BtvYcMr4u7iNwDHDbGrPby9b1FJF1E0jMzM72NoQromWeegSuwLS2l/KUEdjGPTYB+nEOCLy33ZkB7EdkNzARaisjbwCERqQDguc7zi5oxZpIxJsEYk1C2bFkfYqj8WrlyJfPmzbN97Ve4TqMiTjPsyKvFroMo8KG4G2OeNMZUMsZUBToDnxtjugCpQKLnaYnAXJ9TKr8YNmwYZcqUgUauk6iIVAL7u7UJe3KccioQ49xHAq1FZCfQ2nNfObZs2TI+/fRTBg0aZEc2KBUITbHHcrT17pxfirsxZpEx5h7P7SPGmFbGmHjP9VF/vIfyzdNPP025cuXo27ev6ygqkhXHHqzfgo57d0zPUI0CixcvJi0tjcGDB1OiRAnXcVSka4w9prPIcY4op8U9whljGDx4MNdddx1JSUmu46hoUAzbPbMd+M5xliimxT3CzZ07l6+++orhw4dTrFgx13FUtGiMPcC6EDtJnQo6Le4R7OzZszz55JPUrFmTv/3tb67jqGhSBLgd23Lf4TZKtNL5ACNYcnIy27Zt44MPPqBQIf1RqyBrgF3CMQ2IR5uSQaa7O0KdPn2aYcOG0bhxYzp21Ol9lAOxQCvsaYwbHGeJQtqci1Djxo1j//79vPPOO4h4v6alUj6pBVwHfA7URmchDSJtuUeggwcP8vzzz9O+fXuaN2/uOo6KZgLcBZwAljnOEmW0uEegoUOHkpWVxZgxY1xHUQqqAjdh11o94TZKNNHiHmHWrl3L1KlTGTBgAPHx8a7jKGW1xk4MnuY6SPTQ4h5BjDEMHDiQMmXKMHToUNdxlPpVaezY9/XAfsdZooQW9wgyZ84cFi9ezIgRI7jmmmtcx1Hqt/6EPbHpU/TEpiDQ4h4hTp06xcCBA6lTpw49e/Z0HUepC12BHRq5F9uCVwGlQyEjxLPPPst3333HF198oScsqdB1C7AGmA/UwM5DowJCW+4RYMuWLYwZM4bu3btz2223uY6j1MXFAP8FnEEPrgaYFvcwZ4yhb9++lCxZkhdffNF1HKUurwJ2Dd909OBqAGlxD3MpKSksWrSIkSNHomvRqrBxB3Al8BF2iKTyOy3uYSwzM5OBAwfSqFEjPYiqwssVQBvsak1fOc4SobS4h7EBAwZw/PhxpkyZQkyM/ihVmLkZqI6dd+aI4ywRSCtCmEpNTWXmzJkMHTqU2rVru46jVMEJcA929sh5aPeMn2lxD0PHjh0jKSmJOnXqMHjwYNdxlPLeVdiJxXZjh0gqv9HiHob+93//l0OHDjF16lSKFCniOo5SvmkAVAPmw549e1yniRh6tkuYmTNnDtOmTYPboOHHDeFj14mU8pEA7YHxkJiYSFpaGrGxsa5ThT1tuYeRjIwMevXqZccJ3+46jVJ+VApoC4sXL+all15ynSYiaHEPE8YYHnroIU6fPg3/jX7nUpHnFrj33nsZOnQoa9eudZ0m7GlxDxOvv/46n376qV2AQ89VUpFIYOLEicTFxfHggw9y5swZ14nCmtfFXUSuF5H/iMhWEdksIo96tpcWkQUistNzXcp/caPTmjVrePzxx2nXrh0PP/yw6zhKBUyZMmVITk5m69at9O/f33WcsOZLy/0s8Lgx5ibsNPx9RaQWMBhIM8bEY6cG0rF6Pjh27Bh/+ctfKFeuHMnJybrYtYp4rVu3ZsiQIUyZMoXk5GTXccKW18XdGJNhjFnjuX0S2ApUBDoA538iyUBHX0NGK2MM3bt3Z+/evcyaNYu4uDjXkZQKiuHDh3PHHXfw8MMPs3HjRtdxwpIY4/uSKCJSFViCPaH4O2PMNbke+8EYc0HXjIj0BnoDVK5c+Y/RNr5VhuejBb4MO+91G6BJgAMpFQLMsF/r0cGDB6lfvz5XXXUV6enplCxZ0mGy0CQiq40xCXk95vMBVRG5EngfeMwYk++1zY0xk4wxCcaYBJ3NMA9fAwuwq8Y3dpxFKQeuvfZaZsyYwddff82DDz5ITk6O60hhxafiLiKFsYU9xRjzgWfzIRGp4Hm8AnDYt4hR6HvgPeyomI7YkzyUikK33347Y8eOZd68ebroewH5MlpGgCnAVmPMy7keSgUSPbcTgbnex4tCZ4AZ2J/MX4GibuMo5Vrfvn3p3bs3I0eO5O2333YdJ2z40nJvBnQFWorIOs+lHTASaC0iO4HWnvsqP3KwLfYfgP/BnrWnVJQTEV599VVatGhBz549WbZsmetIYcHr8xyNMUu5eIdBK2//3ah1Dvsd5xvsPBtVnaZRKqQUKVKE2bNn06RJE/785z+zdOlSbrrpJtexQpqeoRoqFgIbsMuPNXCcRakQFBcXx2effUbhwoVp06YN+/btcx0ppOkMJaFgmefSEGjuOItSDuVriHBHYBpc3+B66A4Ut5tzD6NU2nJ3bxV2LHstoC06Mkapy6mAHWxwFHgbOwhBXUCLu0vp2PnYq2NnetSfhlL5Uw076OAgtsD/5DZOKNJy4spq4CMgHvtLqh1kShVMDeB+IAN4C44fP+44UGjRkuKDfPUP5mU58Bm2sN+P/hSU8lYNbONoFrRq1YpPPvmEcuXKuU4VErTlHkwGOyrmM2wfuxZ2pXxXE+gMW7Zs4bbbbmP37t2uE4UELe7BkoM9d3cp8EfgL2hhV8pfqkNaWhrff/89TZs2ZcOGDa4TOafFPRhOYw/6rMUOdbwH3fNK+VmTJk344osviImJoVmzZqSmprqO5JSWmEA7DEwGvsOOz22JDndUKkBq167NihUrqFmzJh07duSFF17AH9OahyMt7oG0GXgT+Bn4G3CL0zRKRYWKFSuyZMkSOnfuzFNPPUXnzp05cSLfs5FHDC3ugZCNHeZ4ftre3sD1ThMpFVWKFStGSkoKL774Iu+//z4NGjRg9erVrmMFlRZ3fzuEba2nA02xp0df7TSRUlFJRBg0aBCLFy8mKyuLJk2a8PLLL0fNoh9a3P0lB7vQ4CTgJPAAcBc6IkYpx5o1a8a6deto27Ytjz/+OC1atGDHjh2uYwWcFnd/yMC21j/Hjrnti51SQCkVEsqUKcOHH35IcnIymzdvpl69eowePZrs7GzX0QJGi0uXXr8AAAhpSURBVLsvTgPzgInACeA+z6WEy1BKqbyICN26dWPLli20adOGQYMGUbduXRYsWOA6WkBocffCTz/9xNixY+FVYA12Aet+QG23uZRSl1ehQgU+/PBD5s2bR3Z2NnfddRedOnVi69atrqP5lYTCGNCEhASTnp7uOsZlZWdnk5yczPDhw+1CAdWAu4HyrpMppbySDSyHK1deyenTp+natSvDhg2jWrVqrpPli4isNsYk5PWYttzz4dSpU7z66qvEx8fTq1cvKlasyMKFC+3y31rYlQpfhYHm8M033zBw4EBmzpxJ9erV6datGxs3bnSdzida3C9hz549DBkyhCpVqjBgwAAqVapEamoqy5cvp1UrXSZWqUhRtmxZxowZw65du+jXrx8ffPABdevWpV27dnz88cdhOXwy6rtlLpi29yywCzvf+k7sVAHVgWZA5SCHU0oFxe+X6Dt69ChvvPEGr732GocOHaJKlSr06tWLrl27Urly6BSCS3XLaHEfLnaM+nfAJuyUAT8BV2IXqm4AXOMkmlIqSC62/mp2djZz585lwoQJpKWlAdC8eXMefPBBOnbs6HzueC3uefjhhx/4/PPP+cuIv8AO7DqMhbHj1OsCNwCxQY2klHIkP4trf/vtt7zzzju8/fbbbNu2DRGhadOmdOjQgbZt21K7dm1EgjsroBZ34NixYyxfvpylS5eSlpbGqlWrOHfuHFyB7XapCdwIFA1oDKVUCMpPcf/lucawYcMG5s6dy4cffsjatWsBO8TyzjvvpEWLFjRr1owaNWoEvNhHXXE/efIkmzdvZvXq1axevZpVq1axefNmjDHExsZy6623ctddd9G6dWtu++w2baErpbx3HKbWmcr8+fNJS0sjMzMTgLi4OBo2bEhCQgIJCQnUq1ePypUr+7XgR2xxP3nyJOnp6bR8uSUcAb4HMoFjuZ5UHLgOOytjZaAiUMT3zEopdQGDrUXfeS4HsDXpfJktApQD4oAy9rLj2R3Ex8d79XaXKu4Bm9ZKRO4GxmHbxW8aY0b6+z02b95My5Yt7Z1CQGmgEvYgaDmgAnAVujiGUio4BFu447B1COx6DgexC/ecv+wC1tmHh5ghzJo1y+9RAlLcRSQWeB1oDewDVolIqjFmiz/fp06dOixcuJA7P74TSqKj9pVSoacIttfg9yMos4AjMLTv0IC8baDK4a3ALmPMN8aYn4GZQAd/v0mJEiXsyURXo4VdKRVeigLXQd26dQPyzweqW6YisDfX/X1Ao9xPEJHe2DWKAH4Uke0BynIxcdhe+lAWDhkhPHJqRv8Jh5xhk1Ge8anfuMrFHghUcc8r7W+O3BpjJmGXtnBCRNIvdiAiVIRDRgiPnJrRf8Ihp2YMXGfGPn67amgl7HFjpZRSQRCo4r4KiBeRaiJSBOgMpAbovZRSSv1OQLpljDFnRaQf8Bl2KORUY8zmQLyXD5x1CRVAOGSE8MipGf0nHHJGfcaQOIlJKaWUf+kAQqWUikBa3JVSKgJFXHEXkUdFZJOIbBaRx/J4/AkRWee5bBKRHBEp7XnsbhHZLiK7RGRwCOfcLSIbPY8FbDrNfGS8WkTmich6z3O653osKPvSx4xB2Y/5zFlKROaIyAYRWSkiN+d6LFT25aUyBmxfishUETksIptybSstIgtEZKfnutRFXpvnvsvv6x1nfEZE9ueqA+0KFMoYEzEX4GbskhvFsQeLFwLxl3j+n4HPPbdjga+xM7kXAdYDtUItp+f+biDO9b4EngJe9NwuCxz17Lug7EtfMgZrPxYg52hgmOd2TSAtmL+XvmQM9L4EmmNnatmUa9soYLDn9uDzP+Pfve6i+y4/rw+BjM8Af/c2U6S13G8CvjLGnDbGnAUWA50u8fy/AjM8t4MyZYIfcgZLfjIaoKSICHbtqqPYhQqDtS99yRhM+clZC0gDMMZsA6qKSHlCa19eLGNAGWOWYH9uuXUAkj23k4GOebz0UvsuP693ndEnkVbcNwHNRaSMiBQH2vHbk6l+4Xn8buB9z6a8pkyoGII5wRas+SKyWuw0Dq4yvoYtCgeAjcCjxphzBG9f+pIRgrMf85tzPfDfACJyK/a08kqE1r68WEYI3r48r7wxJgPAc53XeneX2nf5eb3rjAD9PN1gUwvadRSwKX9dMMZsFZEXgQXAj9hfxou10v4MfGmMOf/X9rJTJviLjzkBmhljDohIOWCBiGzztByCnbENduLSlth1rBaIyBcEaV/6ktEYc4Ig7McC5BwJjBORddg/Qms9zwmlfXmxjBCkfVlAQftM++BSGccDIzz3RwAvAQ/l9x+OtJY7xpgpxpgGxpjm2K9JOy/y1M78tqsjqFMm+JATY8wBz/VhYA72q52LjN2BD4y1C/gW2xcbtH3pQ8ag7cf85DTGnDDGdDfG3AJ0wx4f+JYQ2peXyBjUfelxSEQqAHiuD+fxnEvtu/y83mlGY8whY0yO55vmZAq4TyOuuHtaDohIZexXyAv6qkXkaqAFMDfX5qBOmeBtThEpISIlz98G7sJ+pXaR8Tuglec55YEawDcEcV96mzGY+zE/OUXkGs++AugJLPF8uwiZfXmxjMHelx6pQKLndiK//Syfd6l9l5/XO814/g+DRycKuk+9PRIbqhfgC2AL9mtlK8+2JCAp13P+BszM47XtgB3Yo9dDQjEn9qj6es9lcyBzXi4jdgHD+div6JuALsHel95mDOZ+zGfOJtiW8jbgA6BUCO7LPDMGel9i/8hkANnYlm4P7CJ1aZ48aUDpXD/vTy637y72+hDL+Jbn93YDtuBXKEgmnX5AKaUiUMR1yyillNLirpRSEUmLu1JKRSAt7kopFYG0uCulVATS4q6UUhFIi7tSSkWg/w9niHKk5IVBzAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "def Gauss_dist(t,lam,sigma):\n", " return np.exp(-0.5*((t-lam)/sigma)**2)/np.sqrt(2*np.pi)/sigma\n", "\n", "fig1, ax1 = plt.subplots()\n", "\n", "nbins = 20\n", "range1 = (L_mean - 3*L_err,L_mean + 3*L_err)\n", "ax1.hist(L_test,nbins,range1,color='b',label='L');\n", "\n", "fig2, ax2 = plt.subplots()\n", "\n", "nbins = 20\n", "range2 = (T_mean - 3*T_err,T_mean + 3*T_err)\n", "ax2.hist(T_test,nbins,range2,color='r',label='T');\n", "\n", "\n", "fig3, ax3 = plt.subplots()\n", "\n", "nbins = 20\n", "range3 = (g_m1 - 3*g_e1,g_m1 + 3*g_e1)\n", "binwidth = 6*g_e1/nbins\n", "\n", "ax3.hist(g_test,nbins,range3,color='g',label='g');\n", "\n", "t = np.linspace(g_m1 - 3*g_e1,g_m1 + 3*g_e1,100)\n", "yval = Gauss_dist(t,g_m1,g_e1)*npoints*binwidth\n", "#yval = np.array([np.exp(-0.5*((t1-g_m1)/g_e1)**2)/np.sqrt(2*np.pi)/g_e1 for t1 in t])\n", "#print(yval)\n", "\n", "ax3.plot(t,yval,c='k')\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 11.4 Simulazione di eventi complessi con metodo MC" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.12" } }, "nbformat": 4, "nbformat_minor": 4 }