{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "1c7d465d", "metadata": { "execution": { "iopub.execute_input": "2026-06-19T07:57:01.838570Z", "iopub.status.busy": "2026-06-19T07:57:01.838393Z", "iopub.status.idle": "2026-06-19T07:57:04.871688Z", "shell.execute_reply": "2026-06-19T07:57:04.870403Z" }, "tags": [ "remove-cell" ] }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "An NVIDIA GPU may be present on this machine, but a CUDA-enabled jaxlib is not installed. Falling back to cpu.\n" ] } ], "source": [ "%matplotlib inline\n", "import brainmass\n", "import brainstate\n", "import brainunit as u\n", "import jax.numpy as jnp\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "brainstate.random.seed(0)\n", "brainstate.environ.set(dt=0.1 * u.ms)" ] }, { "cell_type": "markdown", "id": "25ceb375", "metadata": {}, "source": [ "# From Activity to Signals\n", "\n", "Neural mass models produce *neural activity* — firing rates, membrane potentials,\n", "synaptic currents. But experiments do not measure activity directly; they measure\n", "**signals**: a BOLD time course from fMRI, voltages from EEG, magnetic fields from\n", "MEG. A **forward model** is the biophysics that turns simulated activity into a\n", "predicted signal, so a model can be compared to (and fit against) real data.\n", "\n", "This page explains the theory of the two forward models ``brainmass`` ships — the\n", "hemodynamic (BOLD) model and the electromagnetic (EEG/MEG) lead-field model — what\n", "each observable measures, and the assumptions baked into each. The hands-on version\n", "is {doc}`/tutorials/05_forward_models`; the architectural role of this layer is in\n", "{doc}`/concepts/architecture_overview`." ] }, { "cell_type": "markdown", "id": "02fe290a", "metadata": {}, "source": [ "## The forward-model idea\n", "\n", "Every forward model is a map\n", "\n", "$$\n", "\\text{neural activity } a(t) \\;\\xrightarrow{\\;\\text{biophysics}\\;}\\; \\text{signal } y(t).\n", "$$\n", "\n", "What distinguishes the modalities is *what physics does the mapping* and *what is\n", "actually measured*:\n", "\n", "| Modality | Measures | Forward physics | Timescale |\n", "|---|---|---|---|\n", "| **fMRI / BOLD** | blood-oxygenation changes | neurovascular coupling (hemodynamics) | seconds (slow) |\n", "| **EEG** | scalp electric potential | volume conduction of dipolar currents | milliseconds (fast) |\n", "| **MEG** | extracranial magnetic field | magnetic field of the same currents | milliseconds (fast) |\n", "\n", "In ``brainmass`` every forward model is itself a **differentiable** module, so a\n", "loss defined in *signal space* (BOLD, EEG, or MEG) backpropagates all the way to\n", "the underlying neural parameters — the end-to-end differentiability discussed in\n", "{doc}`/concepts/why_differentiable`." ] }, { "cell_type": "markdown", "id": "1c841af3", "metadata": {}, "source": [ "## BOLD: the hemodynamic forward model\n", "\n", "fMRI does not see neurons; it sees **blood**. Local neural activity raises metabolic\n", "demand, which (via *neurovascular coupling*) increases blood flow, blood volume,\n", "and the local ratio of oxygenated to deoxygenated haemoglobin. Because oxygenated\n", "and deoxygenated blood have different magnetic susceptibilities, this changes the\n", "MR signal — the **blood-oxygen-level-dependent (BOLD)** contrast. The forward model\n", "is the chain *activity → flow → volume & deoxyhaemoglobin → BOLD*.\n", "\n", "### The Balloon–Windkessel model\n", "\n", "``brainmass.BOLDSignal`` implements the Balloon–Windkessel model of Friston et al.\n", "(2003). For each region, synaptic activity $z_i$ drives a vasodilatory signal\n", "$x_i$, which changes inflow $f_i$, blood volume $v_i$, and deoxyhaemoglobin content\n", "$q_i$:\n", "\n", "$$\n", "\\begin{aligned}\n", "\\dot{x}_i &= z_i - \\kappa_i x_i - \\gamma_i (f_i - 1), \\\\\n", "\\dot{f}_i &= x_i, \\\\\n", "\\tau_i \\dot{v}_i &= f_i - v_i^{1/\\alpha}, \\\\\n", "\\tau_i \\dot{q}_i &= \\frac{f_i}{\\rho}\\bigl[1 - (1-\\rho)^{1/f_i}\\bigr] - q_i\\,v_i^{\\,1/\\alpha - 1},\n", "\\end{aligned}\n", "$$\n", "\n", "with the measured BOLD a static nonlinear read-out of volume and deoxyhaemoglobin,\n", "\n", "$$\n", "\\mathrm{BOLD}_i = V_0\\bigl[k_1(1 - q_i) + k_2(1 - q_i/v_i) + k_3(1 - v_i)\\bigr],\n", "$$\n", "\n", "where $\\rho$ is the resting oxygen-extraction fraction, $\\alpha$ is Grubb's\n", "exponent, $\\tau$ the hemodynamic transit time, and $k_1, k_2, k_3$ are fixed\n", "biophysical constants. The key assumptions: the response is **regional and local**\n", "(no hemodynamic coupling between regions), and it is a **slow temporal low-pass** —\n", "the haemodynamics blur fast neural dynamics into a response that peaks several\n", "seconds after activity.\n", "\n", "### The convolution (HRF) shortcut\n", "\n", "The full ODE is detailed but stiff and slow. For fitting, ``brainmass`` also offers\n", "``HRFBold``: a **single linear convolution** of (downsampled) activity with a\n", "closed-form **hemodynamic response function** $h(t)$,\n", "\n", "$$\n", "\\mathrm{BOLD}(t) = k_1 V_0\\,\\bigl(h * y_{\\downarrow}\\bigr)(t) - 1,\n", "$$\n", "\n", "using kernels such as the canonical first-order Volterra (damped oscillator) or a\n", "mixture of gammas. This is fast, simple, and **differentiable in its handful of\n", "scalar parameters** — the right tool when the BOLD model sits inside an\n", "optimisation loop. Empirically the two routes track each other closely (best-lag\n", "correlation $\\approx 0.98$ on a slow drive), differing mainly in their intrinsic\n", "latency. So: ``BOLDSignal`` when biophysical realism matters, ``HRFBold`` when you\n", "are fitting." ] }, { "cell_type": "code", "execution_count": 2, "id": "328b5eb3", "metadata": { "execution": { "iopub.execute_input": "2026-06-19T07:57:04.874123Z", "iopub.status.busy": "2026-06-19T07:57:04.873679Z", "iopub.status.idle": "2026-06-19T07:57:05.307334Z", "shell.execute_reply": "2026-06-19T07:57:05.306002Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAArIAAAGGCAYAAACHemKmAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAgbpJREFUeJzt3Xd4U+X7BvA7SZvuQfegdLIKtIUCZW8pyEYEAVkCIoqAfHHgAHHhAhRRUBRBf2xEFFRWZe9VNoWWLtrSSfdu3t8fpZHQAm2acpr2/lxXLpqTc07urPL0zXueIxNCCBARERER6Rm51AGIiIiIiLTBQpaIiIiI9BILWSIiIiLSSyxkiYiIiEgvsZAlIiIiIr3EQpaIiIiI9BILWSIiIiLSSyxkiYiIiEgvsZAlIiIiIr3EQpaoBkycOBEeHh41tv8ePXqgR48eNbb/ukzb5+7AgQOQyWQ4cOCAzjM9SCaT4f3336/2ft5//33IZDKkpKQ8dl0PDw9MnDix2vf5pJQ9NqlV57OuD5/jyrwXyz4bW7dufTKhiO7DQpZqjYiICEybNg1eXl4wNjaGpaUlOnfujK+//hp5eXlSx9NLHh4eGDhwYIW3VfSfz5o1ayCTydQXAwMDuLq6YuLEiYiLiyu3jx49emisf//l+vXrNfa4Hufq1at4//33ERUVVaP3s379enz11Vc1eh9UfSqVCl9++SUaN24MExMTeHt7Y/r06cjOzpY6Gj0BD/5eu/9y584dqeNRNRlIHYAIAP766y88++yzMDIywvjx49GyZUsUFhbiyJEjeP3113HlyhX88MMPUsestFWrVkGlUkkdQ2sffPABPD09kZ+fjxMnTmDNmjU4cuQILl++DGNjY411GzZsiEWLFpXbh4uLy5OKW87Vq1excOFC9OjRo9xo2Z49e7TaZ7du3ZCXlwelUqletn79ely+fBmzZ8+uRtry8vLyYGDAX8+68vXXX+P111/H0KFD8frrryM6OhobNmzAm2++CXNzc6nj0RNS9nvtftbW1tKEIZ3hb0qSXGRkJJ577jm4u7vj33//hbOzs/q2V155BeHh4fjrr78kTFh1hoaGUkeolv79+6Nt27YAgClTpsDOzg6fffYZ/vzzT4wcOVJjXSsrKzz//PNSxNTK/YVoVcjl8nJFfE15UvdTU4QQyM/Ph4mJSbnb8vPzoVQqIZc/uS8EN27ciBYtWmDbtm3q6QgffvihXv+xSVV3/+81qjs4tYAk9/nnnyM7Oxs//fSTRhFbxsfHB7NmzVJf//nnn9GrVy84ODjAyMgIvr6+WLFiRbntyr5WP3LkCNq3bw9jY2N4eXnhl19+KbfurVu38Oyzz8LGxgampqbo0KFDueK57Kv4zZs34+OPP0bDhg1hbGyM3r17Izw8XGPdiubNqVQqfP3112jVqhWMjY1hb2+Pfv364cyZM1V+bE9a165dAZRO/6iKoqIiXL9+HQkJCY9d9+LFi5g4caJ6aomTkxNeeOEFpKamlls3Li4OkydPhouLC4yMjODp6Ynp06ejsLAQa9aswbPPPgsA6Nmzp/orxLK5rffPS0xMTISBgQEWLlxY7j7CwsIgk8mwfPlyAOXnyPbo0QN//fUXoqOj1ffh4eGB7OxsmJmZabxny9y+fRsKhaLCEez7PTgvMSsrC7Nnz4aHhweMjIzg4OCAp556CufOnXvc0woASElJwciRI2FpaQlbW1vMmjUL+fn5j9zmYXNQy76mvX/aRtlnbffu3Wjbti1MTEzw/fffq5+zjRs34t1334WrqytMTU2RmZkJADh58iT69esHKysrmJqaonv37jh69Gi5+zxy5AjatWsHY2NjeHt74/vvv6/U4y4jl8uhUqk0Ho9cLq/yqPf27dvRsmVLGBsbo2XLlvj9998rXE+lUuGrr75CixYtYGxsDEdHR0ybNg1379595P4LCwsxf/58BAYGwsrKCmZmZujatSv279+vXkcIAQ8PDwwZMqTc9vn5+bCyssK0adPUywoKCrBgwQL4+PjAyMgIbm5ueOONN1BQUKCxbUFBAV577TXY29vDwsICgwcPxu3bt6vy9KCkpARvv/02nJycYGZmhsGDByM2NlZ9+4IFC2BoaIjk5ORy27744ouwtrZ+7PuyurKyslBSUlKlbc6cOYPg4GDY2dnBxMQEnp6eeOGFF2ooIVUVR2RJcjt27ICXlxc6depUqfVXrFiBFi1aYPDgwTAwMMCOHTvw8ssvQ6VS4ZVXXtFYNzw8HCNGjMDkyZMxYcIErF69GhMnTkRgYCBatGgBoLSY6dSpE3JzczFz5kzY2tpi7dq1GDx4MLZu3Yphw4Zp7PPTTz+FXC7H3LlzkZGRgc8//xxjx47FyZMnH5l78uTJWLNmDfr3748pU6aguLgYhw8fxokTJ9SjBFV5bJVVVFRU4cE+GRkZld5HWdHSoEGDcreVlJSU27+xsTHMzc0RFxeH5s2bY8KECVizZs0j72Pv3r24desWJk2aBCcnJ/V0kitXruDEiRPqIiQ+Ph7t27dHeno6XnzxRTRr1gxxcXHYunUrcnNz0a1bN8ycORPLli3D22+/jebNmwOA+t/7OTo6onv37ti8eTMWLFigcdumTZugUCjURfGD3nnnHWRkZOD27dtYunQpAMDc3Bzm5uYYNmwYNm3ahCVLlkChUKi32bBhA4QQGDt27COfiwe99NJL2Lp1K2bMmAFfX1+kpqbiyJEjuHbtGtq0afPY7UeOHAkPDw8sWrQIJ06cwLJly3D37t0K/6jTVlhYGEaPHo1p06Zh6tSpaNq0qfq2Dz/8EEqlEnPnzkVBQQGUSiX+/fdf9O/fH4GBgViwYAHkcrn6D7nDhw+jffv2AIBLly6hb9++sLe3x/vvv4/i4mIsWLAAjo6Olc42adIkTJs2Dd9//71GkVcVe/bswTPPPANfX18sWrQIqampmDRpEho2bFhu3WnTpmHNmjWYNGkSZs6cicjISCxfvhznz5/H0aNHH/qNTWZmJn788UeMHj0aU6dORVZWFn766ScEBwfj1KlTCAgIgEwmw/PPP4/PP/8caWlpsLGxUW+/Y8cOZGZmqr8hUalUGDx4MI4cOYIXX3wRzZs3x6VLl7B06VLcuHED27dvV287ZcoU/N///R/GjBmDTp064d9//8WAAQOq9Bx9/PHHkMlkePPNN5GUlISvvvoKffr0QWhoKExMTDBu3Dh88MEH2LRpE2bMmKHerrCwEFu3bsUzzzxTo99G9OzZE9nZ2VAqlQgODsbixYvRuHHjR26TlJSkfv+99dZbsLa2RlRUFLZt21ZjOamKBJGEMjIyBAAxZMiQSm+Tm5tbbllwcLDw8vLSWObu7i4AiEOHDqmXJSUlCSMjI/G///1PvWz27NkCgDh8+LB6WVZWlvD09BQeHh6ipKRECCHE/v37BQDRvHlzUVBQoF7366+/FgDEpUuX1MsmTJgg3N3d1df//fdfAUDMnDmzXHaVSlXlx9a9e3fRvXv3cus+qOw5eNRly5Yt6vV//vlnAUDs27dPJCcni9jYWLF161Zhb28vjIyMRGxsbLkcFe1zwoQJQgghIiMjNa4/SkWPfcOGDeVew/Hjxwu5XC5Onz5dbv2y53LLli0CgNi/f3+5dR587r7//vtyr58QQvj6+opevXqpr5e9/vfvc8CAARqvc5ndu3cLAOKff/7RWO7n51ep1w2AWLBggfq6lZWVeOWVVx673YMWLFggAIjBgwdrLH/55ZcFAHHhwgX1Mnd3d43XqWzbB5W9RyIjIzW2BSB27dqlsW7Zc+bl5aXx+qpUKtG4cWMRHBxc7v3v6ekpnnrqKfWyoUOHCmNjYxEdHa1edvXqVaFQKCrMV5G33npLKJVKoVAoxLZt2yq1zYMCAgKEs7OzSE9PVy/bs2ePAKDxHjh8+LAAINatW6ex/a5du8otf/C9WFxcrPG7RQgh7t69KxwdHcULL7ygXhYWFiYAiBUrVmisO3jwYOHh4aF+Tn/99Vchl8s1frcJIcTKlSsFAHH06FEhhBChoaECgHj55Zc11hszZky592JFyl5nV1dXkZmZqV6+efNmAUB8/fXX6mUdO3YUQUFBGttv27btoZ/Xh8nMzFT/bn5QTk6OKCoqUl/ftGmTmDhxoli7dq34/fffxbvvvitMTU2FnZ2diImJeeT9/P777wJAhb9vqHbg1AKSVNlXjBYWFpXe5v55dxkZGUhJSUH37t1x69atcqOMvr6+6q/FAcDe3h5NmzbFrVu31Mv+/vtvtG/fHl26dFEvMzc3x4svvoioqChcvXpVY5+TJk3SmGdZtv/79/mg3377DTKZrNyoHwCNrzur8tgqKygoCHv37i13+fLLLx+6TZ8+fWBvbw83NzeMGDECZmZm+PPPPyscffLw8Ci37zfeeEN9mxDisaOxgOZjz8/PR0pKCjp06AAA6q/QVSoVtm/fjkGDBlU4102bdkzDhw+HgYEBNm3apF52+fJlXL16FaNGjary/oDS58/FxQXr1q3T2OfFixe1mk9sbW2NkydPIj4+Xqs8D47mv/rqqwBK3/u64unpieDg4ApvmzBhgsbrGxoaips3b2LMmDFITU1FSkoKUlJSkJOTg969e+PQoUNQqVQoKSnB7t27MXToUDRq1Ei9ffPmzR96Xw9atmwZlixZgqNHj2L06NF47rnnyh3wZ2RkhPfee++h+0hISEBoaCgmTJgAKysr9fKnnnoKvr6+Gutu2bIFVlZWeOqpp9SPKyUlBYGBgTA3N9eYJvAghUKh/t2iUqmQlpaG4uJitG3bVmMaSZMmTRAUFKTx/kpLS8M///yDsWPHqj8HW7ZsQfPmzdGsWTONLL169QIAdZay98HMmTM18lT1IMbx48dr/C4fMWIEnJ2dNd5n48ePx8mTJzWmKa1btw5ubm7o3r37I/dfUFCATz/9FI0aNYKlpSVMTEzw1FNPYfny5bh06RJu3bqF1atXo1WrVhodKUaOHImff/4Z48ePx9ChQ/Hhhx9i9+7dSE1Nxccff/zI+yw7GGznzp0oKiqqytNBTwinFpCkLC0tAZTOW6qso0ePYsGCBTh+/Dhyc3M1bsvIyND4j+b+//zKNGjQQGOuWnR0NIKCgsqtV/ZVdHR0NFq2bPnQfZZ93f6o+W8RERFwcXHR+BqwIlV5bJVlZ2eHPn36lFv+qPmB3377LZo0aYKMjAysXr0ahw4dgpGRUYXrmpmZVbj/qkpLS8PChQuxceNGJCUladxWVsQnJycjMzNT4/WoLjs7O/Tu3RubN2/Ghx9+CKB0WoGBgQGGDx+u1T7lcjnGjh2LFStWIDc3F6ampli3bh2MjY0fOlXhUT7//HNMmDABbm5uCAwMxNNPP43x48fDy8urUts/+PWpt7c35HK5TtuTPXg0+KNuu3nzJoDSAvdhMjIyUFBQgLy8vAq//m3atOljC/G8vDwsWLAAU6ZMQdu2bfHzzz8jJSUFw4YNw+7du9GlSxfcvHkThYWFFf4OKBMdHQ2g/PNYluP+IvPmzZvIyMiAg4NDhft68L39oLVr12Lx4sW4fv26RuH04HM4fvx4zJgxA9HR0XB3d8eWLVtQVFSEcePGaWS5du0a7O3tH5klOjoacrkc3t7e5R5bVTz4/MhkMvj4+Gi8z0aNGoXZs2dj3bp1mD9/PjIyMrBz50689tprj/1DdOPGjfjhhx8wa9YsNGvWDFFRUfjrr7/UU1aA0sLz9ddfV//f8jBdunRBUFAQ9u3b98j1unfvjmeeeQYLFy7E0qVL0aNHDwwdOhRjxox56O9EerJYyJKkLC0t4eLigsuXL1dq/YiICPTu3RvNmjXDkiVL4ObmBqVSib///htLly4tdxTy/fMT7yeE0DpzTewTqPpjq0nt27dXj3gOHToUXbp0wZgxYxAWFlZj7YpGjhyJY8eO4fXXX0dAQADMzc2hUqnQr1+/Gn/szz33HCZNmoTQ0FAEBARg8+bN6N27N+zs7LTe5/jx4/HFF19g+/btGD16NNavX4+BAwdq9cfIyJEj0bVrV/z+++/Ys2cPvvjiC3z22WfYtm0b+vfvX+X9VWbk+mHrPOxAmYo6FDzstrLX84svvkBAQECF25ibm5c7IKmqrl27hvT0dPXIvoGBAbZu3YpevXphwIAB2L9/PzZs2KA+eE4XVCoVHBwcNEZL7/ewohIA/u///g8TJ05UtwlzcHBQHxz44IGWzz33HF577TWsW7cOb7/9Nv7v//4Pbdu21Sg+VSoVWrVqhSVLllR4f25ublo8wupp0KABBg4cqC5kt27dioKCgkp9U9GtWzdcvXpVYx7tK6+8gpycHFy5cgUA4O/vX+kC083NDWFhYY9cp6zX9okTJ7Bjxw7s3r0bL7zwAhYvXowTJ06wfVstwEKWJDdw4ED88MMPOH78ODp27PjIdXfs2IGCggL8+eefGiOjj/q67nHc3d0r/GVW1tDf3d1d632X8fb2xu7du8sdnHG/mnhsulD2H2nPnj2xfPlyvPXWWzq/j7t37yIkJAQLFy7E/Pnz1cvLRu7K2Nvbw9LS8rF/+FR1isHQoUMxbdo09fSCGzduYN68eY/d7lH307JlS7Ru3Rrr1q1Dw4YNERMTg2+++aZKue7n7OyMl19+GS+//DKSkpLQpk0bfPzxx5UqZG/evKkxohceHg6VSvXIM1KVfdOQnp6u0WuzbHSyOspG/iwtLR85mm9vbw8TE5Ny7wMAjy1AgP9en/uPnDczM8Pff/+NLl26IDg4GPn5+fjoo48eWfyU/Q6oTA5vb2/s27cPnTt3fmRxX5GtW7fCy8tLo00YgAqnJNnY2GDAgAFYt24dxo4di6NHj5Y7OYe3tzcuXLiA3r17P/K96u7uDpVKhYiICI1CuDLP8f0efH6EEAgPD4efn5/G8vHjx2PIkCE4ffo01q1bh9atW6sPvn2Uh436m5mZqQ8OrIpbt2498g+L+3Xo0AEdOnTAxx9/jPXr12Ps2LHYuHEjpkyZUuX7Jd3iHFmS3BtvvAEzMzNMmTIFiYmJ5W6PiIjA119/DeC/0dD7Rz8zMjLw888/a33/Tz/9NE6dOoXjx4+rl+Xk5OCHH36Ah4dHuTlw2njmmWcghKiwzVPZY6mJx6YrPXr0QPv27fHVV19VqT1OZdtvVfTYAZT7j1kul2Po0KHYsWOHRtuyMmXbm5mZASgtwirD2toawcHB2Lx5MzZu3AilUomhQ4c+djszM7NHzl0eN24c9uzZg6+++gq2trZajZ6WlJSUuw8HBwe4uLhUesTy22+/1bheVlA/Kk9ZsXno0CH1spycHKxdu7ZS9/kogYGB8Pb2xpdfflnh2bXK2jMpFAoEBwdj+/btiImJUd9+7do17N69+7H306pVKzg6OmL58uUaX+nb2tqqpxnk5eVh0KBBj9yPs7MzAgICsHbtWo3XYu/eveXm0I8cORIlJSXqaSr3Ky4ufuR7sqLPwcmTJzV+N91v3LhxuHr1Kl5//XUoFAo899xz5bLExcVh1apV5bbNy8tDTk4OgP/eB8uWLdNYp6pnrfvll180polt3boVCQkJ5d5n/fv3V/emPnjwYI33oa6o3dfff/+Ns2fPol+/fo/c9u7du+V+L5V9i1DdbwxINzgiS5Lz9vbG+vXrMWrUKDRv3lzjzF7Hjh3Dli1b1OeA79u3L5RKJQYNGoRp06YhOzsbq1atgoODQ6V6lVbkrbfewoYNG9C/f3/MnDkTNjY2WLt2LSIjI/Hbb7/ppHF7z549MW7cOCxbtgw3b95Uf11++PBh9OzZEzNmzKiRx6ZLr7/+Op599lmsWbMGL730UqW2qWz7LUtLS3Tr1g2ff/45ioqK4Orqij179iAyMrLcup988gn27NmD7t27q1sKJSQkYMuWLThy5Aisra0REBAAhUKBzz77DBkZGTAyMlL3532YUaNG4fnnn8d3332H4ODgSp3xJzAwEJs2bcKcOXPQrl07mJubaxRFY8aMwRtvvIHff/8d06dP1+pEGVlZWWjYsCFGjBgBf39/mJubY9++fTh9+jQWL15cqX1ERkZi8ODB6NevH44fP65us+Tv7//Qbfr27YtGjRph8uTJ6kJp9erVsLe31ygqtSGXy/Hjjz+if//+aNGiBSZNmgRXV1fExcVh//79sLS0xI4dOwAACxcuxK5du9C1a1e8/PLLKC4uxjfffIMWLVrg4sWLj7wfAwMDLF++HKNGjUKrVq0wbdo0uLu749q1a+qDgm7fvo0hQ4bg6NGjj5xXuWjRIgwYMABdunTBCy+8gLS0NHWO+4vx7t27Y9q0aVi0aBFCQ0PRt29fGBoa4ubNm9iyZQu+/vprjBgxosL7GDhwILZt24Zhw4ZhwIABiIyMxMqVK+Hr61thwT9gwADY2tpiy5Yt6N+/f7n397hx47B582a89NJL2L9/Pzp37oySkhJcv34dmzdvVvf9DQgIwOjRo/Hdd98hIyMDnTp1QkhISLn+2I9jY2ODLl26YNKkSUhMTMRXX30FHx8fTJ06VWM9Q0NDPPfcc1i+fDkUCgVGjx5dpfupqk6dOqF169Zo27YtrKyscO7cOaxevRpubm54++23H7nt2rVr8d1332HYsGHw9vZGVlYWVq1aBUtLSzz99NM1mpsqSZpmCUTl3bhxQ0ydOlV4eHgIpVIpLCwsROfOncU333wj8vPz1ev9+eefws/PTxgbGwsPDw/x2WefidWrV1fYEmjAgAHl7qei1lURERFixIgRwtraWhgbG4v27duLnTt3aqxT1mLm/nZVQvzXYurnn39WL3uw/ZYQpa11vvjiC9GsWTOhVCqFvb296N+/vzh79myVH1tV2m9V9Bw87PGUtVaqqNVMSUmJ8Pb2Ft7e3qK4uFido0WLFg+9/6q037p9+7YYNmyYsLa2FlZWVuLZZ58V8fHxFbb/iY6OFuPHj1e3BfPy8hKvvPKKRuuiVatWCS8vL3WbprLWPg977jIzM4WJiYkAIP7v//6v3O0Vtd/Kzs4WY8aMEdbW1uXaMJV5+umnBQBx7Nixxz4HZe5/zAUFBeL1118X/v7+wsLCQpiZmQl/f3/x3XffPXY/ZS20rl69KkaMGCEsLCxEgwYNxIwZM0ReXp7Gug+23xJCiLNnz4qgoCChVCpFo0aNxJIlSx7afqui99nDPjNlzp8/L4YPHy5sbW2FkZGRcHd3FyNHjhQhISEa6x08eFAEBgYKpVIpvLy8xMqVKx/aHqwihw4dEsHBwcLS0lIYGRmJli1bikWLFonc3Fzxzz//CLlcLvr27avRsqkiv/32m2jevLkwMjISvr6+Ytu2bRV+1oUQ4ocffhCBgYHCxMREWFhYiFatWok33nhDxMfHq9d58L2oUqnEJ598Itzd3YWRkZFo3bq12Llz50PvQ4j/WqmtX7++wtsLCwvFZ599Jlq0aCGMjIxEgwYNRGBgoFi4cKHIyMhQr5eXlydmzpwpbG1thZmZmRg0aJCIjY2tUvutDRs2iHnz5gkHBwdhYmIiBgwYoNE27X6nTp0SAETfvn0fuW9deOedd0RAQICwsrIShoaGolGjRmL69Onizp07j9323LlzYvTo0aJRo0bCyMhIODg4iIEDB4ozZ87UeG6qHJkQ1TxChYiIHmrYsGG4dOlSlUe3iCrjtddew08//YQ7d+7A1NRU6jiVduHCBQQEBOCXX37R6LRAVFWcI0tEVEMSEhLw119/8T9qqhH5+fn4v//7PzzzzDN6VcQCwKpVq2Bubq51izuiMpwjS0SkY5GRkTh69Ch+/PFHGBoaan1aVKKKJCUlYd++fdi6dStSU1Mxa9YsqSNV2o4dO3D16lX88MMPmDFjhvrATCJtsZAlItKxgwcPYtKkSWjUqBHWrl0LJycnqSNRHXL16lWMHTsWDg4OWLZs2UN78dZGr776KhITE/H0009X2MWFqKo4R5aIiIiI9BLnyBIRERGRXmIhS0RERER6iXNkK0GlUiE+Ph4WFhZVPvUlEREREVWeEAJZWVlwcXF57EmJWMhWQnx8PNzc3KSOQURERFRvxMbGomHDho9ch4VsJVhYWAAofUIfdQpDIiIiIqqezMxMuLm5qeuvR2EhWwll0wksLS1ZyBIRERE9AZWZzsmDvYiIiIhIL7GQJSIiIiK9xEKWiIiIiPQSC1kiIiIi0kssZImIiIhIL7GQJSIiIiK9xEKWiIiIiPQSC1kiIiIi0kssZImIiIhIL7GQJSIiIiK9xFPUEmkhNbsA/1y+g9NRaYhIzkZWfjFMDBVo2MAEAW7W6NvCCU0cH3+OaCIiItKeTAghpA5R22VmZsLKygoZGRmwtLSUOg5JKDYtF8tCbuL383EoVj36oxPgZo1Xe/mgVzOHSp0vmoiIiKpWd3FElqgSVCqBn49F4Yvd15FfpAIAtHK1Qp/mjmjhYgkrU0PkFpYgIikbxyJScPBGMkJj0zF57Rl0bWyHT4a1gpuNqcSPgoiIqG7hiGwlcES2fssuKMacTaHYczURABDkaYM3+zdDm0YNHrpNclYBfjx8Cz8fjUJhiQoWRgb49Bk/DPBzflKxiYiI9FJV6i4WspXAQrb+uptTiPGrT+FSXAaUCjnmD/LF2KBGlZ4qEJmSg9e3XMCZ6LsAgFm9G2N2n8acakBERPQQVam72LWA6CEycoswetUJXIrLgI2ZEpumdcDzHdyrVIR62plh44sdMK2bFwDg65CbeHf7ZageM7+WiIiIHo+FLFEF8otKMPWXM7h+Jwv2FkbY9GIHtH7EVIJHMVDIMe/p5lg0vBVkMmDdyRi8+8dl8MsQIiKi6mEhS/QAIQTe/v0STkWlwcLIAL9Obo/GOmilNbp9I3w1KgByGbD+ZAyW/xuug7RERET1FwtZogf838kYbDsXB7kMWDkuEM2cdDcvekiAKxYObgEAWLz3BracidXZvomIiOobFrJE9wm7k4UPd14FAMzr3xydfex0fh/jOnrgpe7epfex7RKOhafo/D6IiIjqAxayRPcUFqswa+N5FBar0KuZA6Z09ayx+3ojuCkG+7ugWCUwc2MokrLya+y+iIiI6ioWskT3rDp8C9fvZMHWTInPnvGr0RZZcrkMnz3jh6aOFkjJLsBrm0JRwk4GREREVcJClghAXHoevvn3JgDg3YHNYW9hVOP3aaJUYPmY1jAxVOBoeCpWHODBX0RERFXBQpYIwAc7riC/SIX2njYYGuD6xO63saMFPhhSevDXkr03cC7m7hO7byIiIn3HQpbqvf1hSdh9JREKuQwfDmn5xM+69WxbNwwNcIFKAG9uvYiC4pInev9ERET6ioUs1WtFJSp8sKO0S8GkTh5o6lT9frHaWDCoBezMlbiZlI1v90dIkoGIiEjfsJClem3buduITMmBrZkSs59qIlmOBmZKLBzcEgDw3f5wXEvIlCwLERGRvmAhS/VWYbEKy0JKD7Ca3sMb5kYGkuZ5upUTgls4olgl8OZvF9nFgIiI6DFYyFK9teVsLOLS82BvYYSxQe5Sx4FMVjpH18LYABdvZ2Azz/pFRET0SCxkqV4qKC7B8n9LR2Nf7uENE6VC4kSlHCyNMbtP6RSHL3aHISOvSOJEREREtRcLWaqXNp2ORUJGPpwsjTG6fSOp42gY39EdPg7mSMspxNf7bkodh4iIqNZiIUv1TnGJCt8fvAUAeKWnN4wNa8dobBlDhRzzB/oCANYej8LNxCyJExEREdVOLGSp3tlzNRFx6XmwMVPi2bZuUsepULcm9njK1xElKoEP/7omdRwiIqJaiYUs1Ts/HYkEAIwNalTrRmPv9+6A5jCQy3DoRjKORaRIHYeIiKjWYSFL9UpobDrORt+FoUKGcR2k71TwKO62ZhgTVDp/97NdYRCC7biIiIjux0KW6pWy0dhB/i5wsDSWOM3jzejlAxNDBS7EpmP3lUSp4xAREdUqLGSp3ohPz8PflxIAAJO7eEqcpnIcLIwxpWtp1i/3hKG4RCVxIiIiotqDhSzVG+tORqNEJdDBywYtXKykjlNpU7t5oYGpIcKTsrHtfJzUcYiIiGqNWlfIfvvtt/Dw8ICxsTGCgoJw6tSph6575coVPPPMM/Dw8IBMJsNXX31Vbp33338fMplM49KsWbMafARUGxWXqLDlzG0AwLgOHtKGqSJLY0NM7+ENAFj+bziKOCpLREQEoJYVsps2bcKcOXOwYMECnDt3Dv7+/ggODkZSUlKF6+fm5sLLywuffvopnJycHrrfFi1aICEhQX05cuRITT0EqqUOhCUjKasANmZKPOXrKHWcKnu+gztszZSIScvFdo7KEhERAahlheySJUswdepUTJo0Cb6+vli5ciVMTU2xevXqCtdv164dvvjiCzz33HMwMjJ66H4NDAzg5OSkvtjZ2dXUQ6BaauPpWADAM21coTSoVW/7SjFVGuDFbl4AgG/3h3OuLBEREWpRIVtYWIizZ8+iT58+6mVyuRx9+vTB8ePHq7XvmzdvwsXFBV5eXhg7dixiYmKqG5f0SGJmPvaHlY7qj2pXO0+AUBnPd3CHjZkSUam5+PNCvNRxiIiIJFdrCtmUlBSUlJTA0VHza19HR0fcuXNH6/0GBQVhzZo12LVrF1asWIHIyEh07doVWVkPP+1nQUEBMjMzNS6kv7aevY0SlUBb9wbwcbCQOo7WzIwM1B0Mlv8bjhIV+8oSEVH9VmsK2ZrSv39/PPvss/Dz80NwcDD+/vtvpKenY/PmzQ/dZtGiRbCyslJf3Nz0dxSvvhNCYPOZ0mkF+jwaW2Z8Rw9YmxriVkoO/rrXSoyIiKi+qjWFrJ2dHRQKBRITNZu+JyYmPvJArqqytrZGkyZNEB4e/tB15s2bh4yMDPUlNjZWZ/dPT9apyDREp+bC3MgAA/ycpY5TbeZGBpjYyQMA8MOhCJ7ti4iI6rVaU8gqlUoEBgYiJCREvUylUiEkJAQdO3bU2f1kZ2cjIiICzs4PL2qMjIxgaWmpcSH9tD20dC5p/5ZOMFUaSJxGN8Z39ICxoRyX4zJxPCJV6jhERESSqTWFLADMmTMHq1atwtq1a3Ht2jVMnz4dOTk5mDRpEgBg/PjxmDdvnnr9wsJChIaGIjQ0FIWFhYiLi0NoaKjGaOvcuXNx8OBBREVF4dixYxg2bBgUCgVGjx79xB8fPVmFxSr1mbyGtnaVOI3u2Jgp8Wxg6TSJ7w/dkjgNERGRdGrVENWoUaOQnJyM+fPn486dOwgICMCuXbvUB4DFxMRALv+v9o6Pj0fr1q3V17/88kt8+eWX6N69Ow4cOAAAuH37NkaPHo3U1FTY29ujS5cuOHHiBOzt7Z/oY6Mn79CNZGTkFcHBwggdvGyljqNTU7p6Yt3JaBy8kYxrCZlo7sxvDYiIqP6RCU6ye6zMzExYWVkhIyOD0wz0yIz157DzYgImd/HEewN9pY6jc6+sO4e/LiVgeGtXLBkVIHUcIiIinahK3VWrphYQ6Up2QTH2XSs9cHBIgIvEaWpG2QkS/rwQj/j0PInTEBERPXksZKlO2nv1DvKLVPC0M0MrVyup49QIfzdrdPCyQbFKYPWRSKnjEBERPXEsZKlO+uNet4LB/i6QyWQSp6k507p5AwA2nIpBRl6RxGmIiIieLBayVOfczSnE4ZspAIDBdXRaQZkeTe3RxNEcOYUl2HiKp14mIqL6pVqFbFFREWJjYxEWFoa0tDRdZSKqlr3XElGiEmjmZAFve3Op49QomUyGKV1K58r+cjyap60lIqJ6pcqFbFZWFlasWIHu3bvD0tISHh4eaN68Oezt7eHu7o6pU6fi9OnTNZGVqFJ2Xb4DAOjfUv/P5FUZgwNc0MDUEHHpeeoD3IiIiOqDKhWyS5YsgYeHB37++Wf06dMH27dvR2hoKG7cuIHjx49jwYIFKC4uRt++fdGvXz/cvHmzpnITVSgrvwhH7k0r6N9Kd6c2rs2MDRUY1a4RAGDtsShpwxARET1BVTohwunTp3Ho0CG0aNGiwtvbt2+PF154AStWrMCaNWtw+PBhNG7cWCdBiSrj3+tJKCxRwcvODI0d6va0gvuN6+iOHw5F4FhEKsLuZKGpk4XUkYiIiGpclUZkN2zYoC5iO3XqhMzMzArXMzY2xksvvYQXXnih+gmJqmD3ldJpBf1aOtXpbgUPcrU2QV/f0hHotcejpA1DRET0hGh9sNeJEyeQn59fbnlmZibefPPNaoUi0kZeYQn2X08GUFrI1jcTOnkAAH4/F4eMXLbiIiKiuq/KheyIESPw6aefQiaTISkpqdztOTk5+PLLL3USjqgqDt5IRl5RCVytTersSRAepYOXDZo5WSCvqASbz8RKHYeIiKjGVWmOLAA0atQIO3fuhBAC/v7+sLW1hb+/P/z9/REQEICwsDA4O9ePo8WpdimbVhDcon5NKygjk8kwoZMH5m27hF9OROGFLp5QyOvf80BERPVHlQvZJUuWAACUSiWOHj2K+Ph4nD9/HqGhofj999+hUqnw+eef6zwo0aMUl6jw7/XSbwiCWzhKnEY6QwNc8ek/1xGblod/ryfhKd/6+1wQEVHdV+VCtkxOTg4MDQ0BAEOGDNFZICJtnI2+i4y8IlibGiLQvYHUcSRjolRgVDs3/HDoFn49Ec1CloiI6jStD/YqK2KJaoOQe6OxPZs6wEBRv8+8PDaotKfs4ZvJiEnNlTgNERFRzanS//gxMVU7l3tcXFyV1ifSVtkZrXo3d5A4ifTcbc3QtbEdhAA2nK7aZ5aIiEifVKmQbdeuHaZNm/bIU9BmZGRg1apVaNmyJX777bdqByR6nMiUHNxKzoGBXIZuTeyljlMrjA1yBwBsPh2LwmKVxGmIiIhqRpXmyF69ehUff/wxnnrqKRgbGyMwMBAuLi4wNjbG3bt3cfXqVVy5cgVt2rTB559/jqeffrqmchOphdwbjW3vaQNLY055AUpHph0sjJCUVYDdV+5gkL+L1JGIiIh0rkojsra2tliyZAkSEhKwfPlyNG7cGCkpKbh58yYAYOzYsTh79iyOHz/OIpaemLJuBb2b88CmMoYKOZ5r5wYAWH+S0wuIiKhu0qprgYmJCUaMGIERI0boOg9RlWTmF+FUZBoAoA/nx2oY1b4Rlu8Px/FbqQhPyoaPg7nUkYiIiHSqfh/eTXrvYFgyilUC3vZmcLc1kzpOreJqbYJezUqL+w2nOCpLRER1j9aF7OnTp9G7d2/4+flh+PDh+OCDD/Dnn39WubMBUXXsD+O0gkcZc68V19azt5FfVCJxGiIiIt3SupAdN24cFAoFXnzxRXh6euLgwYOYNGkSPDw8YGtrq8uMRBVSqQQO3UgBAPRoym4FFenexAGu1ibIyCvC35cSpI5DRESkU1qf2Ss2NhZ//fUXvL29NZZHR0cjNDS0urmIHutqQiZSsgtgqlSgrbuN1HFqJYVchtHt3fDlnhtYfzIGw9s0lDoSERGRzmg9Itu5c2fcvn273HJ3d3eespaeiIM3kgEAnbztoDTgdO+HebatG+Qy4Ez0XYQnZUsdh4iISGeq9L//8OHD8f777+P333/HSy+9hA8//BB3796tqWxEj1RWyHbntIJHcrQ0Vh/0tflMrMRpiIiIdKdKUwu8vb1x9OhRfPfdd0hJKZ2b2KRJEwwZMgQdOnRA69at0apVKyiVyhoJS1QmM78I56JL/4jqwbN5PdbItm7Ydy0J287dxuvBTWGo4Ag2ERHpvyoVsl988YX657i4OISGhqovn332GW7dugUDAwM0bdoUFy9e1HlYojLHwlNRrBLwsjeDm42p1HFqvZ7NHGBnboSU7AKEXEtCv5ZOUkciIiKqNq0P9nJ1dYWrqysGDBigXpadnY3Q0FBcuHBBJ+GIHubgjdK2W905Glsphgo5RgQ2xMqDEdh0OoaFLBER1Qk6/X7R3NwcXbp0wSuvvKLL3RJpEELgYNi9+bEsZCttZNvSjgUHbyQjISNP4jRERETVx4lypHfCk7IRn5EPIwM5OnixZ3Fledmbo72nDVQC2HqmfMcRIiIifcNClvTO4ZulBxq297SBsaFC4jT6ZVRbNwDA5rOxUKmExGmIiIiqh4Us6Z2j4aWFbNfGdhIn0T9Pt3KGhZEBYtPycOJWqtRxiIiIqoWFLOmVohKVugDr5M1CtqpMlAoMDnABAGw8zZ6yRESk36pVyB4+fBjPP/88OnbsiLi4OADAr7/+iiNHjugkHNGDLsSmI6ewBA1MDeHrbCl1HL00ql3p9IJdV+4gI7dI4jRERETa07qQ/e233xAcHAwTExOcP38eBQUFAICMjAx88sknOgtIdL8j96YVdPKxg1wukziNfmrlaoXmzpYoLFZhe2ic1HGIiIi0pnUh+9FHH2HlypVYtWoVDA0N1cs7d+6Mc+fO6SQc0YOOhZdOK+jiw2kF2pLJZBh1rxXXxtOxEIIHfRERkX7SupANCwtDt27dyi23srJCenp6dTIRVSinoBjnYkpPS8tCtnqGtnaF0kCOawmZuByXKXUcIiIirWhdyDo5OSE8PLzc8iNHjsDLy6taoYgqcioyDcUqgUY2pjwtbTVZmyoR3KL07F5bz/KgLyIi0k9aF7JTp07FrFmzcPLkSchkMsTHx2PdunWYO3cupk+frsuMRAD+mx/bmaOxOjEisHR6wR8X4lFQXCJxGiIioqoz0HbDt956CyqVCr1790Zubi66desGIyMjzJ07F6+++qouMxIB+K9/LKcV6EYXHzs4WhohMbMA/15LQv9WzlJHIiIiqhKtR2RlMhneeecdpKWl4fLlyzhx4gSSk5Px4Ycf6jIfEQAgKSsf1+9kQSYDOnrztLS6oJDLMLxN6ajs1rM8ZS0REekfrQvZvLw85ObmQqlUwtfXF46Ojvjxxx+xZ88eXeYjAgAcjyjtVuDrbAkbM6XEaeqOsukFB24kIykrX+I0REREVaN1ITtkyBD88ssvAID09HQEBQVh8eLFGDJkCFasWKGzgEQAcOQmpxXUBG97c7RpZI0SlcAf5+OljkNERFQlWhey586dQ9euXQEAW7duhaOjI6Kjo/HLL79g2bJlOgtIJIRQz4/lgV66NyKw9ExfW86ypywREekXrQvZ3NxcWFhYAAD27NmD4cOHQy6Xo0OHDoiOjtZZQKKo1FzEZ+RDqZCjnYeN1HHqnAF+zjAykONGYjYuxWVIHYeIiKjStC5kfXx8sH37dsTGxmL37t3o27cvACApKQmWlpY6C0hU1narjbs1TJQKidPUPVYmhvf1lOVBX0REpD+0LmTnz5+PuXPnwsPDA0FBQejYsSOA0tHZ1q1b6ywg0YlbpQd6dfbmtIKaou4pG8qeskREpD+07iM7YsQIdOnSBQkJCfD391cv7927N4YNG6aTcERCCJy8lQYACPJi262a0tnHDk6WxriTmY+Qa0l4mj1liYhID2g9IguUnqa2devWkMv/20379u3RrFmzagcjAoCI5BykZBfAyEAOfzcrqePUWaU9ZV0BcHoBERHpD61HZAEgJCQEISEhSEpKgkql0rht9erV1QpGBPw3raBNowYwMuD82Jo0IrAhvjsQgYM3kpGUmQ8HS2OpIxERET2S1iOyCxcuRN++fRESEoKUlBTcvXtX40KkCycjS6cVdOC0ghrnZW+OQPcGKFEJbA+NkzoOERHRY2k9Irty5UqsWbMG48aN02UeIjUhhHpEtoMX2249CSMCG+Js9F1sPXsbU7t6QSaTSR2JiIjoobQekS0sLESnTp10mYVIw62UHCRnFUBpIIe/m7XUceoF9pQlIiJ9onUhO2XKFKxfv16XWYg0lHUraNPIGsaGnB/7JFgaG6Jfy9KeslvO8KAvIiKq3bQuZPPz87FkyRJ0794dr776KubMmaNx0da3334LDw8PGBsbIygoCKdOnXrouleuXMEzzzwDDw8PyGQyfPXVV9XeJ9Ue/00r4PzYJ6msp+yfF+KRX8SeskREVHtpXchevHgRAQEBkMvluHz5Ms6fP6++hIaGarXPTZs2Yc6cOViwYAHOnTsHf39/BAcHIykpqcL1c3Nz4eXlhU8//RROTk462SfVDkIInIwsLWSDPFnIPkmdvO3gbGWMjLwihFzj54SIiGovmRBCSB2iTFBQENq1a4fly5cDAFQqFdzc3PDqq6/irbfeeuS2Hh4emD17NmbPnq2zfZbJzMyElZUVMjIyePrdJyQyJQc9vzwApYEcFxf05dSCJ+yL3dfx7f4I9Gxqj58ntZc6DhER1SNVqbuqdUKE9PR0LF68GFOmTMGUKVOwdOlSZGRod4BIYWEhzp49iz59+vwXTi5Hnz59cPz48Se6z4KCAmRmZmpc6Mkqm1YQ4Mb5sVJ4pk3p9IKynrJERES1kdaF7JkzZ+Dt7Y2lS5ciLS0NaWlpWLJkCby9vXHu3Lkq7y8lJQUlJSVwdHTUWO7o6Ig7d+5olVHbfS5atAhWVlbqi5ubm1b3T9o7yfmxkirrKasSwLbz7ClLRES1k9aF7GuvvYbBgwcjKioK27Ztw7Zt2xAZGYmBAweW+3pf38ybNw8ZGRnqS2xsrNSR6pXS/rFlJ0Jg/1iplB309dvZ26hFM5CIiIjUqjUi++abb8LA4L9zKhgYGOCNN97AmTNnqrw/Ozs7KBQKJCYmaixPTEx86IFcNbVPIyMjWFpaalzoyYlOzcWdzHwoFXK0adRA6jj1VllP2ZtJ2bh4mz1liYio9tG6kLW0tERMTEy55bGxsbCwsKjy/pRKJQIDAxESEqJeplKpEBISgo4dO2qVsSb2STWvrFsB58dK6/6eslvPsqcsERHVPloXsqNGjcLkyZOxadMmxMbGIjY2Fhs3bsSUKVMwevRorfY5Z84crFq1CmvXrsW1a9cwffp05OTkYNKkSQCA8ePHY968eer1CwsLERoaitDQUBQWFiIuLg6hoaEIDw+v9D6p9uG0gtrj/p6yBcXsKUtERLWLweNXqdiXX34JmUyG8ePHo7i4GABgaGiI6dOn49NPP9Vqn6NGjUJycjLmz5+PO3fuICAgALt27VIfrBUTEwO5/L/aOz4+Hq1bt9bI9OWXX6J79+44cOBApfZJtYsQQn2gVxAP9JJcWU/ZhIx8hFxLwtOtnKWOREREpFbtPrK5ubmIiIgAAHh7e8PU1FQnwWoT9pF9cmJSc9Hti/0wVMhwcUEwTJScWiC1sp6yvZo5YPXEdlLHISKiOu6J9ZEFAFNTU7Rs2RItW7ask0UsPVll/WP9G1qziK0lhrOnLBER1VLVKmR/+ukntGzZEsbGxjA2NkbLli3x448/6iob1UMnItk/trbxtjdHm0bWKFEJbA9lT1kiIqo9tC5k58+fj1mzZmHQoEHYsmULtmzZgkGDBuG1117D/PnzdZmR6onS+bFlB3qxkK1NRgSWnhRkK3vKEhFRLaL1wV4rVqzAqlWrNDoUDB48GH5+fnj11VfxwQcf6CQg1R+37+YhLj0PBnIZ2rhbSx2H7jPAzxkLd1zBjcRsXIrLgF9Da6kjERERaT8iW1RUhLZt25ZbHhgYqO5iQFQV6vmxbtYwVWr9NxbVACsTQwS3KO0p+xt7yhIRUS2hdSE7btw4rFixotzyH374AWPHjq1WKKqf2D+2divrKfsHe8oSEVEtUa1hr59++gl79uxBhw4dAAAnT55ETEwMxo8fjzlz5qjXW7JkSfVSUr1QNiIb5Mn5sbVRZx87OFka405mPv69loT+7ClLREQS07qQvXz5Mtq0aQMA6j6ydnZ2sLOzw+XLl9XryWSyakak+iA2LVc9PzbQvYHUcagCCrkMw9q4YsWBCGw9e5uFLBERSU7rQnb//v26zEH13MnI0mkFrRpawcyI82Nrq2faNMSKAxE4cCMZSVn5cLAwljoSERHVY1rPkc3Ly0Nubq76enR0NL766ivs2bNHJ8Gofik7LS3bbtVuPg7maH2vp+wf5+OljkNERPWc1oXskCFD8MsvvwAA0tPT0b59eyxevBhDhgyp8CAwokfhiRD0R9lBX+wpS0REUtO6kD137hy6du0KANi6dSucnJwQHR2NX375BcuWLdNZQKr74tLzEJuWBwXnx+qFgX4uUBrIEZaYhSvxmVLHISKiekzrQjY3NxcWFhYAgD179mD48OGQy+Xo0KEDoqOjdRaQ6r6yaQWtXK1gzvmxtZ6ViSH6+joCKB2VJSIikorWhayPjw+2b9+O2NhY7N69G3379gUAJCUlwdLSUmcBqe47wfmxeqdsesH20Dj2lCUiIsloXcjOnz8fc+fOhYeHB4KCgtCxY0cApaOzrVu31llAqvvKToQQxBMh6I2uje3haGmE9Nwi7L+eJHUcIiKqp7QuZEeMGIGYmBicOXMGu3btUi/v3bs3li5dqpNwVPfFp+chJi0XCrkMbTk/Vm8o5DIMa/3fQV9ERERS0LqQBQAnJye0bt0acvl/u2nfvj2aNWtW7WBUP5y8162gpYslLIwNJU5DVTEi0BUAsD8sGclZBRKnISKi+qhahezhw4fx/PPPo2PHjoiLiwMA/Prrrzhy5IhOwlHdd/LetALOj9U/Pg4WCHC711M2NE7qOEREVA9pXcj+9ttvCA4OhomJCc6fP4+CgtIRmYyMDHzyySc6C0h1Gw/00m/PsKcsERFJSOtC9qOPPsLKlSuxatUqGBr+95Vw586dce7cOZ2Eo7rtTkY+olJzIZcBbT04P1YfDfZzgVIhx/U77ClLRERPntaFbFhYGLp161ZuuZWVFdLT06uTieoJ9fxYVyvOj9VTVqaGeKoFe8oSEZE0tC5knZycEB4eXm75kSNH4OXlVa1QVD9wWkHdUNZT9o/QOBQWqyROQ0RE9YnWhezUqVMxa9YsnDx5EjKZDPHx8Vi3bh3mzp2L6dOn6zIj1VHq/rGe7B+rz7r62MHBwgh3c4vwL3vKEhHRE6T1+UDfeustqFQq9O7dG7m5uejWrRuMjIwwd+5cvPrqq7rMSHVQYmY+IlNy7s2PZSGrzwwUcgxr44rvD97C1rO30a+lk9SRiIiontBqRLaoqAh9+vTByJEjkZaWhsuXL+PEiRNITk7Ghx9+qOuMVAeVTSvwdbGElQnnx+q7EW1KpxccCEtCSjZ7yhIR0ZOhVSFraGiIixcvAgCUSiV8fX3Rvn17mJub6zQc1V1l0wo6eHJ+bF3Q2NEC/g2tUKwS+CM0Xuo4RERUT2g9R/b555/HTz/9pMssVI+Ujch29GYhW1eUHfS15Uwse8oSEdETofUc2eLiYqxevRr79u1DYGAgzMzMNG5fsmRJtcNR3XQn47/5se14oFedMcjfBR/+dQ3X72ThclwmWjW0kjoSERHVcVoXspcvX0abNm0AADdu3NC4TSaTVS8V1Wllo7EtXa1gyf6xdYa1qRLBLZyw40I8Np2JQauGraSOREREdZzWhez+/ft1mYPqEfaPrbtGtXXDjgvx+CM0Hu887QsTpULqSEREVIdpPUeWSFvHy+bHspCtczp526JhAxNk5Rfjn8sJUschIqI6joUsPVHx6XmITs2FQi5DW48GUschHZPLZRjZ1g0AsOl0rMRpiIiormMhS0/U/fNjLTg/tk4aEdgQMhlwMjINUSk5UschIqI6jIUsPVHHI8rmx7JbQV3lYm2Cbo3tAQCbz3BUloiIag4LWXqiTkRyfmx98Fy70ukFW8/eRnGJSuI0RERUV2ldyC5atAirV68ut3z16tX47LPPqhWK6qbbd3MRm5Z3b34sR2Trst7NHWFrpkRSVgEOhCVLHYeIiOoorQvZ77//Hs2aNSu3vEWLFli5cmW1QlHdVHZaWr+GVjA30rrzG+kBpYEcw1q7AgA2cXoBERHVEK0L2Tt37sDZ2bnccnt7eyQksO0Olfff/FhOK6gPRt2bXvDv9SQkZeVLnIaIiOoirQtZNzc3HD16tNzyo0ePwsXFpVqhqG46wf6x9UpjRwu0aWSNEpXAtnNxUschIqI6SOvvd6dOnYrZs2ejqKgIvXr1AgCEhITgjTfewP/+9z+dBaS6ITYtF3HpeTCQyxDozv6x9cWodm44F5OOzadjMa2bF09fTUREOqV1Ifv6668jNTUVL7/8MgoLCwEAxsbGePPNNzFv3jydBaS6oexsXv5u1jDj/Nh6Y4CfCxbuuIpbKTk4HXUX7T15kB8REemO1lMLZDIZPvvsMyQnJ+PEiRO4cOEC0tLSMH/+fF3mozriBPvH1kvmRgYY6Fc6l55n+iIiIl2rdh9Zc3NztGvXDi1btoSRkZEuMlEdI4RQj8h29LKTOA09aaPaNQIA/H0pAZn5RRKnISKiuqRK3/HOmTMHH374IczMzDBnzpxHrrtkyZJqBaO6IyYtFwkZ+TBUcH5sfdSmkTV8HMwRnpSNHRfiMTbIXepIRERUR1SpkD1//jyKiorUPz8MD+ig+5W13Qpws4aJUiFxGnrSZDIZRrV1w8d/X8Om07EsZImISGeqVMju37+/wp+JHqWs7Rb7x9Zfw9u44vPd13HxdgYux2WgpauV1JGIiKgOqPYcWaJH0Zwfy0K2vrI1N0K/lqUHfa07GSNxGiIiqiuq1QcpPz8fFy9eRFJSElQqlcZtgwcPrlYwqhsiknOQmFkApUKO1o04P7Y+G9O+EXZciMefoXF4Z0BznqaYiIiqTev/SXbt2oXx48cjJSWl3G0ymQwlJSXVCkZ1w5GbyQCAth4NOD+2nuvgZQMvezPcSs7BH6FxnCtLRETVpvXUgldffRXPPvssEhISoFKpNC4sYqnMkfDSaQVdGrPtVn0nk8kwpn1pK651J2IghJA4ERER6TutC9nExETMmTMHjo6OusxDdUhRiUp9oFdXH3uJ01BtMCKwIZQGclxNyMSF2xlSxyEiIj2ndSE7YsQIHDhwQIdRqK65EJuO7IJiWJsawtfFUuo4VAtYmyoxoFXpQV/rT0ZLnIaIiPSd1nNkly9fjmeffRaHDx9Gq1atYGhoqHH7zJkzqx2O9NuR8NL505297aCQs7cwlRob1Ai/n4/Dnxfi8c4AX1iZGD5+IyIiogpoXchu2LABe/bsgbGxMQ4cOKBxEgSZTMZClnDkZmkhy/mxdL9A9wZo4miOG4nZ2H4+DhM6eUgdiYiI9JTWUwveeecdLFy4EBkZGYiKikJkZKT6cuvWrWqF+vbbb+Hh4QFjY2MEBQXh1KlTj1x/y5YtaNasGYyNjdGqVSv8/fffGrdPnDgRMplM49KvX79qZaRHy8ovwvnYdABAFx8WsvSf+w/6Wn+SB30REZH2tC5kCwsLMWrUKMjluj2nwqZNmzBnzhwsWLAA586dg7+/P4KDg5GUlFTh+seOHcPo0aMxefJknD9/HkOHDsXQoUNx+fJljfX69euHhIQE9WXDhg06zU2aTt5KQ4lKwN3WFG42plLHoVpmWJuGMDaUIywxC2ej70odh4iI9JTWVeiECROwadMmXWYBACxZsgRTp07FpEmT4Ovri5UrV8LU1BSrV6+ucP2vv/4a/fr1w+uvv47mzZvjww8/RJs2bbB8+XKN9YyMjODk5KS+NGjA5vw1qWx+LEdjqSJWJoYY5OcCoHRUloiISBtaz5EtKSnB559/jt27d8PPz6/cwV5Lliyp8j4LCwtx9uxZzJs3T71MLpejT58+OH78eIXbHD9+HHPmzNFYFhwcjO3bt2ssO3DgABwcHNCgQQP06tULH330EWxtecrUmsJClh5nTFAjbDl7GzsvJWD+IF9YmyqljkRERHpG60L20qVLaN26NQCU+xr//gO/qiIlJQUlJSXletM6Ojri+vXrFW5z586dCte/c+eO+nq/fv0wfPhweHp6IiIiAm+//Tb69++P48ePQ6Eof7apgoICFBQUqK9nZmZq9Xjqq4SMPIQnZUMmAzp5s5CligW4WcPX2RJXEzKx5cxtTO3mJXUkIiLSM1oXsvv379dljhr13HPPqX9u1aoV/Pz84O3tjQMHDqB3797l1l+0aBEWLlz4JCPWKWXdCvxcrWBlytZKVDGZTIZxHd0xb9sl/HoiGpO7eELONm1ERFQFWheyZa5evYqYmBgUFhaql8lkMgwaNKjK+7Kzs4NCoUBiYqLG8sTERDg5OVW4jZOTU5XWBwAvLy/Y2dkhPDy8wkJ23rx5GtMVMjMz4ebmVpWHUq8dDWfbLaqcIQEuWPT3NcSk5eLAjST0asYzBRIRUeVpXcjeunULw4YNw6VLlyCTydQtdMqmFZSUlFR5n0qlEoGBgQgJCcHQoUMBACqVCiEhIZgxY0aF23Ts2BEhISGYPXu2etnevXvRsWPHh97P7du3kZqaCmdn5wpvNzIygpGRUZXzEyCEwJHw0tPSduFpaekxTJUGGNnWDT8eicTaY9EsZImIqEq07lowa9YseHp6IikpCaamprhy5QoOHTqEtm3bVuvUtXPmzMGqVauwdu1aXLt2DdOnT0dOTg4mTZoEABg/frzGwWCzZs3Crl27sHjxYly/fh3vv/8+zpw5oy58s7Oz8frrr+PEiROIiopCSEgIhgwZAh8fHwQHB2udkyp2LSELKdkFMDFUoI27tdRxSA+M6+gOmQw4eCMZkSk5UschIiI9onUhe/z4cXzwwQews7ODXC6HXC5Hly5dsGjRomqd1WvUqFH48ssvMX/+fAQEBCA0NBS7du1SH9AVExODhIQE9fqdOnXC+vXr8cMPP8Df3x9bt27F9u3b0bJlSwCAQqHAxYsXMXjwYDRp0gSTJ09GYGAgDh8+zFHXGnDgRmm/307etjAyKH8gHdGD3G3N0LOpAwDg1+PREqchIiJ9IhNanlanQYMGOHfuHDw9PeHt7Y0ff/wRPXv2REREBFq1aoXc3FxdZ5VMZmYmrKyskJGRAUtLS6nj1GojVx7Hqag0fDikBcZ19JA6DumJA2FJmPjzaVgYGeDE271hZlTt6ftERKSnqlJ3aT0i27JlS1y4cAEAEBQUhM8//xxHjx7FBx98AC8vttGpjzLyinA2pvQsTT3ujbARVUa3xvbwtDNDVkExfj8fJ3UcIiLSE1oXsu+++y5UKhUA4IMPPkBkZCS6du2Kv//+G8uWLdNZQNIfR26moEQl4G1vxtPSUpXI5TKM6+AOAPjleBS0/KKIiIjqGa2/v7v/QCkfHx9cv34daWlpaNCggdYnRCD9diCsdH5sT47GkhaeCWyIL/eE4UZiNk7cSkNHb555j4iIHk3rEdmK2NjYsIitp1QqgQM3kgEAPZuxkKWqszIxxLDWrgCAtceipA1DRER6oVpHVISEhCAkJARJSUnqaQZlVq9eXa1gpF+uJmQiOasApkoF2no0kDoO6anxHT2w7mQM9ly9g7j0PLham0gdiYiIajGtR2QXLlyIvn37IiQkBCkpKbh7967GheqXsmkFnX3s2HaLtNbUyQIdvWyhEsD6k2zFRUREj6b1iOzKlSuxZs0ajBs3Tpd5SE/tDyudVtCjKc/mRdUzoZM7jt9KxYZTsXi1V2MYG/IPIyIiqpjWI7KFhYXo1KmTLrOQnkrNLsB5tt0iHenT3BGu1iZIyynEtnNsxUVERA+ndSE7ZcoUrF+/XpdZSE+FXE+CSgAtXCw5p5GqzUAhx6TOHgCAn47cgkrFVlxERFQxracW5Ofn44cffsC+ffvg5+cHQ0NDjduXLFlS7XCkH/ZcSQQA9PV1kjgJ1RWj2rnh6303EZGcg4M3ktkJg4iIKqR1IXvx4kUEBAQAAC5fvqxxG1tw1R95hSU4El46P/YpX0eJ01BdYWFsiOfau2HV4UisOnyLhSwREVVI60J2//79usxBeurwzWTkF6ngam2C5s4WUsehOmRiZ0+sPhqFYxGpuBKfgRYuVlJHIiKiWkbrObITJkzAoUOHdJmF9NCeq/emFbRw5Eg86ZSrtQmebuUMAPjpSKTEaYiIqDbSupDNyMhAnz590LhxY3zyySeIi+PRxfVNiUrg3+ul/WM5rYBqwtSungCAHRfikZiZL3EaIiKqbbQuZLdv3464uDhMnz4dmzZtgoeHB/r374+tW7eiqKhIlxmpljobfRdpOYWwMjFEew8bqeNQHeTX0BrtPWxQVCJ42loiIipH60IWAOzt7TFnzhxcuHABJ0+ehI+PD8aNGwcXFxe89tpruHnzpq5yUi2058odAEDvZg4wUFTrrUT0UJPvjcquOxmD3MJiidMQEVFtopPqIyEhAXv37sXevXuhUCjw9NNP49KlS/D19cXSpUt1cRdUywgh8M/l0kK2bwtOK6Ca06e5IzxsTZGRV4StZ29LHYeIiGoRrQvZoqIi/Pbbbxg4cCDc3d2xZcsWzJ49G/Hx8Vi7di327duHzZs344MPPtBlXqolzsWkIy49D2ZKBc/mRTVKIZfhhS6lo7Krj0SihCdIICKie7Ruv+Xs7AyVSoXRo0fj1KlT6p6y9+vZsyesra2rEY9qq50X4wGUHuRlbKiQOA3VdSMCG2LxnhuISs3F3quJ6NeSJ98gIqJqjMguXboU8fHx+PbbbyssYgHA2toakZFsm1PXqFQCf19KAAAM9HOROA3VB6ZKAzzfoREAYMXBCAjBUVkiIqpGITtu3DgYGxvrMgvpidNRaUjMLICFsQG6NrGTOg7VE5M6e8LIQI4Lsek4HpEqdRwiIqoFqjS1YM6cOZVed8mSJVUOQ/ph58XS0djgFk4wMuC0Anoy7MyNMKqdG345Ho3vDkSgkw//iCIiqu+qVMieP39e4/q5c+dQXFyMpk2bAgBu3LgBhUKBwMBA3SWkWqW4RIV/LpdNK3CWOA3VN1O7emHdyRgcCU/Bxdvp8GtoLXUkIiKSUJUK2f3796t/XrJkCSwsLLB27Vo0aNAAAHD37l1MmjQJXbt21W1KqjVORqYhJbsQDUwN0ZkjYvSEudmYYoi/C7adj8N3+yOwchz/aCYiqs+0niO7ePFiLFq0SF3EAkCDBg3w0UcfYfHixToJR7XPb+dK+3j2a+kMQ54EgSTwUg9vAMDuq3cQnpQtcRoiIpKS1pVIZmYmkpOTyy1PTk5GVlZWtUJR7ZRdUIx/LpWeBGFEYEOJ01B91cTRAk/5OkIIYOXBCKnjEBGRhLQuZIcNG4ZJkyZh27ZtuH37Nm7fvo3ffvsNkydPxvDhw3WZkWqJvy7GI6+oBF72ZmjTyFrqOFSPvXxvVHb7+TjEpedJnIaIiKSidSG7cuVK9O/fH2PGjIG7uzvc3d0xZswY9OvXD999950uM1ItseVM6bSCZwPdIJPJJE5D9VnrRg3Q0csWxSqBVYduSR2HiIgkonUha2pqiu+++w6pqak4f/48zp8/j7S0NHz33XcwMzPTZUaqBSKSs3Em+i7kMmB4G1ep4xDh1V4+AID1p2KQmJkvcRoiIpJCtY/WMTMzg5+fH/z8/FjA1mFbz5aOxnZvYg9HS54Ig6TX0dsW7TwaoLBYhRUHOFeWiKg+qlYhm56ejsWLF2PKlCmYMmUKlixZgoyMDF1lo1qiRCWw7V63gmfbukmchqiUTCbD7D5NAHBUloiovtK6kD1z5gy8vb2xdOlSpKWlIS0tDUuXLoW3tzfOnTuny4wksZBriUjMLIC1qSF6N3eQOg6RWqf7RmXZwYCIqP7RupB97bXXMHjwYERFRWHbtm3Ytm0bIiMjMXDgQMyePVuHEUlqPx+NAgCMaufGU9JSrSKTyTCr971R2ZMxSOKoLBFRvVKtEdk333wTBgb/nRzMwMAAb7zxBs6cOaOTcCS9awmZOH4rFQq5DOM7ekgdh6iczj62aOveAAXFKqzgqCwRUb2idSFraWmJmJiYcstjY2NhYWFRrVBUe6y5Nxob3MIRrtYm0oYhqoDGXFmOyhIR1StaF7KjRo3C5MmTsWnTJsTGxiI2NhYbN27ElClTMHr0aF1mJImk5RRie2gcAGBSZ0+J0xA9XGcfWwRyVJaIqN4xePwqFfvyyy8hk8kwfvx4FBcXQwgBpVKJ6dOn49NPP9VlRpLIhlMxKChWoaWrJdq6N5A6DtFDlY7KNsa4n05h3YkYTOnqxW8QiIjqAa1HZJVKJb7++mvcvXsXoaGhuHDhgrpzgZGRkS4zkgTyi0rwy/EoAMCkTp48kxfVel187NDRyxaFJSp8ve+G1HGIiOgJ0HpEFgDy8/Nx+fJlJCUlQaVSITIyUn3b4MGDqx2OpLPpdCwSMwvgbGWMgf7OUscheiyZTIbX+zXF8O+OYevZ23ixmxd8HDhfn4ioLtO6kN21axfGjRuH1NTUcrfJZDKUlJRUKxhJJ7+oBN/uDwcAvNLThy23SG+0adQAT/k6Yu/VRCzecwMrng+UOhIREdUgracWvPrqqxg5ciQSEhKgUqk0Lixi9du6kzFIyiqAq7UJRvJMXqRn5vZtCpkM+OfyHVy8nS51HCIiqkFaF7KJiYmYM2cOHB0ddZmHJJZbWIwVB0pHY1/t5QOlQbXOYkz0xDV1ssCw1q4AgC92h0mchoiIapLWVcqIESNw4MABHUah2uDX49FIyS6Em40JnglsKHUcIq281qcJDBUyHL6ZgmPhKVLHISKiGqL1HNnly5fj2WefxeHDh9GqVSsYGhpq3D5z5sxqh6MnKzW7AMvvzY2d2asxDBUcjSX95GZjijHtG2Ht8Wh8tus6fn+5M+Rydt4gIqprtC5kN2zYgD179sDY2BgHDhzQaM8kk8lYyOqhL/eEISu/GC1cLDG8DUdjSb/N6NUYW8/exoXbGdhxMR5DAlyljkRERDqm9ZDbO++8g4ULFyIjIwNRUVGIjIxUX27duqXLjPQEnIu5i42nYwEA7w9uAQVHr0jP2VsY4eWePgCAz/65jvwiHoRKRFTXaF3IFhYWYtSoUZDL+fWzvisoLsGbWy9CCGB4G1e087CROhKRTkzu4gkXK2PEZ+TjpyORj9+AiIj0itZV6IQJE7Bp0yZdZiGJfL3vJm4mZcPWTIn3BvhKHYdIZ4wNFXizfzMAwHf7w5GUlS9xIiIi0iWt58iWlJTg888/x+7du+Hn51fuYK8lS5ZUOxzVvGMRKVhxMAIA8NHQlmhgppQ4EZFuDfJzweqjUbgQm46le29g0XA/qSMREZGOaF3IXrp0Ca1btwYAXL58WeO2+w/8otorMTMfszaGQgjguXZu6N+Kp6Klukcul+G9Ac0xYuVxbDodi/EdPdDc2VLqWEREpANaF7L79+/XZQ56wvKLSvDir2eRnFWApo4WmD+IUwqo7mrrYYMBrZzx16UEfPzXNfw6uT3/4CYiqgN4pFY9VFSiwsvrzuFCbDqsTQ2xanxbmCq1/puGSC+81b8ZlAo5joSnYNflO1LHISIiHahWIXv48GE8//zz6NixI+Li4gAAv/76K44cOaKTcKR7BcUlmLH+HP69ngRjQzlWjW+LRramUsciqnFuNqZ4qbsXAOCDnVeRW1gscSIiIqourQvZ3377DcHBwTAxMcH58+dRUFAAAMjIyMAnn3yis4CkO+m5hXhhzWnsvpIIpUKOFWMD2WqL6pXpPXzQsIEJEjLy8c2/4VLHISKiatK6kP3oo4+wcuVKrFq1SqNjQefOnXHu3DmdhCPdORdzF4OWH8HR8FSYKhVYPbEdejZzkDoW0RNlolTg/UEtAAA/Hr6F8KRsiRMREVF1aF3IhoWFoVu3buWWW1lZIT09vTqZ8O2338LDwwPGxsYICgrCqVOnHrn+li1b0KxZMxgbG6NVq1b4+++/NW4XQmD+/PlwdnaGiYkJ+vTpg5s3b1Yro75Izy3Ewh1XMGLFMcSm5cHNxgRbX+qELo3tpI5GJIk+vo7o3cwBRSUCC/68DCGE1JGIiEhLWheyTk5OCA8v/9XckSNH4OXlpXWgTZs2Yc6cOViwYAHOnTsHf39/BAcHIykpqcL1jx07htGjR2Py5Mk4f/48hg4diqFDh2q0BPv888+xbNkyrFy5EidPnoSZmRmCg4ORn193m6NHJGfjk7+voetn+/Hz0SioBDCstSt2zugKXxe2HqL67f3BLWBkIMfR8FTsvJggdRwiItKSTGg5HLFo0SL83//9H1avXo2nnnoKf//9N6Kjo/Haa6/hvffew6uvvqpVoKCgILRr1w7Lly8HAKhUKri5ueHVV1/FW2+9VW79UaNGIScnBzt37lQv69ChAwICArBy5UoIIeDi4oL//e9/mDt3LoDSebyOjo5Ys2YNnnvuucdmyszMhJWVFTIyMmBpWfuKwOyCYty+m4sbidk4FZmKk7fScPO+r0ybOVng7aebo1sTewlTEtUuX++7iaX7bsDR0gj75nSHhbHh4zciIqIaV5W6S+ueS2+99RZUKhV69+6N3NxcdOvWDUZGRpg7d67WRWxhYSHOnj2LefPmqZfJ5XL06dMHx48fr3Cb48ePY86cORrLgoODsX37dgBAZGQk7ty5gz59+qhvt7KyQlBQEI4fP15hIVtQUKA+eA0ofUIBoKioCEVFRVo9tsradfkOtofGQSUEhAAEBFSq0ukRpdcBlQDyCkuQXViErLxiZOaXz2RqKENnb1s8174RunjbQS6X1Xh2In0yubMb/roQi5i0XHy56yre5emZiYhqharUK1oXsjKZDO+88w5ef/11hIeHIzs7G76+vjA3N9d2l0hJSUFJSQkcHR01ljs6OuL69esVbnPnzp0K179z54769rJlD1vnQYsWLcLChQvLLQ8JCYGpac22qpIBGKaTgVMBIAV5t1Kw95Yu9kdU98xsWvZTLHbvjpUyChER3ZObm1vpdavVBT8kJAQhISFISkqCSqXSuG316tXV2bWk5s2bpzHKm5mZCTc3N/Tu3bvGpxbcTMrCtfgsyOWADDLIZIBcJoNMVvYzAMhgYiiHmZEBzI0N4GRpDHMjntCASBsLd1zF1rOxaGRjim0vd4KRgULqSERUxwghkF+sQlZeMbILi5FXWIycgmLkFqqQW1iMnIIS5BUVI7ugGHmFKuQUFiO3oBi5hSWlP99bL6+wBIXFJcgvUiG/uAQlqid7sOpL3b3xSk+fGr+fsm/CK0Pr6mfhwoX44IMP0LZtWzg7O+vkdI92dnZQKBRITEzUWJ6YmAgnJ6cKt3Fycnrk+mX/JiYmwtnZWWOdgICACvdpZGQEIyOjcssNDQ01Wo3VBF9XG/i6srcr0ZPyxtO+2Hs9GTeS8/DNgUjM699c6khEVEvlFhYjNbsQaTmFyMgrQmZ+Uem/ecUPXC9CZn5x6b/3lheV1ETR+V/tZWQgh7GhAsaG9/41KP3ZyFBx73rpz0qFHEoDGZQKOQwVcigN/vtX+cB1Q4UMRvddd7cxq/E6CECV7kPrQnblypVYs2YNxo0bp+0uylEqlQgMDERISAiGDh0KoPRgr5CQEMyYMaPCbTp27IiQkBDMnj1bvWzv3r3o2LEjAMDT0xNOTk4ICQlRF66ZmZk4efIkpk+frrPsRKSfLI0N8fHQVpjyyxmsOnQLA1o5w6+htdSxiOgJyC8qQXJWAZKyCpCaXYC7uYVIzSlE2r1iNTWn9N/SnwuQX6R6/E4fQSGXwdzIAGZKBcyMDGBqZABzIwVMlQYwNzKA6b3lZkoDmBndW0epuHebAUyU9wpVA4VG0apUyCGXV39AUR9pXcgWFhaiU6dOuswCAJgzZw4mTJiAtm3bon379vjqq6+Qk5ODSZMmAQDGjx8PV1dXLFq0CAAwa9YsdO/eHYsXL8aAAQOwceNGnDlzBj/88AOA0rm8s2fPxkcffYTGjRvD09MT7733HlxcXNTFMhHVb318HTHI3wU7LsTjja0X8eeMLlAaVOsM3kQkocJiFRIz85GUlY/EzAIkZeYjMasAiZn5SL73b2JmATLyqn4QtNJADlszJaxMDGFlYghLE0NYGpf9bFD6r/F/t92/3MRQoZNvsOk/WheyU6ZMwfr16/Hee+/pMg9GjRqF5ORkzJ8/H3fu3EFAQAB27dqlPlgrJiYGcvl//8F06tQJ69evx7vvvou3334bjRs3xvbt29GyZUv1Om+88QZycnLw4osvIj09HV26dMGuXbtgbGys0+xEpL/eH+SLIzeTcf1OFpaF3MTc4KaP34iIJJFdUIz49DzE3c3D7Xv/xqXnIe5uLuLS85CUVYDKNhc1MpDD3sIItuZGsDVTwuaBi636XyPYmCthpmQxWptUqY/s/QdAqVQqrF27Fn5+fvDz8ys3n2HJkiW6Symx2t5Hloh04+9LCXh53TnIZcCmaR3RzoPz1YmkIIRAYmYBolJzEJWSg6jUXESl5CAmrbRQrcxIqtJADkdLIzhaGMPB0ggOFsZwtDSGg4URHC2N4XhvmaWJAQvTWqbG+sieP39e43rZnNP7z6IFgG8IItJLT7dyxjNtGuK3c7fx2qZQ/DOrK0+UQFSDMnKLcDMpC+FJ2YhMzUF0Sm5p8Zqa89j5qJbGBnBtYApXaxM0bGACV2sTuN73r62ZkvVIPaD1mb3qE47IEtUfWflF6P/1Ydy+m4dn2jTE4pH+Ukci0nup2QW4mZSNm0nZCE/MUv+cnFXw0G0UchkaNjCBh60ZPGxN4W5rBndbU3Wxyj8y664ncmYvIqK6yMLYEEtHBWDU98fx27nb6NXMAQP8nB+/IRGhsFiF8KRsXInPwNWETFyNz8TNpGyk5RQ+dBtXaxN4O5jDy660UPWwM4OHrRkaNjCBoYIHXdKjsZAlInpAOw8bvNzDB8v3h+Pt3y+hjbs1nK1MpI5FVKtk5hfhWnwmrsRn3le0Zj20X6qbjQkaO1igsYM5fBzM0djRAj4O5jyhD1UL3z1ERBWY1acxDt1MxsXbGXh1/XlseLEDR4eo3iooLsHV+EyExqbjQmw6QmPTEZVa8WlELYwN4OtsCV8XS/g6W6KZkyW8HcxgqmTJQbrHdxURUQUMFXJ8M7o1Bn5zBGei7+Kzf67j3YG+UsciqnEqlUBUag5C7xWsF2LTcTUhs8KRVldrEzR3tkQLl/8K14YNTHiQFT0xLGSJiB7C3dYMXz7rj2m/nsWPRyLR1qMB+rXkfFmqWwqLVbgUl4HTUWk4FZmGs9F3K2xvZWOmRICbNfwbWiOgkTX8XK3QwEwpQWKi/7CQJSJ6hOAWTnixmxd+OHQLr2+5iKZOlvC0M5M6FpHWcguLcT4mHScj03A6Mg3nY++Wa3WlNJCjpYslAtwaIKCRNQIaWsPNhiOtVPuwkCUieozXg5vifMxdnI66i+n/dxbbX+kMY0OF1LGIKiWvsASno9JwNDwFJyPTcDkuA8UqzWkCDUwN0c7DBu09bdDOwwbNnS15mmbSCyxkiYgew1Ahx/IxbTBg2WFcv5OFt367iKWjAjg6RbVSiUrgUlwGjoan4MjNFJyNvovCEs0RVxcrY7TzLC1c23vYwNveHHI538+kf1jIEhFVgqOlMZaNbo1xP53C9tB4+DiYY0avxlLHIoIQApEpOaWFa3gKjkekIjO/WGMdFytjdPaxQ0dvW7T3tEHDBqYSpSXSLRayRESV1MnbDh8MaYF3fr+ML/fcgKedOU+WQJLIKyzBiVup+Pd6EvaHJeH23TyN2y2MDdDJ2xZdfOzQ2ccOnnZm/AaB6iQWskREVTA2yB3hSdn4+WgU/rclFG42JvBraC11LKoHYtNysT8sCfuvJ+FYRCoKiv+bLqBUyBHo3gBdGpcWri1dLGHAvsdUD7CQJSKqoncH+CIyJQcHwpIxZe0Z/DGjM8/8RTpXVKLC6ag0HAhLxr/XkxCelK1xu4uVMXo2c0DPpg7o5GPLEw5QvSQTQlR8LjlSy8zMhJWVFTIyMmBpaSl1HCKqBbLyi/DMimO4kZiNpo4W2DStA6xN2VOTqievsASHbiZj95U7CLmWpNHPVSGXIdC9AXrdK16bOJpzugDVSVWpu1jIVgILWSKqSGxaLp5ZcQxJWQVo3cga66YEcVSMqiw9txAh15Kw+8odHLqZrNHT1cZMiR5N7dGrmQO6NraHlYmhhEmJngwWsjrGQpaIHibsThZGfn8cGXlF6NrYDj9OaAsjA/aYpUdLyMjDniuJ2HP1Dk7cSkPJfX1dXa1NENzCCcEtHNHWwwYKtsWieoaFrI6xkCWiRzkXcxdjV51EXlEJBrRyxrLRrVl8UDnx6XnYeTEef11MwIXbGRq3NXOyQF9fR/Rt4YQWLpacMkD1WlXqLn4HRkRUTW0aNcAP4wPxwprT+OtSAsyMFFg03I/FLCEluwD/XErAnxficTrqrnq5TFb6vglu4Yi+vk7w4GmPibTCQpaISAe6NrbHsuda45X157D5zG0UlQh8McKPLZDqoYy8Iuy+cgc7LsTjWESqxrSB9h42GOTvjOCWTnCwMJYwJVHdwEKWiEhH+rdyxtfPtcbsTaH4/XwcCotV+Oq5ABiymK3z8gpLsO9aInZciMeBsGSNU8L6NbTCID8XDPR3Zps2Ih1jIUtEpEOD/F2gNJBjxvpz+OtSAgqKVfh2bGseAFYHFRarcOhGMv68EI991xKRW1iivq2xgzkG+7tgoL8LPDltgKjG8GCvSuDBXkRUVfvDkvDSr2dRUKxC18Z2+G5sG1gYs3WSvitRCZy4lYo/Q+Pxz+UEZOYXq29zszHBID8XDA5wQVNHCx6wRaQldi3QMRayRKSNY+EpmLz2DPKKStDMyQKrJ7aDizW/WtY3Qgici7mLHRcSsPNiAlKyC9S3OVgYYaCfCwb5OyPAzZrFK5EOsJDVMRayRKSti7fTMXntGSRnFcDBwgg/TWiHVg2tpI5FjyGEwNWETOy4kIAdF+IRl56nvs3a1BD9WzpjsL8L2nuyzyuRrrGQ1TEWskRUHbfv5mLymjMIS8yCiaECy0a3xlO+jlLHogrcSs7GnxfiseNCPCKSc9TLzZQK9G3hhMH+LujsYwelAQ/gI6opLGR1jIUsEVVXZn4RXll3DodvpkAmA2b09MHsPk04mlcLxKXnYeeFePx5IR5X4jPVy5UGcvRu5oBB/i7o1cwBxoY8YI/oSWAhq2MsZIlIF4pKVPhgx1X8eiIaANDBywbLnmsNB0v2E33SkrMK8M/lBPwZGo8z0f+dqEAhl6FrYzsM8nNB3xaOPECPSAIsZHWMhSwR6dIfoXGYt+0ScgtLYGduhGXPBaCTj53Useq89NxC7L5yBzsvJuBoeArKzlMgk5WeqGBwgAv6t3SGjZlS2qBE9RwLWR1jIUtEuhaelI1X1p1DWGIWZDLghc6emNu3KUyU/PpalzLzi7D3SiJ2XozHkfAUFJX891+ev5s1Bvk5Y6CfC5ysOCpOVFuwkNUxFrJEVBPyCkuwcMcVbDwdCwBwtzXFp8P90NHbVuJk+i23sBj7riVh54V4HLiRjMLi/86y1czJAoP8XTDQzxnutjxRAVFtxEJWx1jIElFN2h+WhLe3XUJCRj4AYGxQI7wR3AxWppyfWVn5RSXYfz0JOy8mIOR6IvKL/iteve3N7hWvLvBxMJcwJRFVBgtZHWMhS0Q1LSu/CIv+uY71J2MAlPYqnd27McZ2cIehgq2eKpKZX4T915Ow50oiDoQlIee+U8S625pi4L1pA82ceJYtIn3CQlbHWMgS0ZNyLDwFC/68gptJ2QAAL3szvN2/OXo3d2AxBiApMx97ryVi95VEHI/QnPPqam2iLl5bulry+SLSUyxkdYyFLBE9ScUlKmw8HYule28gNacQQOmBSS/38MZTzR0hr0e9Z4UQCE/Kxr/Xk7D7yh2cj03H/f9redubIbiFE4JbOKGVq1W9em6I6ioWsjrGQpaIpJCZX4QVByKw+kgkCu4dsOTjYI6XuntjsL9LnT27VFZ+EY6Gp+LgjWQcupGscXpYAAhws0bfFo7o6+vEOa9EdRALWR1jIUtEUkrOKsCaY5H45Xg0svKLAQC2ZkoMb+OKUe3c4ONgIXHC6ikuUeFyfCaORaTgYFgyzkbfRbHqv/+alAZyBHnaoK+vI57ydWKrLKI6joWsjrGQJaLaIDO/COtPxmD1kUgkZRWol7dpZI2hrV3Rp7kjXKxNJExYOUUlKlyKy8CJW6k4eSsNZ6LSNA7UAgBPOzN0b2KP7k3t0cHTlv11ieoRFrI6xkKWiGqT4hIVDoQlY9OZWPx7PQkl941etnCxRJ/mjuje1B4tXawkn36gUglEpebg4u0MhMam4+LtdFyJz1RPlShjaWyA9p626NbEDt2b2LPHK1E9xkJWx1jIElFtlZSVjz/Ox2PP1Ts4E31X40AoIwM5/N2s0da9AfzdrNHYwRyNbExhUAPtvIQQSMjIx82kbIQnZSM8KQs3E7NxIzELmfemQ9zP2tQQQZ42CPK0RZCXDZo5WULBA7WICCxkdY6FLBHpg9TsAvx7PQn7riXiVGQa7uYWlVtHqZDDw84UnnZmcLAwhr2FERwsjGBrbgQTQwWMDeUwNlRAaSBHiUqguESgSKVCUbEKmfnFyMgrQnpuIdJzi5CYmY/4jDzEp+cjPj2v3Cir+j4N5GjpYgm/htYIcLOGX0MreNiascMAEVWIhayOsZAlIn0jhEBEcg7ORqfhTNRdXE3IRERytsYZr3TNQC6Dh50ZGjuYw+e+SxNHC57UgYgqrSp1l8ETykRERE+QTCZTF5Kj2jUCUDpfNS49D+HJ2YhNy0VyVgGSMguQnF2A1JxC5BeWIL+4BPlFJSgsVkEhl8FQIYeBQgZDuRwWJoawNjGEtakhrEwM4WBhBBdrEzhbmcDV2gTO1sYsWInoiWIhS0RUT8jlMrjZmMLNxlTqKEREOsE/nYmIiIhIL7GQJSIiIiK9xEKWiIiIiPQSC1kiIiIi0kssZImIiIhIL7GQJSIiIiK9xEKWiIiIiPQSC1kiIiIi0kssZImIiIhIL7GQJSIiIiK9xEKWiIiIiPSSgdQB9IEQAgCQmZkpcRIiIiKiuq2s3iqrvx6FhWwlZGVlAQDc3NwkTkJERERUP2RlZcHKyuqR68hEZcrdek6lUiE+Ph4WFhaQyWQ1fn+ZmZlwc3NDbGwsLC0ta/z+SDf4uuknvm76ia+bfuLrpp+e9OsmhEBWVhZcXFwglz96FixHZCtBLpejYcOGT/x+LS0t+UHXQ3zd9BNfN/3E100/8XXTT0/ydXvcSGwZHuxFRERERHqJhSwRERER6SUWsrWQkZERFixYACMjI6mjUBXwddNPfN30E183/cTXTT/V5teNB3sRERERkV7iiCwRERER6SUWskRERESkl1jIEhEREZFeYiFby3z77bfw8PCAsbExgoKCcOrUKakj0SO8//77kMlkGpdmzZpJHYsqcOjQIQwaNAguLi6QyWTYvn27xu1CCMyfPx/Ozs4wMTFBnz59cPPmTWnCktrjXreJEyeW+wz269dPmrCktmjRIrRr1w4WFhZwcHDA0KFDERYWprFOfn4+XnnlFdja2sLc3BzPPPMMEhMTJUpMQOVetx49epT7zL300ksSJWYhW6ts2rQJc+bMwYIFC3Du3Dn4+/sjODgYSUlJUkejR2jRogUSEhLUlyNHjkgdiSqQk5MDf39/fPvttxXe/vnnn2PZsmVYuXIlTp48CTMzMwQHByM/P/8JJ6X7Pe51A4B+/fppfAY3bNjwBBNSRQ4ePIhXXnkFJ06cwN69e1FUVIS+ffsiJydHvc5rr72GHTt2YMuWLTh48CDi4+MxfPhwCVNTZV43AJg6darGZ+7zzz+XKDEAQbVG+/btxSuvvKK+XlJSIlxcXMSiRYskTEWPsmDBAuHv7y91DKoiAOL3339XX1epVMLJyUl88cUX6mXp6enCyMhIbNiwQYKEVJEHXzchhJgwYYIYMmSIJHmo8pKSkgQAcfDgQSFE6efL0NBQbNmyRb3OtWvXBABx/PhxqWLSAx583YQQonv37mLWrFnShXoAR2RricLCQpw9exZ9+vRRL5PL5ejTpw+OHz8uYTJ6nJs3b8LFxQVeXl4YO3YsYmJipI5EVRQZGYk7d+5ofP6srKwQFBTEz58eOHDgABwcHNC0aVNMnz4dqampUkeiB2RkZAAAbGxsAABnz55FUVGRxmeuWbNmaNSoET9ztciDr1uZdevWwc7ODi1btsS8efOQm5srRTwAgIFk90waUlJSUFJSAkdHR43ljo6OuH79ukSp6HGCgoKwZs0aNG3aFAkJCVi4cCG6du2Ky5cvw8LCQup4VEl37twBgAo/f2W3Ue3Ur18/DB8+HJ6enoiIiMDbb7+N/v374/jx41AoFFLHIwAqlQqzZ89G586d0bJlSwClnzmlUglra2uNdfmZqz0qet0AYMyYMXB3d4eLiwsuXryIN998E2FhYdi2bZskOVnIElVD//791T/7+fkhKCgI7u7u2Lx5MyZPnixhMqL64bnnnlP/3KpVK/j5+cHb2xsHDhxA7969JUxGZV555RVcvnyZxw/omYe9bi+++KL651atWsHZ2Rm9e/dGREQEvL29n3RMHuxVW9jZ2UGhUJQ7YjMxMRFOTk4SpaKqsra2RpMmTRAeHi51FKqCss8YP3/6z8vLC3Z2dvwM1hIzZszAzp07sX//fjRs2FC93MnJCYWFhUhPT9dYn5+52uFhr1tFgoKCAECyzxwL2VpCqVQiMDAQISEh6mUqlQohISHo2LGjhMmoKrKzsxEREQFnZ2epo1AVeHp6wsnJSePzl5mZiZMnT/Lzp2du376N1NRUfgYlJoTAjBkz8Pvvv+Pff/+Fp6enxu2BgYEwNDTU+MyFhYUhJiaGnzkJPe51q0hoaCgASPaZ49SCWmTOnDmYMGEC2rZti/bt2+Orr75CTk4OJk2aJHU0eoi5c+di0KBBcHd3R3x8PBYsWACFQoHRo0dLHY0ekJ2drTFiEBkZidDQUNjY2KBRo0aYPXs2PvroIzRu3Bienp5477334OLigqFDh0oXmh75utnY2GDhwoV45pln4OTkhIiICLzxxhvw8fFBcHCwhKnplVdewfr16/HHH3/AwsJCPe/VysoKJiYmsLKywuTJkzFnzhzY2NjA0tISr776Kjp27IgOHTpInL7+etzrFhERgfXr1+Ppp5+Gra0tLl68iNdeew3dunWDn5+fNKGlbptAmr755hvRqFEjoVQqRfv27cWJEyekjkSPMGrUKOHs7CyUSqVwdXUVo0aNEuHh4VLHogrs379fACh3mTBhghCitAXXe++9JxwdHYWRkZHo3bu3CAsLkzY0PfJ1y83NFX379hX29vbC0NBQuLu7i6lTp4o7d+5IHbveq+g1AyB+/vln9Tp5eXni5ZdfFg0aNBCmpqZi2LBhIiEhQbrQ9NjXLSYmRnTr1k3Y2NgIIyMj4ePjI15//XWRkZEhWWbZveBERERERHqFc2SJiIiISC+xkCUiIiIivcRCloiIiIj0EgtZIiIiItJLLGSJiIiISC+xkCUiIiIivcRCloiIiIj0EgtZIiIiItJLLGSJ6LEOHDgAmUyG9PR0qaNUaPv27fDx8YFCocDs2bOljlOh999/HwEBAVpt6+Hhga+++qpK2/zwww9wc3ODXC6v8ra1RY8ePWrt66mtwsJC+Pj44NixY5Va18PDA2fOnHkCyYj0EwtZItJQUfHQqVMnJCQkwMrKSppQjzFt2jSMGDECsbGx+PDDDzFx4kQMHTpU6lg6c/r0abz44ouVXj8zMxMzZszAm2++ibi4uCptS1VXlT9SVq5cCU9PT3Tq1Omx6yqVSsydOxdvvvlmNRMS1V0sZInosZRKJZycnCCTyaSOUk52djaSkpIQHBwMFxcXWFhY6GzfhYWFOttXde7f3t4epqamld4uJiYGRUVFGDBgAJydnau07ZNQVFQkdQRJCCGwfPlyTJ48udLbjB07FkeOHMGVK1dqMBmR/mIhS0RqEydOxMGDB/H1119DJpNBJpMhKiqq3NSCNWvWwNraGjt37kTTpk1hamqKESNGIDc3F2vXroWHhwcaNGiAmTNnoqSkRL3/goICzJ07F66urjAzM0NQUBAOHDjwyExLlixBq1atYGZmBjc3N7z88svIzs4GUDrloaxw7dWrF2QyGXr06IG1a9fijz/+UD+GsvuIjY3FyJEjYW1tDRsbGwwZMgRRUVEaj3/o0KH4+OOP4eLigqZNm1aYKSIiAkOGDIGjoyPMzc3Rrl077Nu3r1LP8ffffw83NzeYmppi5MiRyMjIeOz9Pzi1ID09HVOmTIG9vT0sLS3Rq1cvXLhwQf3atGrVCgDg5eWlfg0fFBUVBZlMhm3btqFnz54wNTWFv78/jh8/rrHekSNH0LVrV5iYmMDNzQ0zZ85ETk6O+naZTIbt27drbGNtbY01a9Zo3M+mTZvQvXt3GBsbY926dUhNTcXo0aPh6uoKU1NTtGrVChs2bKjUc3i/P/74A23atIGxsTG8vLywcOFCFBcXAwDGjBmDUaNGaaxfVFQEOzs7/PLLLwAAlUqFRYsWwdPTEyYmJvD398fWrVvV65e990NCQtC2bVuYmpqiU6dOCAsLUz/fCxcuxIULF9Tvt7LH/qCzZ88iIiICAwYMUC8rLCzEjBkz4OzsDGNjY7i7u2PRokXq2xs0aIDOnTtj48aNVX5uiOoFQUR0T3p6uujYsaOYOnWqSEhIEAkJCaK4uFjs379fABB3794VQgjx888/C0NDQ/HUU0+Jc+fOiYMHDwpbW1vRt29fMXLkSHHlyhWxY8cOoVQqxcaNG9X7nzJliujUqZM4dOiQCA8PF1988YUwMjISN27ceGimpUuXin///VdERkaKkJAQ0bRpUzF9+nQhhBAFBQUiLCxMABC//fabSEhIEBkZGWLkyJGiX79+6sdQUFAgCgsLRfPmzcULL7wgLl68KK5evSrGjBkjmjZtKgoKCoQQQkyYMEGYm5uLcePGicuXL4vLly9XmCk0NFSsXLlSXLp0Sdy4cUO8++67wtjYWERHRz/0cSxYsECYmZmJXr16ifPnz4uDBw8KHx8fMWbMGPU6D7t/d3d3sXTpUvV6ffr0EYMGDRKnT58WN27cEP/73/+Era2tSE1NFbm5uWLfvn0CgDh16pT6NXxQZGSkACCaNWsmdu7cKcLCwsSIESOEu7u7KCoqEkIIER4eLszMzMTSpUvFjRs3xNGjR0Xr1q3FxIkT1fsBIH7//XeNfVtZWYmff/5Z4348PDzEb7/9Jm7duiXi4+PF7du3xRdffCHOnz8vIiIixLJly4RCoRAnT55U76d79+5i1qxZD31ODx06JCwtLcWaNWtERESE2LNnj/Dw8BDvv/++EEKInTt3ChMTE5GVlaXeZseOHcLExERkZmYKIYT46KOPRLNmzcSuXbtERESE+Pnnn4WRkZE4cOCAEEKo3/tBQUHiwIED4sqVK6Jr166iU6dOQgghcnNzxf/+9z/RokUL9fstNze3wrxLliwRzZo101j2xRdfCDc3N3Ho0CERFRUlDh8+LNavX6+xzptvvim6d+/+0OeBqD5jIUtEGioqHioqZAGI8PBw9TrTpk0TpqamGkVDcHCwmDZtmhBCiOjoaKFQKERcXJzGvnv37i3mzZtX6XxbtmwRtra26ut3794VAMT+/fvVyyZMmCCGDBmisd2vv/4qmjZtKlQqlXpZQUGBMDExEbt371Zv5+joqC5sq6JFixbim2++eejtCxYsEAqFQty+fVu97J9//hFyuVwkJCQ88v7vL2QPHz4sLC0tRX5+vsY63t7e4vvvvxdCCHH+/HkBQERGRj40T1mB+eOPP6qXXblyRQAQ165dE0IIMXnyZPHiiy9qbHf48GEhl8tFXl6eEKLyhexXX3310CxlBgwYIP73v/+prz+ukO3du7f45JNPNJb9+uuvwtnZWQghRFFRkbCzsxO//PKL+vbRo0eLUaNGCSGEyM/PF6ampuLYsWMa+5g8ebIYPXq0EOK/9/6+ffvUt//1118CgPo5WLBggfD393/s45s1a5bo1auXxrJXX31V9OrVS+N9+aCvv/5aeHh4PHb/RPWRwRMfAiaiOsHU1BTe3t7q646OjvDw8IC5ubnGsqSkJADApUuXUFJSgiZNmmjsp6CgALa2tg+9n3379mHRokW4fv06MjMzUVxcjPz8fOTm5lZp7ueFCxcQHh5ebg5tfn4+IiIi1NdbtWoFpVL5yH1lZ2fj/fffx19//YWEhAQUFxcjLy8PMTExj9yuUaNGcHV1VV/v2LEjVCoVwsLC4OTkVKn7v3DhArKzs8s9Z3l5eRqPo7L8/PzUPzs7OwMAkpKS0KxZM1y4cAEXL17EunXr1OsIIaBSqRAZGYnmzZtX+n7atm2rcb2kpASffPIJNm/ejLi4OBQWFqKgoKDKr+nRo0fx8ccfa+z3/vfHyJEjsW7dOowbNw45OTn4448/1F/Th4eHIzc3F0899ZTGfgsLC9G6dWuNZQ97nho1alTpvHl5eTA2NtZYNnHiRDz11FNo2rQp+vXrh4EDB6Jv374a65iYmCA3N7fS90NUn7CQJSKtGBoaalyXyWQVLlOpVABKiz+FQoGzZ89CoVBorHd/8Xu/qKgoDBw4ENOnT8fHH38MGxsbHDlyBJMnT0ZhYWGVip7s7GwEBgZqFGVl7O3t1T+bmZk9dl9z587F3r178eWXX8LHxwcmJiYYMWKETg4Oe9z9Z2dnw9nZucK5xdbW1lW+v/tfs7KD+e5/zaZNm4aZM2eW266sgJPJZBBCaNxW0cFcDz6uL774Al9//TW++uor9Rzo2bNnV+k5zM7OxsKFCzF8+PByt5UVjGPHjkX37t2RlJSEvXv3wsTEBP369VNvDwB//fWXxh8YAGBkZKRx/VHPU2XZ2dnh0qVLGsvatGmDyMhI/PPPP9i3bx9GjhyJPn36aMzTTUtL03iPEtF/WMgSkQalUqlxgJautG7dGiUlJUhKSkLXrl0rtc3Zs2ehUqmwePFiyOWlx6Zu3rz5sdtV9BjatGmDTZs2wcHBAZaWllV/APc5evQoJk6ciGHDhgEoLYgqOqDqQTExMYiPj4eLiwsA4MSJE5DL5Q89qKwibdq0wZ07d2BgYAAPDw9t4lfpvq5evQofH5+HrmNvb4+EhAT19Zs3b1Zq9PDo0aMYMmQInn/+eQClReGNGzfg6+tbpXxhYWGPzNepUye4ublh06ZN+Oeff/Dss8+qi1JfX18YGRkhJiYG3bt3r/T9Pqiyn5nWrVtjxYoVEEJodACxtLTEqFGjMGrUKIwYMQL9+vVDWloabGxsAACXL18uN0JMRKXYtYCINHh4eODkyZOIiopCSkpKlUedHqZJkyYYO3Ysxo8fj23btiEyMhKnTp3CokWL8Ndff1W4jY+PD4qKivDNN9/g1q1b+PXXX7Fy5cpKPYaLFy8iLCwMKSkpKCoqwtixY2FnZ4chQ4bg8OHDiIyMxIEDBzBz5kzcvn27So+lcePG2LZtG0JDQ3HhwgWMGTOmUs+TsbExJkyYgAsXLuDw4cOYOXMmRo4cqZ5WUBl9+vRBx44dMXToUOzZswdRUVE4duwY3nnnHZ03zn/zzTdx7NgxzJgxA6Ghobh58yb++OMPzJgxQ71Or169sHz5cpw/fx5nzpzBSy+9VG5kviKNGzfG3r17cezYMVy7dg3Tpk1DYmJilfLNnz8fv/zyCxYuXIgrV67g2rVr2LhxI959912N9caMGYOVK1di7969GDt2rHq5hYUF5s6di9deew1r165FREQEzp07h2+++QZr166tdA4PDw9ERkYiNDQUKSkpKCgoqHC9nj17Ijs7W6OV1pIlS7BhwwZcv34dN27cwJYtW+Dk5KQxun748OFy0w2IqBQLWSLSMHfuXCgUCvj6+sLe3v6x8z6r4ueff8b48ePxv//9D02bNsXQoUNx+vTph84z9Pf3x5IlS/DZZ5+hZcuWWLdunUZrooeZOnUqmjZtirZt28Le3h5Hjx6FqakpDh06hEaNGmH48OFo3rw5Jk+ejPz8/CqP0C5ZsgQNGjRAp06dMGjQIAQHB6NNmzaP3c7HxwfDhw/H008/jb59+8LPzw/fffddle5bJpPh77//Rrdu3TBp0iQ0adIEzz33HKKjo+Ho6FilfT2On58fDh48iBs3bqBr165o3bo15s+frx5RBoDFixfDzc0NXbt2xZgxYzB37txKTfl499130aZNGwQHB6NHjx5wcnKq8kksgoODsXPnTuzZswft2rVDhw4dsHTpUri7u2usN3bsWFy9ehWurq7o3Lmzxm0ffvgh3nvvPSxatAjNmzdHv3798Ndff8HT07PSOZ555hn069cPPXv2hL29/UPbiNna2mLYsGEa01ssLCzw+eefo23btmjXrh2ioqLw999/q7+BOH78ODIyMjBixIhK5yGqT2TiwclNREREVCMuXryIp556ChEREQ+dG36/UaNGwd/fH2+//fYTSEekfzgiS0RE9IT4+fnhs88+Q2Rk5GPXLSwsRKtWrfDaa689gWRE+okjskRERESklzgiS0RERER6iYUsEREREeklFrJEREREpJdYyBIRERGRXmIhS0RERER6iYUsEREREeklFrJEREREpJdYyBIRERGRXmIhS0RERER6iYUsEREREeml/weBqToFdqOo2wAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import brainmass\n", "import brainunit as u\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "# Visualise a canonical hemodynamic response function: the impulse response that\n", "# convolution-BOLD applies to neural activity. Note the ~5-6 s peak and post-undershoot.\n", "kernel = brainmass.MixtureOfGammasHRFKernel()\n", "t = np.linspace(0.0, 25.0, 400) * u.second\n", "h = kernel(t)\n", "\n", "fig, ax = plt.subplots(figsize=(7, 4))\n", "ax.plot(np.asarray(t.to_decimal(u.second)), np.asarray(h))\n", "ax.axhline(0.0, color='0.7', lw=0.8)\n", "ax.set_xlabel('time after a brief neural event (s)')\n", "ax.set_ylabel('hemodynamic response $h(t)$')\n", "ax.set_title('Canonical HRF: activity is blurred & delayed by ~5 s')\n", "fig.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "0116fb86", "metadata": {}, "source": [ "The HRF makes the BOLD assumptions visceral: a brief burst of neural activity\n", "produces a response that **peaks several seconds later** and undershoots before\n", "returning to baseline. fMRI therefore reports a slow, lagged, smoothed image of the\n", "underlying fast dynamics — which is exactly why whole-brain BOLD models compare\n", "*functional connectivity* (correlation structure over minutes) rather than\n", "millisecond waveforms." ] }, { "cell_type": "markdown", "id": "5b098313", "metadata": {}, "source": [ "## EEG / MEG: the lead-field forward model\n", "\n", "EEG and MEG see **electromagnetism**. Synchronous post-synaptic currents in\n", "populations of pyramidal neurons act, at a distance, like **equivalent current\n", "dipoles** (ECDs). Each region's dipole sets up an electric potential (measured by\n", "EEG electrodes) and a magnetic field (measured by MEG sensors) throughout the head.\n", "Because the quasi-static Maxwell equations are **linear**, the mapping from source\n", "currents to sensors is a single matrix multiply.\n", "\n", "### The lead-field equation\n", "\n", "``brainmass``'s ``LeadFieldModel`` implements\n", "\n", "$$\n", "\\mathbf{y}(t) \\;=\\; \\mathbf{s}(t)\\,\\mathbf{L} \\;+\\; \\boldsymbol{\\varepsilon}(t),\n", "$$\n", "\n", "where\n", "\n", "- $\\mathbf{s}(t) \\in \\mathbb{R}^{R}$ is the **equivalent current dipole moment** per\n", " region (typically in $\\mathrm{nA\\cdot m}$), obtained from the NMM observable via a\n", " scale factor.\n", "- $\\mathbf{L} \\in \\mathbb{R}^{R \\times M}$ is the **lead-field** (gain) matrix — the\n", " signal each unit dipole produces at each of $M$ sensors. Its units are\n", " $\\mathrm{V/(nA\\cdot m)}$ for EEG and $\\mathrm{T/(nA\\cdot m)}$ for MEG.\n", "- $\\mathbf{y}(t) \\in \\mathbb{R}^{M}$ are the sensor measurements (volts for EEG,\n", " tesla for MEG).\n", "- $\\boldsymbol{\\varepsilon}(t) \\sim \\mathcal{N}(\\mathbf{0}, \\boldsymbol{\\Sigma})$ is\n", " additive sensor noise.\n", "\n", "The lead field $\\mathbf{L}$ encodes all the head geometry and physics — the\n", "conductivity of brain, skull, and scalp; sensor positions; source orientations — and\n", "is normally precomputed by a **head model** (boundary- or finite-element method) in\n", "a tool such as MNE, FieldTrip, or Brainstorm. ``brainmass`` consumes that\n", "$\\mathbf{L}$ and provides the differentiable, **unit-safe** projection, aggregating\n", "vertex-level lead fields to region level when needed." ] }, { "cell_type": "markdown", "id": "839f844e", "metadata": {}, "source": [ "### What each modality assumes — and how they differ\n", "\n", "The *same* dipolar sources feed both EEG and MEG, but the physics treats them\n", "differently:\n", "\n", "- **EEG** measures scalp **potential**. The currents must cross the **skull**,\n", " whose low conductivity *smears* the potential — EEG has lower spatial resolution\n", " but sees all source orientations, including radial (pointing out of the head).\n", "- **MEG** measures the **magnetic field**, which the skull barely perturbs (magnetic\n", " permeability is nearly uniform), so MEG is less smeared. But a radially-oriented\n", " dipole in a spherically symmetric head produces *no external magnetic field* —\n", " MEG is largely **blind to radial sources** and most sensitive to tangential ones.\n", "\n", "Shared assumptions: the **quasi-static** approximation (no wave propagation at these\n", "frequencies, so the map is an instantaneous linear projection), a **fixed source\n", "space** (one dipole per region, with a chosen orientation), and a head conductivity\n", "model. These assumptions are what make the forward map a constant matrix $\\mathbf{L}$\n", "— and what make the *inverse* problem (recovering sources from sensors)\n", "ill-posed, which is precisely why a *generative* forward model like this one is\n", "useful.\n", "\n", "> A units caveat worth flagging (see {doc}`/tutorials/05_forward_models`): the\n", "> ``EEGLeadFieldModel`` helper is tuned for an mV neural source. For MEG from an mV\n", "> source, drive the base ``LeadFieldModel`` with an explicit ``scale`` and\n", "> ``sensor_unit=u.tesla`` so the output lands in tesla." ] }, { "cell_type": "code", "execution_count": 3, "id": "dcfe3b58", "metadata": { "execution": { "iopub.execute_input": "2026-06-19T07:57:05.309598Z", "iopub.status.busy": "2026-06-19T07:57:05.309359Z", "iopub.status.idle": "2026-06-19T07:57:05.420789Z", "shell.execute_reply": "2026-06-19T07:57:05.419809Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "source activity (T x R): (50, 4) -> sensors (T x M): (50, 6)\n", "sensor output unit: mV\n" ] } ], "source": [ "import brainmass\n", "import brainstate\n", "import brainunit as u\n", "import numpy as np\n", "import jax.numpy as jnp\n", "\n", "brainstate.random.seed(0)\n", "\n", "# A tiny EEG forward projection: R=4 source regions -> M=6 sensors.\n", "# The lead field carries the EEG sensor unit per dipole moment, mV / (nA*m).\n", "R, M = 4, 6\n", "L = jnp.asarray(np.random.randn(R, M)) * (u.mV / (u.nA * u.meter))\n", "eeg = brainmass.EEGLeadFieldModel(in_size=(R,), out_size=(M,), L=L)\n", "\n", "# A short window of regional dipole activity (mV NMM observable), (T, R).\n", "T = 50\n", "sources = (jnp.array([1.0, -0.5, 0.8, 0.2]) * jnp.ones((T, R))) * u.mV\n", "sensors = eeg.update(sources) # linear lead-field projection -> sensor space\n", "\n", "print('source activity (T x R):', sources.shape, ' -> sensors (T x M):', sensors.shape)\n", "print('sensor output unit:', u.get_unit(sensors)) # mV, via s @ L" ] }, { "cell_type": "markdown", "id": "492bbb77", "metadata": {}, "source": [ "The four regional sources are linearly mixed by the lead field into six sensor\n", "measurements, unit-correctly (volts out). The same $\\mathbf{L}$-style multiply,\n", "with a tesla-valued lead field, gives MEG. Because the projection is differentiable,\n", "a loss on the *sensor* signal flows back to the neural model — closing the loop from\n", "parameters all the way to the measured signal." ] }, { "cell_type": "markdown", "id": "0338aeb8", "metadata": {}, "source": [ "## Choosing a forward model\n", "\n", "| You are modelling… | Use | Why |\n", "|---|---|---|\n", "| Resting-state fMRI / BOLD FC | ``BOLDSignal`` or ``HRFBold`` | hemodynamics; FC over minutes |\n", "| **Fitting** in BOLD space | ``HRFBold`` | fast, differentiable convolution |\n", "| Biophysically realistic BOLD | ``BOLDSignal`` | full Balloon–Windkessel ODE |\n", "| EEG (scalp potentials) | ``EEGLeadFieldModel`` | volume conduction, all orientations |\n", "| MEG (magnetic fields) | ``LeadFieldModel`` / ``MEGLeadFieldModel`` | tangential-source sensitivity |\n", "\n", "In all cases the forward model is the *generative* link between a neural mass model\n", "and data — and in ``brainmass`` it is one more differentiable module in the chain,\n", "not a separate post-processing step." ] }, { "cell_type": "markdown", "id": "c92cce8e", "metadata": {}, "source": [ "## Key takeaways\n", "\n", "- A **forward model** maps simulated neural activity to a measured signal; the\n", " modality determines the physics and what is observed.\n", "- **BOLD** is hemodynamic, regional, and slow: a Balloon–Windkessel ODE\n", " (``BOLDSignal``) or an HRF convolution (``HRFBold``) that blurs and delays\n", " activity by seconds.\n", "- **EEG/MEG** are electromagnetic and fast: a linear **lead-field** projection of\n", " dipolar source currents, with EEG smeared by the skull and MEG blind to radial\n", " sources.\n", "- Every forward model in ``brainmass`` is **differentiable and unit-safe**, so\n", " signal-space losses backpropagate to neural parameters.\n", "\n", "## See also\n", "\n", "- {doc}`/tutorials/05_forward_models` — running BOLD and EEG/MEG forward models.\n", "- {doc}`/concepts/architecture_overview` — the observation layer's place in the stack.\n", "- {doc}`/concepts/why_differentiable` — why end-to-end differentiability matters.\n", "- {doc}`/reference/index` — the forward-model and observation API.\n", "\n", "## References\n", "\n", "- Friston, K. J., Harrison, L., & Penny, W. (2003). Dynamic causal modelling.\n", " *NeuroImage*, 19(4), 1273–1302.\n", "- Buxton, R. B., Wong, E. C., & Frank, L. R. (1998). Dynamics of blood flow and\n", " oxygenation changes during brain activation: the balloon model. *Magnetic\n", " Resonance in Medicine*, 39(6), 855–864.\n", "- Hämäläinen, M., Hari, R., Ilmoniemi, R. J., Knuutila, J., & Lounasmaa, O. V.\n", " (1993). Magnetoencephalography — theory, instrumentation, and applications to\n", " noninvasive studies of the working human brain. *Reviews of Modern Physics*,\n", " 65(2), 413–497.\n", "- Nunez, P. L., & Srinivasan, R. (2006). *Electric Fields of the Brain: The\n", " Neurophysics of EEG* (2nd ed.). Oxford University Press." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.11" } }, "nbformat": 4, "nbformat_minor": 5 }