{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Optical Crosstalk in SPAD arrays\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Theory\n",
"\n",
"Optical crosstalk is a interaction between pixels on the same \n",
"array. When an avalanche is triggered on one pixel it emits \n",
"secondary photons that can trigger another pixel on the same chip. \n",
"This secondary photon emission is proportional to the avalanche current \n",
"whose duration is typically < 20ns in SPADs using AQCs. Hence, when a\n",
"crosstalk event occurs, the second pixels is triggered with a delay\n",
"< 20ns. We can estimate the crosstalk by \"coincidence counting\",\n",
"that is counting the number of photons $C$ arriving in two pixels \n",
"within a short time window $\\Delta t$ (e.g. 25-50ns). \n",
"The number of coincident events due to uncorrelated dark counts \n",
"can be computed from Poisson statistics. Then, the number\n",
"of crosstalk events is simply the total coincidence counts minus\n",
"the coincidence due to uncorrelated counts. The crosstalk probability\n",
"is the number of coincidences divided by the total counts.\n",
"\n",
"**Reference**:\n",
"- [Restelli JMO 2006](http://dx.doi.org/10.1080/09500340600790121).\n",
"\n",
"## Poisson statistics\n",
"\n",
"Given the random variable $X \\sim {\\rm Poiss}\\{\\Lambda\\}$, \n",
"where $\\Lambda = \\lambda\\Delta t$. We ask, what is the probability\n",
"of having at least one count in a time window $\\delta t$?\n",
"\n",
"$$P\\{X \\ge 1\\} = 1 - P\\{X = 0\\} = 1 - e^\\Lambda \n",
"\\xrightarrow[\\scriptstyle\\Lambda\\to0]{} \\Lambda$$\n",
"\n",
"Hence, for two random variables $X_A \\sim {\\rm Poiss}\\{\\Lambda_A\\}$\n",
"and $X_B \\sim {\\rm Poiss}\\{\\Lambda_B\\}$, the propbability of having\n",
"at least one count in each variable in a time $\\Delta t$ is:\n",
"\n",
"$$P\\{X_A \\ge 1\\}P\\{X_B \\ge 1\\} = (1 - e^{-\\Lambda_A})(1 - e^{-\\Lambda_B})\n",
"\\approx (\\lambda_A\\Delta t)(\\lambda_B\\Delta t) \n",
"\\quad {\\rm for}\\quad\\Delta t \\ll \\lambda^{-1}$$\n",
"\n",
"In a measurement of duration $T$, number of \"coincidences\" $C_u$ is the number\n",
"of times both variables have at least one count in a window $\\Delta t$\n",
"\n",
"$$C_u = P\\{X_A \\ge 1\\}P\\{X_B \\ge 1\\} \\frac{T}{\\Delta t}\n",
"\\approx \\lambda_A\\,\\lambda_B\\, T\\,\\,\\Delta t \n",
"\\quad {\\rm for}\\quad\\Delta t \\ll \\lambda^{-1}$$\n",
"\n",
"## Crosstalk probability\n",
"\n",
"Given a measurement of duration $T$, with total counts $N_A$ and $N_B$ in\n",
"pixels $A$ and $B$. If $C$ are the total number of coincidence in windows of \n",
"duration $\\Delta t$, then the coincidence $C_c$ due to crosstalk are:\n",
"\n",
"\\begin{equation}\n",
"C_c = C - C_u \\qquad\\qquad (1)\n",
"\\end{equation}\n",
"\n",
"where $C_u$ are the number of coincidences due to Poisson statistics, \n",
"i.e. the coincidences we would have if the two pixels were uncorrelated.\n",
" \n",
"\\begin{equation}\n",
"C_u = \\left[1 - \\exp\\left(-\\frac{N_A - C_c}{T}\\right)\\right]\n",
" \\left[1 - \\exp\\left(-\\frac{N_B - C_c}{T}\\right)\\right]\n",
" \\frac{T}{\\Delta t}\n",
"\\end{equation}\n",
"\n",
"The expression of $C_u$ can be replaced eq. (1) and then solved for $C_c$ \n",
"(a numerical iterative solution is straightforward).\n",
"For simplicity, we could assume $C_c \\ll N_A,N_B$ obtaining:\n",
"\n",
"\\begin{equation}\n",
"C_u = \\left[1 - \\exp\\left(-\\frac{N_A}{T}\\right)\\right]\n",
" \\left[1 - \\exp\\left(-\\frac{N_B}{T}\\right)\\right]\n",
" \\frac{T}{\\Delta t}\n",
"\\end{equation}\n",
"\n",
"In either cases, the probability of crosstalk is:\n",
"\n",
"\\begin{equation}\n",
"P_c = \\frac{C_c}{N_A + N_B - C_c}\n",
"\\end{equation}\n",
"\n",
"The counts $C_c$ are independetent events and thus are Poisson distributed.\n",
"The standard deviation of $C_c$ is then\n",
"${\\rm Std}\\{C_c\\} = \\sqrt{C_c}$, and the standard deviation of $P_c$ is:\n",
"\n",
"\\begin{equation}\n",
"{\\rm Std}\\{P_c\\} = \\frac{\\sqrt{C_c}}{N_A + N_B - C_c}\n",
"\\end{equation}\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Define functions"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"def coincidence_py(timestamps1, timestamps2):\n",
" \"\"\"Pure python coincidence counting.\"\"\"\n",
" coinc = 0\n",
" i1, i2 = 0, 0\n",
" while (i1 < timestamps1.size) and (i2 < timestamps2.size):\n",
" if timestamps1[i1] == timestamps2[i2]:\n",
" coinc += 1\n",
" i1 += 1\n",
" i2 += 1\n",
" elif timestamps1[i1] > timestamps2[i2]:\n",
" i2 += 1\n",
" elif timestamps1[i1] < timestamps2[i2]:\n",
" i1 += 1\n",
" return coinc"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"%load_ext Cython"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"%%cython\n",
"cimport numpy as np\n",
"\n",
"def coincidence(np.int64_t[:] timestamps1, np.int64_t[:] timestamps2):\n",
" \"\"\"Cython coincidence counting.\"\"\"\n",
" cdef np.int64_t coinc, i1, i2, size1, size2\n",
" size1 = timestamps1.size\n",
" size2 = timestamps2.size\n",
" coinc = 0\n",
" i1, i2 = 0, 0\n",
" while i1 < size1 and i2 < size2:\n",
" if timestamps1[i1] == timestamps2[i2]:\n",
" coinc += 1\n",
" i1 += 1\n",
" i2 += 1\n",
" elif timestamps1[i1] > timestamps2[i2]:\n",
" i2 += 1\n",
" elif timestamps1[i1] < timestamps2[i2]:\n",
" i1 += 1\n",
" return coinc"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"def crosstalk_probability(t1, t2, tol=1e-6, max_iter=100):\n",
" \"\"\"Estimate crosstalk probability between two pixels in a SPAD array.\n",
" \n",
" Given two input arrays of timestamps `t1` and `t2`, estimate \n",
" the crosstalk probability using Poisson statistics without \n",
" approximations.\n",
" \n",
" Arguments:\n",
" t1, t2 (array of ints): arrays of timestamps from DCR measurements\n",
" for the two pixels to be measured. Timestamps need to be \n",
" integers and coincidences are detected when values in the two \n",
" arrays are equal. These arrays need to be rescaled, if \n",
" coincidence need to be computed on a delta t larger than \\\n",
" a single timestamp unit.\n",
" tol (float): tollerance for iterative equasion solution\n",
" max_iter (int): max iterations used to solve the equation\n",
" \n",
" Returns:\n",
" A 3-element tuple:\n",
" - crosstalk probability\n",
" - crosstalk probability standard deviation\n",
" - number of iterations used for the estimation.\n",
" \"\"\"\n",
" T = (max((t1.max(), t2.max())) - min((t1.min(), t2.min())))\n",
" C = coincidence(t1, t2)\n",
" \n",
" # Find C_c by solving eq. (1) iteratively\n",
" C_c, C_u_prev = 0, 0\n",
" for i in range(max_iter):\n",
" C_u = ((1 - np.exp(-(t1.size - C_c)/T)) * \n",
" (1 - np.exp(-(t2.size - C_c)/T)) * T)\n",
" C_c = C - C_u\n",
" if np.abs(C_u - C_u_prev) < tol:\n",
" break\n",
" C_u_prev = C_u\n",
" \n",
" P_c = C_c / (t1.size + t2.size - C_c)\n",
" sigma = np.sqrt(C_c) / (t1.size + t2.size - C_c)\n",
" return P_c, sigma, i"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Simulation\n",
"\n",
"## Poisson processes\n",
"\n",
"Here we simulate two poisson processes to check the the coincidences \n",
"are, as predicted, equal to $C_u$."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"import numpy as np\n",
"%matplotlib inline\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Simulation parameters:"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"λ_A = 1000 # cps\n",
"λ_B = 2000 # cps\n",
"T = 600 # s\n",
"delta_t = 50e-9 # s\n",
"P_c = 0.05"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"From theory, the number of coincidences are:"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"60.0"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"C_u = λ_A * λ_B * T * delta_t\n",
"C_u"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"7.745966692414834"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"np.sqrt(C_u)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's simulate timestamps for two uncorrelated Poisson processes:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(0.0010001265879312606, 1000.1265879312606)"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"dt = np.random.exponential(scale=1/λ_A, size=10**6)\n",
"dt.mean(), dt.sum()"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"np.random.seed(1)\n",
"\n",
"t_A = np.cumsum(np.random.exponential(scale=1/λ_A, size=λ_A * T * 2))\n",
"assert t_A.max() > T\n",
"t_A = t_A[t_A < 600]\n",
"t_A = (t_A / delta_t).astype('int64')\n",
"\n",
"t_B = np.cumsum(np.random.exponential(scale=1/λ_B, size=λ_B * T * 2))\n",
"assert t_B.max() > T\n",
"t_B = t_B[t_B < 600]\n",
"t_B = (t_B / delta_t).astype('int64')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The empirical coincidences are:"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"62"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"C_u_sim = coincidence(t_A, t_B)\n",
"C_u_sim"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"assert C_u - 3*np.sqrt(C_u) < C_u_sim < C_u + 3*np.sqrt(C_u)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's repeat the simulation $N$ times to be sure that the \n",
"estimated coincidence falls within the espected error range\n",
"all the times:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
">>> In 500 simulations, the number of coincincence was outside the error range 0 times\n"
]
}
],
"source": [
"N = 500\n",
"fail = 0\n",
"for _ in range(N):\n",
" t_A = np.cumsum(np.random.exponential(scale=1/λ_A, size=λ_A * T * 2))\n",
" assert t_A.max() > T\n",
" t_A = t_A[t_A < T]\n",
" t_A = (t_A / delta_t).astype('int64')\n",
"\n",
" t_B = np.cumsum(np.random.exponential(scale=1/λ_B, size=λ_B * T * 2))\n",
" assert t_B.max() > T\n",
" t_B = t_B[t_B < T]\n",
" t_B = (t_B / delta_t).astype('int64')\n",
" \n",
" C_u_sim = coincidence(t_A, t_B)\n",
" if C_u < C_u_sim - 3*np.sqrt(C_u_sim) or C_u > C_u_sim + 3*np.sqrt(C_u_sim):\n",
" fail += 1\n",
"\n",
"print('>>> In %d simulations, the number of coincincence was outside the error range %d times' \n",
" % (N, fail/N))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For further information of confidence intervals for Poisson distribution\n",
"estimators see:\n",
"\n",
"- Patil & Kulkarni. *Statistical Journal* 10(2) pp211–227(2012). \n",
" COMPARISON OF CONFIDENCE INTERVALS FOR THE POISSON MEAN: SOME NEW ASPECTS. \n",
" [[PDF]](https://www.ine.pt/revstat/pdf/rs120203.pdf) "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Optical crosstalk\n",
"\n",
"Here we simulate two Poisson processes plus crosstalk, to check \n",
"that the estimated crosstalk is consistent with the ground truth."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"λ_A = 1000 # cps\n",
"λ_B = 1550 # cps\n",
"T = 1200 # s\n",
"delta_t = 50e-9 # s\n",
"P_c = 0.1"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"np.random.seed(1)"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"t_A = np.cumsum(np.random.exponential(scale=1/λ_A, size=λ_A * T * 2))\n",
"assert t_A.max() > T\n",
"t_A = t_A[t_A < 600]\n",
"t_A = (t_A / delta_t).astype('int64')\n",
"t_ct_AB = t_A[np.random.rand(t_A.size) <= P_c]\n",
"\n",
"t_B = np.cumsum(np.random.exponential(scale=1/λ_B, size=λ_B * T * 2))\n",
"assert t_B.max() > T\n",
"t_B = t_B[t_B < 600]\n",
"t_B = (t_B / delta_t).astype('int64')\n",
"t_ct_BA = t_B[np.random.rand(t_B.size) <= P_c]\n",
"\n",
"t_B = np.hstack([t_B, t_ct_AB])\n",
"t_B.sort()\n",
"t_A = np.hstack([t_A, t_ct_BA])\n",
"t_A.sort()"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(9.9967084810755544, 0.076659160092492767, 4)"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P_c_est, P_c_err, i = crosstalk_probability(t_A, t_B, delta_t)\n",
"P_c_est*100, P_c_err*300, i"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's repeat the simulation $N$ times to obtain a distribution\n",
"of estimated crosstalk values:"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"λ_A = 1000 # cps\n",
"λ_B = 300 # cps\n",
"T = 1200 # s\n",
"delta_t = 50e-9 # s\n",
"P_c = 0.05"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"3.0"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"N = 100\n",
"P = []\n",
"I = 0\n",
"for _ in range(N):\n",
" t_A = np.cumsum(np.random.exponential(scale=1/λ_A, size=λ_A * T * 2))\n",
" assert t_A.max() > T\n",
" t_A = t_A[t_A < T]\n",
" t_A = (t_A / delta_t).astype('int64')\n",
" t_ct_AB = t_A[np.random.rand(t_A.size) <= P_c]\n",
"\n",
" t_B = np.cumsum(np.random.exponential(scale=1/λ_B, size=λ_B * T * 2))\n",
" assert t_B.max() > T\n",
" t_B = t_B[t_B < T]\n",
" t_B = (t_B / delta_t).astype('int64')\n",
" t_ct_BA = t_B[np.random.rand(t_B.size) <= P_c]\n",
"\n",
" t_B = np.hstack([t_B, t_ct_AB])\n",
" t_B.sort()\n",
" t_A = np.hstack([t_A, t_ct_BA])\n",
" t_A.sort()\n",
" \n",
" P_c_est, P_c_err, i = crosstalk_probability(t_A, t_B)\n",
" I += i\n",
" P.append(P_c_est*100)\n",
"I/N"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEWCAYAAAB2X2wCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xm8XHV9//HXmwQIiQkQclkkhCug\nKCCKXMFCEUSsCMjihggKak21lYrF0ri04I7VWqz9tTUiDQpCBQGV4MJiDLIkJhDWgGxBwpYbQshC\nAgl8fn98v1eGYeZu59x7bmbez8fjPu6Zs8x5z5kzn/nO95w5o4jAzMw2fBtVHcDMzMrhgm5m1iJc\n0M3MWoQLuplZi3BBNzNrES7oZmYtwgW9nyQdL+k3Q3TfMyR9pcDyqyTtVGamfqxzM0m/kPSUpIuG\nc91DQdIsSX+dh0+S9PsKMgz782itxQW9hqS/lHR9LlLLJF0n6Y0AEXF+RPzVCMj458LTIyJeFhH3\nD3OU9wDbAFtFxHuHed19khSSdqk6RzPD+TxW9QZlw2901QFGCkkTgMuBTwA/ATYBDgCeqTLXCLYj\n8MeIWF/kTiSNLnofNrRG4nM0EjONCBHhv/Rt2S5geS/TTwJ+X3M7gL8F7gFWAl8GdgZuAFaQ3xQa\nLVuz/C55eAbwlTy8JemNpRt4Mg9PztO+CjwHrAVWAf/Z4L42B36Yl38Q+AKwUW0O4Fv5vh8A3tHL\nY34NMAtYDtwBHJnHfxF4FliXc3y0wbKjgM8B9+XtMx/YoSbv3+Vt90Aetx/wB+Cp/H+/um1/f76f\nB4Dj8/hdgN/lZZYC/5fHz87rWJ3zHdvbds3LzAL+uslz/c283TZv8Dg3Aqblx/lEft4n5mljgPPy\n+OX5cW3Tz+dxBvBfwC/zPNcB2wJn5fx3AXvV5OjJsBK4Ezim5jlcm9e3iryP0/d+ch3w78Ay8r5Z\n97j3Ie3ry4FHgf8k7+95+u7AlXn5x4HP9bZfAJ358Y/u5Tl5USbS6+2avH2XAucDW9QsvwNwSX6M\nT+SMm+blX1sz39bAGqCj6jpUuI5VHWCk/AET8pN+LvAOYMu66fUv8gB+npfbndSSvxrYKb9Y7gRO\nbLRszfKNCvpWwLuBscB44CLgsprl/ryTN7mvHwI/y8t2An8kF9ycYx3wsfzC+gTwCKAG22Nj4N78\n4tsEODi/AHfN088Azutle/4jcBuwKyDgdaTumZ68VwITgc3y/yeBD5I+NR6Xb28FjCO9Qfasdztg\n9zx8AfB5UlEdA/xlo20y0O3KC298GwHfB34NjG3yOE8BbgQmk4rF94AL8rS/AX6R1zkK2BuY0M/n\ncQapSO2dH9s1pDezD+X7+grw25pl3wu8PGc+lvRmtl0v+19f+8l64OT8fGzW4HHvDbwpT+8EFgKn\n5GnjSUX+1Jx9PLBvb/sF/SvoL8pEekN/W97uHaQ38rPy/KOAW0hvAOOo2T9Ib5TfqFnPp4BfVF2D\nSqljVQcYSX+k1swMYHHeeX4ObFOzQ9UX9P1rbs8H/qnm9r/V7FyNXlANC3qDTK8Hnqy5/eedvP6+\n8k78DLBbzbS/AWbV5Li3ZtrYvOy2DdZ7APAYudWWx10AnJGHz6D3gn43cFSTaQEcXHP7g8Dcunlu\nyHnHkVqB76ausJCK0nRqWtqNtu9At2te7xzg/4CfUtPybHA/C4G31tzejvSmORr4CHA9sGeD5Zo+\njzX7xPdrpp0MLKy5/Vp6/0S5oGf7N9h3+7Of/GmAr51TgEvz8HHAzQPZL+hfQe81E3B0z3qBvyC1\nzEc3mG9f4CFe+EQyD3jfQB7vSP3zQdEaEbEwIk6KiMnAHqQWz1m9LPJ4zfCaBrdfNtAMksZK+p6k\nByWtILU6tpA0qh+LTyK1ph+sGfcgsH3N7cd6BiLi6TzYKOfLgYci4vle7qs3O5A+VjfzUN26Hqyb\n/iCwfUSsJrU4Pw48KmmmpFfneU4jtfLmSrpD0kearWwQ23UX4CjgixHxbC+PY0fgUknLJS0nFfjn\nSF0rPyK17i+U9Iikf5W0cS/3Va/f+5ekD0laUJNjD9L+0Eh/9pOH6IWkV0m6XNJjeXt+rWZ9vT33\nfe0XvXlRJklbS7pQ0sM5w3l1GR6MBv3sETGH9AnmwLwv7UJqvG3wXNCbiIi7SK2kPUq4u9Wk1jAA\nkrbtZd5TSR9H942ICcCbexbridbLsktJrcMda8ZNAR4eaGBSV8wOkmr3kYHc10OkPs5mah/HI7w4\n84vWFRG/joi3kVq/d5G6QYiIxyLiYxHxclIL8796ObOlr+1abyHwYeCXknbt5XE8RDoOsUXN35iI\neDgi1kXEFyNiN9IxgiNIXSb1j78QSTuStsknSd1aWwC303yf6c9+0le+/yY9F6/M2/NzNevr7blv\nNm11/j+2Zlz966Q+09fzuD1zhhPqMkyR1OzEj3Pz/B8ELo6ItU3m26C4oGeSXi3pVEmT8+0dSB8d\nbyzh7m8Bdpf0ekljSN0VzYwntb6WS5oInF43/XFSP/1LRMRzpINyX5U0Pr/Q/4HUchmonlbMaZI2\nlnQQ8E7gwn4ufzbwZUmvVLKnpK2azHsF8CpJH5A0WtKxwG7A5ZK2kXSkpHGkboJVpBYwkt7b83yR\n+tyjZxov3U59bdeXiIgLSIXqKknNCtT/kLb3jjlTh6Sj8vBbJL02fwpYQSqizfIVMY702Lvzej/M\nixsijwOTJW2SH1cZ+8l40mNalVu5n6iZdjmwraRTJG2a17FvntZwv4iIbtIbygmSRuVPW701CHoy\nrCI9p9uT+ud7zCX1458paZykMZL2r5n+I+AYUlH/4QAe94jmgv6ClaS+tTmSVpMK+e2kll0hEfFH\n4EvAVaQzO3o7J/gs0gGfpTnDr+qmfwd4j6QnJf1Hg+VPJhXi+/N6fgycM4jMzwJHkg4QLyUdSPpQ\n/uTSH98mFY3fkF74PyA9rkbreoLUej2VdGD6NOCIiFhK2kdPJbXilwEHks4uAngj6flaRfrI/KmI\neCBPOwM4N3dBvI++t2tDEXEu6bm7RlJng1m+k9f9G0kr8333FK9tgYvz419IOiPnvJrlense+y0i\n7iQds7mBVLxfSzojpMc1pLOUHpO0NI8rup98BvgA6XXzfdLxhp48K0kHK99J6uK7B3hLntzbfvEx\nUlF+gnSiwfV9ZPgi8AbSWU4zSWe09GR4Lq9/F+BPpONix9ZMXwzcRHojvHYAj3tEUz4oYGbWViSd\nAzwSEV+oOktZ/MUiM2s7+dPWu4C9qk1SLne5mFlbkfRlUnfqN2u66FqCu1zMzFqEW+hmZi1iWPvQ\nJ02aFJ2dncO5SrM+3X333QDsumtvp5ubVWf+/PlLI6Kjr/mGtaB3dnYyb9684VylWZ8OOuggAGbN\nmlVpDrNmJNV/k7ohd7mYmbUIn7Zobe8LX2iZ05CtzbmgW9s75JBDqo5gVgp3uVjbW7BgAQsWLKg6\nhllhbqFb2zvllFMAHxS1DZ9b6GZmLcIF3cysRbigm5m1CBd0M7MW4YOi1nY6p8180e21k49oOH7R\nmYcPWyazMrigW9sbM/k1VUcwK4W7XKztrV28kLWLF1Ydw6wwF3Rre8tnn8vy2edWHcOsMBd0M7MW\n4YJuZtYiXNDNzFqEC7qZWYvwaYvW9ia+dWrVEcxK0WcLXdI5kpZIur1u/MmS7pZ0h6R/HbqIZkNr\nk212YpNtdqo6hllh/elymQEcWjtC0luAo4A9I2J34FvlRzMbHmsWLWDNIl8P3TZ8fXa5RMRsSZ11\noz8BnBkRz+R5lpQfzWx4PHX9hQBs1vn6ipOYFTPYg6KvAg6QNEfS7yS9sdmMkqZKmidpXnd39yBX\nZ2ZmfRlsQR8NbAm8CfhH4CeS1GjGiJgeEV0R0dXR0THI1ZmZWV8GW9AXA5dEMhd4HphUXiwzMxuo\nwRb0y4CDASS9CtgEWFpWKDMzG7g+D4pKugA4CJgkaTFwOnAOcE4+lfFZ4MSIiKEMajZUtnr7J6uO\nYFaK/pzlclyTSSeUnMWsEhtvNbnqCGal8Ff/re09fe8cnr53TtUxzArzV/+t7a2YeykAY3fZt+Ik\nZsW4hW5m1iJc0M3MWoQLuplZi3BBNzNrET4oam1v0hGnVh3BrBQu6Nb2Rk/wNYasNbjLxdre6oWz\nWb1wdtUxzApzC93a3sqbrwBg3GveXHESs2LcQjczaxEu6GZmLcIF3cysRbigm5m1CB8UtbbXcfRn\nq45gVoo+W+iSzpG0JP+YRf20z0gKSf75OdtgjRq7OaPGbl51DLPC+tPlMgM4tH6kpB2AtwF/KjmT\n2bBaddtVrLrtqqpjmBXWZ0GPiNnAsgaT/h04DfBPz9kGzQXdWsWgDopKOhJ4OCJu6ce8UyXNkzSv\nu7t7MKszM7N+GHBBlzQW+DzwL/2ZPyKmR0RXRHR1dPiaGWZmQ2UwLfSdgVcAt0haBEwGbpK0bZnB\nzMxsYAZ82mJE3AZs3XM7F/WuiFhaYi4zMxugPgu6pAuAg4BJkhYDp0fED4Y6mNlw2fq9Z1QdwawU\nfRb0iDiuj+mdpaUxq8BGG4+pOoJZKfzVf2t7K2+aycqbZlYdw6wwF3Rre6vvupbVd11bdQyzwlzQ\nzcxahAu6mVmLcEE3M2sRLuhmZi3C10O3ltY5re+zV7b9wJnDkMRs6LmFbmbWIlzQre09NecSnppz\nSdUxzApzQbe2t+a+uay5b27VMcwKc0E3M2sRLuhmZi3CBd3MrEX4tEVrexq9adURzErhgm5tb5v3\nfbHqCGal6LPLRdI5kpZIur1m3Dcl3SXpVkmXStpiaGOamVlf+tOHPgM4tG7clcAeEbEn8EfgsyXn\nMhs2y6+7gOXXXVB1DLPC+izoETEbWFY37jcRsT7fvJH0Q9FmG6S1D97C2gdvqTqGWWFlnOXyEeCX\nJdyPmZkVUKigS/o8sB44v5d5pkqaJ2led3d3kdWZmVkvBl3QJZ0IHAEcHxHRbL6ImB4RXRHR1dHR\nMdjVmZlZHwZ12qKkQ4F/Ag6MiKfLjWQ2vEZtNqHqCGal6LOgS7oAOAiYJGkxcDrprJZNgSslAdwY\nER8fwpxmQ6bjmM9VHcGsFH0W9Ig4rsHoHwxBFjMzK8DXcrG29+TvZvDk72ZUHcOsMH/139reMw/f\nVXUEs1K4hW5m1iJc0M3MWoQLuplZi3AfurW90eMnVR3BrBQu6Nb2Jr3zM1VHMCuFC7pZE53TZr5k\n3KIzD68giVn/uA/d2t6yq6az7KrpVccwK8wtdGt7zy65v+oIZqVwC93MrEW4oJuZtQgXdDOzFuE+\ndGt7G0/cvuoIZqVwQbe2t9WhJ1cdwawU7nIxM2sRfRZ0SedIWiLp9ppxEyVdKeme/H/LoY1pNnSe\n+NV3eeJX3606hllh/WmhzwAOrRs3Dbg6Il4JXJ1vm22Q1i17mHXLHq46hllhfRb0iJgNLKsbfRRw\nbh4+Fzi65FxmZjZAg+1D3yYiHgXI/7duNqOkqZLmSZrX3d09yNWZmVlfhvygaERMj4iuiOjq6OgY\n6tWZmbWtwZ62+Lik7SLiUUnbAUvKDGU2nDbZeqeqI5iVYrAF/efAicCZ+f/PSktkNswmHjK16ghm\npejPaYsXADcAu0paLOmjpEL+Nkn3AG/Lt83MrEJ9ttAj4rgmk95achazSiz9xbcA/3KRbfj81X9r\ne+tXLq06glkp/NV/M7MW4YJuZtYiXNDNzFqE+9Ct7W26/aurjmBWChd0a3tbHnhS1RHMSuGCbi2j\nc9rMqiOYVcp96Nb2ui/9Gt2Xfq3qGGaFuYVube+5NSuqjmBWCrfQzcxahAu6mVmLcEE3M2sR7kO3\ntjdmx9dVHcGsFC7o1va22L/ZBUXNNizucjEzaxGFCrqkT0u6Q9Ltki6QNKasYGbD5fGfnM7jPzm9\n6hhmhQ26oEvaHvh7oCsi9gBGAe8vK5jZcIn1zxDrn6k6hllhRbtcRgObSRoNjAUeKR7JzMwGY9AF\nPSIeBr4F/Al4FHgqIn5TP5+kqZLmSZrX3d09+KRmZtarIl0uWwJHAa8AXg6Mk3RC/XwRMT0iuiKi\nq6OjY/BJzcysV0VOWzwEeCAiugEkXQLsB5xXRjCz4bLZzvtUHcGsFEUK+p+AN0kaC6wB3grMKyWV\n2TDafN93VR3BrBRF+tDnABcDNwG35fuaXlIuMzMboELfFI2I0wGfwGsbtMd+PA2AbT9wZsVJzIrx\nN0XNzFqEC7qZWYtwQTczaxEu6GZmLcKXz7W2N+7VB1QdwawULujW9sa/4fCqI5iVwl0u1vaeX7eW\n59etrTqGWWFuoVvbW3LRGUD/zkPvnDaz4fhFZ760ld9o3kbzmZXFLXQzsxbhgm5m1iJc0M3MWoQL\nuplZi/BBUWt7L3vtIVVHMCuFC7q1PRd0axXucrG299zTT/Hc009VHcOsMBd0a3vdl32d7su+XnUM\ns8IKFXRJW0i6WNJdkhZK+ouygpmZ2cAU7UP/DvCriHiPpE2AsSVkMjOzQRh0QZc0AXgzcBJARDwL\nPFtOLDMzG6giLfSdgG7gfyW9DpgPfCoiVtfOJGkqMBVgypQpBVZnNnI1u8aL2XAq0oc+GngD8N8R\nsRewGphWP1NETI+Irojo6ujoKLA6s6Exfq/DGL/XYVXHMCusSAt9MbA4Iubk2xfToKCbjXTjXvPm\nqiOYlWLQLfSIeAx4SNKuedRbgTtLSWU2jNav6Gb9iu6qY5gVVvQsl5OB8/MZLvcDHy4eyWx4Lb38\n34D+XQ/dbCQrVNAjYgHQVVIWMzMrwN8UNTNrES7oZmYtwgXdzKxF+PK51vYm7HNM1RHMSuGCbm1v\n7C77Vh3BrBTucrG2t+6Jxax7YnHVMcwKcwvdRrxG10lZdObhpd3/E7/+T8DnoduGzy10M7MW4YJu\nZtYiXNDNzFqEC7qZWYvwQVFre5vv9/6qI5iVwgXd2t5mna+vOoJZKdzlYm3v2cfv59nH7686hllh\nLujW9pZdPZ1lV0+vOoZZYYULuqRRkm6WdHkZgczMbHDKaKF/ClhYwv2YmVkBhQq6pMnA4cDZ5cQx\nM7PBKnqWy1nAacD4ZjNImgpMBZgyZUrB1Zklja7vYtbuBl3QJR0BLImI+ZIOajZfREwHpgN0dXXF\nYNdnNlS2ePOJVUcwK0WRFvr+wJGSDgPGABMknRcRJ5QTzWx4jJn8mqojmJVi0H3oEfHZiJgcEZ3A\n+4FrXMxtQ7R28ULWLvZxfdvw+Tx0a3vLZ5/L8tnnVh3DrLBSvvofEbOAWWXcl5mZDY5b6GZmLcIF\n3cysRbigm5m1CF8+19rexLdOrTqCWSlc0K3tbbLNTlVHMCuFu1ys7a1ZtIA1ixZUHcOsMLfQre09\ndf2FgH+5yDZ8bqGbmbUIF3Qzsxbhgm5m1iJc0M3MWoQPilrb2+rtn6w6glkpXNCt7W281eSqI5iV\nwl0u1vaevncOT987p+oYZoW5hW5tb8XcSwEYu8u+FScxK2bQLXRJO0j6raSFku6Q9Kkyg5mZ2cAU\naaGvB06NiJskjQfmS7oyIu4sKZuZmQ1Akd8UfTQibsrDK4GFwPZlBTMzs4EppQ9dUiewF/CSI0uS\npgJTAaZMmVLG6qyFdU6bWXUEsw1W4YIu6WXAT4FTImJF/fSImA5MB+jq6oqi6zMr26QjTq06glkp\nChV0SRuTivn5EXFJOZHMhtfoCR1VRzArRZGzXAT8AFgYEd8uL5LZ8Fq9cDarF86uOoZZYUW+WLQ/\n8EHgYEkL8t9hJeUyGzYrb76ClTdfUXUMs8IG3eUSEb8HVGIWMzMrwF/9NzNrES7oZmYtwgXdzKxF\n+OJc1vY6jv5s1RHMSuGCbm1v1NjNq45gVgp3uVjbW3XbVay67aqqY5gV5oJubc8F3VqFu1xsWPii\nW2ZDzy10M7MW4YJuZtYiXNDNzFqE+9Ct7W393jOqjmBWChd0a3sbbTym6ghmpXCXi7W9lTfNZOVN\nPgvHNnwu6Nb2Vt91LavvurbqGGaFuaCbmbWIQgVd0qGS7pZ0r6RpZYUyM7OBK/KboqOA/we8A9gN\nOE7SbmUFMzOzgSnSQt8HuDci7o+IZ4ELgaPKiWVmZgNV5LTF7YGHam4vBvatn0nSVGBqvrlK0t2D\nXN8kYOkglx1KzjUwIzbXg984Yshz6RsDXmTEbi+ca6CKZNuxPzMVKeiNfiA6XjIiYjowvcB60sqk\neRHRVfR+yuZcA+NcA+NcAzNSc8HwZCvS5bIY2KHm9mTgkWJxzMxssIoU9D8Ar5T0CkmbAO8Hfl5O\nLDMzG6hBd7lExHpJnwR+DYwCzomIO0pL9lKFu22GiHMNjHMNjHMNzEjNBcOQTREv6fY2M7MNkL8p\nambWIlzQzcxaxIgp6JJGSbpZ0uUNpk2R9Ns8/VZJh9VM+2y+9MDdkt4+EnJJ6pS0RtKC/Pc/w5xr\nR0lX50yzJE2umXaipHvy34kjKNdzNdur1IPrkhZJui3f97wG0yXpP/J+dKukN9RMG7LtVTBXldvr\n1ZJukPSMpM/UTRuyy4EUzNXrskOc6/j8/N0q6XpJr6uZVu72iogR8Qf8A/Bj4PIG06YDn8jDuwGL\naoZvATYFXgHcB4waAbk6gdsr3F4XASfm4YOBH+XhicD9+f+WeXjLqnPl26uGcFstAib1Mv0w4Jek\n71a8CZgzHNtrsLlGwPbaGngj8FXgMzXjR+XX4E7AJvm1uVvVufqz7BDn2q9nvyFdKqVn/yp9e42I\nFnpuqR0OnN1klgAm5OHNeeF896OACyPimYh4ALiXdEmCqnMNqX7k2g24Og//lhcuyfB24MqIWBYR\nTwJXAoeOgFxVOwr4YSQ3AltI2o4h3l4FclUqIpZExB+AdXWTKr0cSC+5KhUR1+f9B+BG0nd2YAi2\n14go6MBZwGnA802mnwGcIGkxcAVwch7f6PID24+AXACvyF0Pv5N0QImZ+pPrFuDdefgYYLykrah+\nezXLBTBG0jxJN0o6usRMkN54fyNpvtKlKOo12y5Dvb0Gmwuq3V7NVL29hmrZMu/7o6RPXTAE26vy\ngi7pCGBJRMzvZbbjgBkRMZn0MfRHkjain5cfqCDXo8CUiNiL3AUhaUIv91N2rs8AB0q6GTgQeBhY\nT/Xbq1kuSNurC/gAcJakncvIle0fEW8gfdz9O0lvro/fYJnoZXzVuaDa7dVM1dtrqJYt5b4lvYVU\n0P+pZ1SD2Qptr8oLOrA/cKSkRaSPHAdLOq9uno8CPwGIiBuAMaQL3Qzl5QcGnSt3AT2Rx88n9ZO9\narhyRcQjEfGu/Iby+TzuKSreXr3kIiIeyf/vB2YBe5WUq/a+lwCX8tJuuWbbZUgvb1EgV9Xbq5mq\nt9eQLFvGfUvak9QVeVRPbWAottdQHCQY7B9wEI0Ppv0SOCkPvyY/aAG78+KDovdT8kHRQebq6MlB\nOuDxMDBxGHNNAjbKw18FvpSHJwIPkA7wbZmHR0KuLYFNa+a5h5IOpgHjgPE1w9cDh9bNczgvPvg4\nd6i3V8FclW6vmnnP4MUHRUfn1+AreOEg3+4jIFe/lx2i53EK6fjefnXjS99ehR9QmX+1hQD4EnBk\nHt4NuC4/4AXAX9Us83lSC/hu4B0jIRepn/iOPP4m4J3DnOs9+UX+R1KrYNOaZT6Sd657gQ+PhFyk\nswBuy9vrNuCjJWbZKd/vLfk5+Xwe/3Hg43lYpB9ruS+vv2uot1eRXCNge21Lal2uAJbn4Ql52mH5\n+b2vZ9mqczVbdhhznQ08SaoRC4B5NcuXur381X8zsxYxEvrQzcysBC7oZmYtwgXdzKxFuKCbmbUI\nF3Qzsxbhgm4NSdpW0oWS7pN0p6QrJJX15aj6dR0tabd+zHdGz1X0lK7WOKQ/uFu7vn7O3ynp9ibT\nzu55jPnqfJPy8PU1y35gEBn3knR2Hn63pDskXdtzSQVJO0u6sGb+TSTNllTkB+JthHJBt5eQJNI3\n3mZFxM4RsRvwOWCbuvlGlbTKo0nn9A+74SpsEfHXEXFng/H75cFO0tf4B+pzwHfz8KmkLyD9sOa+\nvgL8c836niVdIO3YQazLRjgXdGvkLcC6iPjzddwjYkFEXCvpIKVrwP+Y9KUWJP2DpNvz3yl53DhJ\nMyXdkscfm8efmVv8t0r6lqT9gCOBb+brSe8s6WOS/pCX/amksc2CStpI0rmSvtJg2iJJ35A0N//t\nksfPkPRtSb8FviFpoqTLcqYb89e0e7xO0jVK10P/WF7+ZUrXdb9J6TrYtVfIG53z3Crp4p7szT5R\nSFqVB88EDsjb4NO5lf36mvmuq8uFpPHAnhFxSx71POlb02OBdUoXhXs0Iu6pW+1lwPHNtqltuPyx\nyxrZA+jtIlv7AHtExAOS9gY+DOxL+mbjHEm/I32D7pGIOBxA0uaSJpKusvjqiAhJW0TEcqUfaLg8\nIi7O8y6PiO/n4a+Qrpnz3foQpP33fNK157/aJOuKiNhH0odIV4M8Io9/FXBIRDwn6bvAzRFxtKSD\nSS3cnmK6J6nVOw64WdJMYAlwTESsyF0nN+qFH5nYlfTNzesknQP8LfCtXrZlj2mkr6sfkR/3MuAk\n4JTc1bVpRNxat0wXUNvF80XSj7Y/ApxAus7Q+xus63bSdcOtxbiFboMxN9L15wH+Erg0IlZHxCrg\nEuAAUuv9kNxCPiDSRbhWAGuBsyW9C3i6yf3vkVuot5Fakrs3me979F7MAS6o+f8XNeMviojnah7D\njwAi4hpgK0mb52k/i4g1EbGUdA33fUhvXF+TdCtwFemSpz3dUQ9FxHV5+Lx834NxEXCEpI1Jlx+Y\n0WCe7YDunhsRcWVE7B0R7yR1Y10B7Jo/KXy/59NCftzP5ha+tRAXdGvkDmDvXqavrhludAlQIuKP\n+T5uA74u6V8iYj2pIP6UVHB+1eT+ZwCfjIjXklqdY5rMdz3wFknNpsOLL0daO9zXY4i6/7Xjjydd\ngG3viHg98HhNxkbzD1hEPE36QY2jgPeRfgWq3hoabJtcuE8E/gv4OukNYT4v7mbZlPTmai3EBd0a\nuQbYtKfPGEDSGyUd2GDe2cDRksZKGkfqUrlW0suBpyPiPFKXwxskvQzYPCKuAE7hhW6NlUBta3E8\n8GhunfbW1/sDUiv0ol4Obh5VrDqoAAABpElEQVRb8/+GJvPM7lmPpIOApRGxIk87StKYfNbIQcAf\nSL9OtSQi1ild43rHmvuaIqnnk8BxwO97yV+rfhtAuqjTfwB/iIhlDZZZCOzSYPxpwHciYh2wGelN\n5XlS3zr5sXTn6dZC3IduL5H7t48h/XDCNFJLbhGpCG9fN+9NkmYAc/OosyPiZqUf7P6mpOdJPwn2\nCVLB+lluUQv4dF7mQuD7kv6edEXGfwbmAA+SWvhNuwYi4tu5e+RHko6PiPpfS9pU0hxS4+W4Jndz\nBvC/uQvlaVLrtsdcYCbpEqhfjohHJJ0P/ELpB4EXAHfVzL8QOFHS90hXlvzvZtnr3Aqsl3QL6UdT\n/j0i5ktaAfxvk8d+Vz42MT4iVgLkN9KuiDgjz/ZvpJ89W076VATpoPcV/cxlGxBfbdFaltKPbXTl\n/u8NTi7Os0gHkRv+rJ+kTwMrI6LZ77g2WuYS4LMRcXcpQW3EcJeL2QiUz8qZQ7pGdrPfaIX0CeCZ\nAdzvJsBlLuatyS10M7MW4Ra6mVmLcEE3M2sRLuhmZi3CBd3MrEW4oJuZtYj/D8zusdxaRDMBAAAA\nAElFTkSuQmCC\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x968f4a8>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.hist(P, range=(96*P_c, 104*P_c), bins=60, label='Estimator');\n",
"plt.axvline(P_c*100, ls='--', color='k', label='True value')\n",
"plt.xlabel('Crosstalk probability (%)')\n",
"plt.title('Simulation of crosstalk estimator accuracy');"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"> **NOTE** The crosstalk estimator is well centered around the true value."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python [conda env:py36]",
"language": "python",
"name": "conda-env-py36-py"
},
"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.6.0"
},
"nav_menu": {},
"toc": {
"navigate_menu": true,
"number_sections": true,
"sideBar": true,
"threshold": 6,
"toc_cell": false,
"toc_section_display": "block",
"toc_window_display": false
}
},
"nbformat": 4,
"nbformat_minor": 1
}