From 8bcb03e00d41e4bdf6f515704e37871db670ec2f Mon Sep 17 00:00:00 2001 From: Francois Lanusse Date: Sun, 10 May 2020 19:41:04 +0200 Subject: [PATCH] Delete DifferentiableLensing.ipynb --- notebooks/DifferentiableLensing.ipynb | 497 -------------------------- 1 file changed, 497 deletions(-) delete mode 100644 notebooks/DifferentiableLensing.ipynb diff --git a/notebooks/DifferentiableLensing.ipynb b/notebooks/DifferentiableLensing.ipynb deleted file mode 100644 index af6cfa3..0000000 --- a/notebooks/DifferentiableLensing.ipynb +++ /dev/null @@ -1,497 +0,0 @@ -{ - "cells": [ - { - "cell_type": "code", - "execution_count": 1, - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Populating the interactive namespace from numpy and matplotlib\n" - ] - } - ], - "source": [ - "%pylab inline\n", - "\n", - "import jax.numpy as np\n", - "from jax import jit, grad, vmap" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## Let's define some redshift distributions" - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "metadata": {}, - "outputs": [], - "source": [ - "from jax_cosmo.redshift import smail_nz" - ] - }, - { - "cell_type": "code", - "execution_count": 21, - "metadata": {}, - "outputs": [], - "source": [ - "nz1 = smail_nz(a=1.0, b=2., z0=0.5, zmax=2.)\n", - "nz2 = smail_nz(a=2.5, b=4., z0=1.5, zmax=3.)" - ] - }, - { - "cell_type": "code", - "execution_count": 4, - "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/home/francois/.local/lib/python3.8/site-packages/jax/lib/xla_bridge.py:116: UserWarning: No GPU/TPU found, falling back to CPU.\n", - " warnings.warn('No GPU/TPU found, falling back to CPU.')\n" - ] - }, - { - "data": { - "text/plain": [ - "Text(0.5, 0, 'redshift z')" - ] - }, - "execution_count": 4, - "metadata": {}, - "output_type": "execute_result" - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEGCAYAAABrQF4qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deXhU5dn48e89k30DsgFJCAQIqwICAoKoiCKuuItW21ot2mrVX19ttbb6tn27WH3tW1uq4lJqtaJ1xcpWN1BQJCD7vhMSSAhZyL49vz/OBIeYZZLM5Mxyf65rLiZnmXOfHHLPM888537EGINSSqng5bA7AKWUUr6liV4ppYKcJnqllApymuiVUirIaaJXSqkgF2Z3AC1JTk42AwYMsDsMpZQKGGvXrj1mjElpaZ1fJvoBAwaQk5NjdxhKKRUwRORAa+u060YppYJcuy16EXkRuAwoMMac1sL6B4Bvub3ecCDFGHNcRPYDJ4AGoN4YM95bgSullPKMJy36+cDM1lYaYx43xowxxowBHgKWG2OOu20yzbVek7xSStmg3Ra9MWaFiAzw8PVuBF7tSkBKKeVNdXV15ObmUl1dbXcoXhEVFUVGRgbh4eEe7+O1L2NFJAar5X+322IDLBMRAzxrjJnXxv5zgDkAmZmZ3gpLKRXicnNziY+PZ8CAAYiI3eF0iTGGoqIicnNzycrK8ng/b34Zezmwslm3zRRjzFjgYuAuETmntZ2NMfOMMeONMeNTUlocIaSUUh1WXV1NUlJSwCd5ABEhKSmpw59OvJnoZ9Os28YYk+f6twB4G5jgxeMppZRHgiHJN+nMuXgl0YtID+Bc4F23ZbEiEt/0HJgBbPbG8fxdbX0ju46eYPGmfJ5bsZcDRRV2h6SUCmGeDK98FTgPSBaRXOBRIBzAGPOMa7OrgGXGGPeM1ht42/XuEwb80xizxHuh+6evDhbznRe/pKy6/uSyP324i8euGcWlo/raGJlSyt+tWLGC++67j40bN7JgwQKuvfZar7yuJ6NubvRgm/lYwzDdl+0FRnc2sEB0sKiS2/+eQ4+YcH45ayTZqfFERzi5/18buOuf61i9rz8PXzqcyDCn3aEqpfxQZmYm8+fP54knnvDq6/plCYRAVFpZx63zv6S+0fC3705gcGrcyXWvzTmLPyzZzvOf7aOmrpHHrh1lY6RKqe60f/9+Lr74Ys4++2xWrVpFeno6f/vb37j00ktPbrNp0yb27t1LU40vh8O7RQs00XtBXUMjd7ycw6HjVfzjtlOTPEBEmIOfXzYCp1N4dvleZp2RxuRByTZFq1To+uV7W9iaV+bV1xyRlsCjl49sc5tdu3bx6quv8txzz3H99dfz8ccfs379egDmzp3L8uXL6d+/v1fjcqe1brzg7XWH+WLvcX579elMHJjU6nb3TR9CZmIMP3trE9V1Dd0YoVLKTllZWYwZMwaAcePGsX//fgBWrlzJ888/z4svvujT42uLvosaGw3zPt3L8L4JXDM2vc1toyOc/Paq07n5hdX8+aNdPHDRsG6KUikFtNvy9pXIyMiTz51OJ1VVVeTn53PbbbexcOFC4uLi2ti767RF30XLdxayu6CcOedkeTS+9ezsZK4em86zy/eyLd+7HyGVUoGhrq6O66+/nscee4whQ4b4/Hia6Lvo2RV76NsjistGpXm8zy8uHUF8VBhPLN3hw8iUUv5q1apVrFmzhkcffZQxY8YwZswY8vLyWLNmDRkZGfzrX//ijjvuYORI73wC0a6bLtiYW8IXe4/z8CXDCXd6/p7ZKzaCb581gKc+2sW+YxVkJcf6MEqllJ0GDBjA5s1f3yt6//33c//997e4bVpaGrm5uV6PQVv0XfDcp/uIjwxj9oR+Hd73W5MyCXc4mL9ynw8iU0qpr2mi76Tc4koWbcrnpomZxEd5Xi60SWp8FJePTuNfa3MprarzQYRKKWXRRN9J72/Mp6HRcPOkzo99vXXKACprG3htzUEvRqaUUqfSRN9JH2w7ysi0BPolxnT6NU5L78HErET+vuoA9Q2NXoxOKaW+pom+E4rKa1h7oJgLhvfu8mt97+wsDpdUsWzrUS9EppRS36SJvhM+2l5Ao4ELR3Q90V8wvDcZvaJ5ZfUBL0SmlFLfpIm+Ez7YdpS+PaIYmZbQ5ddyOoSrzkjn8z1FFJwIjjktlVKd8+STTzJixAhGjRrF9OnTOXDAOw1ATfQdVF3XwIqdx7hgeG+vzVoza0wajcb6glcpFbrOOOMMcnJy2LhxI9deey0/+clPvPK6mug7aNWeY1TVNXCBF7ptmgxOjWd43wQWbsjz2msqpfzD/v37GT58ON///vcZOXIkM2bM4PDhwyfviB0zZgxOp5MDBw4wbdo0YmKsAR6TJk3y2s1TemdsB/1nawGxEU4mDUz06uteMTqNx5Zs59Dxyi6N5FFKtWHxg3Bkk3dfs8/pcPHv29ykM2WKX3jhBS6++GKvhKgt+g5obDR8uO0o5w5N8fosUZePtqYZ1Fa9UsGno2WKX375ZXJycnjggQe8cnxt0XfApsOlFJyo8cpom+YyesUwrn8v3tuQx13TBnv99ZVStNvy9pWOlCn+4IMP+M1vfsPy5ctP2a8rtEXfAZ/tPgbAuUNSffL6V4xOY/uRE+w8esInr6+U8g+tlSn+6quvuOOOO1i4cCGpqd7LM+0mehF5UUQKRGRzK+vPE5FSEVnvejzitm6miOwQkd0i8qDXorbJ2gPFDEqJJTE2wievf8npfXEILFyv3TdKBbPWyhQ/8MADlJeXc9111zFmzBiuuOIKrxzPk66b+cBfgJfa2OZTY8xl7gtExAnMBS4EcoE1IrLQGLO1k7HaqrHRsPZAMTNH9vHZMVLiI5k8KJnFm/O5/6KhPjuOUqr7dKRM8QcffOCTGNpt0RtjVgDHO/HaE4Ddxpi9xphaYAEwqxOv4xf2FJZTWlXHuAG9fHqcC4ansqewgv3HKnx6HKVU6PBWH/1ZIrJBRBaLSNOUKOnAIbdtcl3LWiQic0QkR0RyCgsLvRSW9+QcKAZgfH/fJvrprvo5H24v8OlxlFKhwxuJfh3Q3xgzGvgz8I5reUu3jZrWXsQYM88YM94YMz4lJcULYXlXzv5ikmIjfD4bVL/EGIb0juPDbVrkTClvMabV1BNwOnMuXU70xpgyY0y56/kiIFxEkrFa8O5TL2UAAfst49oDxxnbv5fXyh60Zfrw3ny57zhl1TohiVJdFRUVRVFRUVAke2MMRUVFREVFdWi/Lo+jF5E+wFFjjBGRCVhvHkVACZAtIlnAYWA2cFNXj2eHwhM17C+q5MYJmd1yvOnDUnn6kz2s2FnYoUnHlVLflJGRQW5uLv7YJdwZUVFRZGRkdGifdhO9iLwKnAcki0gu8CgQDmCMeQa4FviBiNQDVcBsY7111ovI3cBSwAm8aIzZ0qHo/MTapv55H38R2+SMzF70ignnw20FmuiV6qLw8HCysrLsDsNW7SZ6Y8yN7az/C9bwy5bWLQIWdS40/7H2wHEiwhyclt6jW47ndAjThqby0Y4C6hsaCXPqfW1Kqc7TDOKBnAPFjErv4fX6Nm2ZPrw3JZV1fHWopNuOqZQKTpro21Fd18Dmw6U+Hz/f3DlDkglzCB/o6BulVBdpom/HxtxS6hoMZ/b3blni9sRHhTNxYCIf63h6pVQXaaJvx7qD1hexY318o1RLzslOYefRco6U6hSDSqnO00Tfjq15ZaT3jPZZIbO2TM22bhxrqpqplFKdoYm+HdvyyxjeN96WYw/rE09yXASf7QqO8b9KKXtoom9DdV0DewrLGdE3wZbjOxzC2YOT+Wz3MRobA/+uPqWUPTTRt2HHkRM0GhiRZk+iBzg7O4Vj5bVsP6KTkSilOkcTfRu25ZcBMNymFj3A1OxkAD7V7hulVCdpom/D1vwy4iLD6NcrxrYYeidEMaR3nH4hq5TqNE30bdiaZ30R63D4vmJlW6Zmp7B633Gq6xpsjUMpFZg00beisdGw/cgJW7ttmpydnUxtfSNr9ndmoi+lVKjTRN+KQ8WVlNfU2zbixt3ErEQinA4+3aXdN0qpjtNE34qtedYXsXaOuGkSExHGuP69NNErpTpFE30rtuaX4RAY0tuem6WaOzs7mW35ZRyvqLU7FKVUgNFE34pt+WUMSokjKrz7ShO3ZdLAJAC+2FtkcyRKqUCjib4VW/PK/KLbpsmojB7ERjj5fI8meqVUx2iib0FJZS15pdV+MeKmSbjTwZlZiazao/30SqmO0UTfgq2uO2L9YcSNu8mDkthTWEFBmZYtVkp5ThN9C5pG3PhTix7grIFWOYTPtZ9eKdUB7SZ6EXlRRApEZHMr678lIhtdj1UiMtpt3X4R2SQi60Ukx5uB+9Kuo+UkxUaQEh9pdyinGJGWQEJUmPbTK6U6xJMW/XxgZhvr9wHnGmNGAb8G5jVbP80YM8YYM75zIXa/PYXlDEqNszuMb3A6hIkDk1iliV4p1QHtJnpjzAqg1XvvjTGrjDHFrh+/ADK8FJstjDHsLixnUIr/JXqAswYmcfB4JbnFlXaHopQKEN7uo78NWOz2swGWichaEZnT1o4iMkdEckQkp7DQvpK8xytqKamsY7AftugBJg+2xtNr941SylNeS/QiMg0r0f/UbfEUY8xY4GLgLhE5p7X9jTHzjDHjjTHjU1JSvBVWh+0uKAdgUEqsbTG0ZUhqPImxEfqFrFLKY15J9CIyCngemGWMOZmBjDF5rn8LgLeBCd44ni/tKawA8NsWvcMhnDUwic/3FGGMTi+olGpflxO9iGQCbwG3GGN2ui2PFZH4pufADKDFkTv+ZHdBOdHhTtJ6RNsdSqsmDUoiv7SaQ8er7A5FKRUAwtrbQEReBc4DkkUkF3gUCAcwxjwDPAIkAX8VEYB61wib3sDbrmVhwD+NMUt8cA5etaewnIEpsbZPNtKWSVmJgFX3JjPJvtmvlFKBod1Eb4y5sZ31twO3t7B8LzD6m3v4t90F5Yzr38vuMNo0ODWOpNgIvthXxPVn9rM7HKWUn9M7Y91U1TZwuKTKb/vnm4gIEwcmsnqvzjillGqfJno3ewqbRtz4d6IHq2zx4ZIqDh3X8fRKqbZponfTlOj9vUUPMDFL69MrpTyjid7NnoJyHAIDkv3/C87s1DgSYyNYvU+7b5RSbdNE72ZPYQWZiTFEhvnHrFJtcTiECQMStUWvlGqXJno3uwv8t8ZNSyYNTCS3uErr3iil2qSJ3qWh0bDvWEVA9M83meiaR1ZH3yil2qKJ3uXQ8UpqGxoDqkU/tHc8PWPCWb1Pu2+UUq3TRO9ycmhlALXov+6n1xa9Uqp1muhdmqpWDg6gFj1Y4+kPHq8kr0Tr3iilWqaJ3mVvYQXJcRH0iAm3O5QOmTjQqnuj3TdKqdZoonfZX1TBgCT/rEHflmF9rHlk9QtZpVRrNNG7HDxeGZCVIJ0OYUJWot44pZRqlSZ6oLqugfzS6oBs0YNVDmHfsQoKyqrtDkUp5Yc00WO15gH6B2CLHqwvZAG+0Fa9UqoFmuiBA0VNiT4wW/Qj0hKIjwxjtZZDUEq1QBM9cKDImid2QIC26J0OYfyAXlr3RinVIk30WCNuEqLC6BkTYXconTZxYBJ7CisoPFFjdyhKKT+jiR6r62ZAcmB22zSZ6JpH9kvtp1dKNdNuoheRF0WkQEQ2t7JeROQpEdktIhtFZKzbupkissO17kFvBu5NB4oqA7Z/vslp6T2IjXDqjVNKqW/wpEU/H5jZxvqLgWzXYw7wNICIOIG5rvUjgBtFZERXgvWFuoZGDpdU0T8xMPvnm4Q7HYwboPPIKqW+qd1Eb4xZAbSVPWYBLxnLF0BPEekLTAB2G2P2GmNqgQWubf3K4eIqGhpNwA6tdDcxK5EdR09wvKLW7lCUUn7EG3306cAht59zXctaW+5X9rtG3AR61w18PZ5eh1kqpdx5I9FLC8tMG8tbfhGROSKSIyI5hYWFXgjLM003SwXq0Ep3ozJ6EB3u1GGWSqlThHnhNXKBfm4/ZwB5QEQry1tkjJkHzAMYP358q28I3rb/WCXR4U5S4iO765A+E+50uMbTaz89ANVlULQbygug/Ag0NkBcb4hLhV5ZEJdid4RKdQtvJPqFwN0isgCYCJQaY/JFpBDIFpEs4DAwG7jJC8fzqgNFFfRPikGkpQ8ggeesQUn8YckOisprSIoL/DevDivYBtv+DXs+hENfgmlofds+oyD7Qhh2KaSP674Ylepm7SZ6EXkVOA9IFpFc4FEgHMAY8wywCLgE2A1UAre61tWLyN3AUsAJvGiM2eKDc+iSA8crGZQS+P3zTU720+87ziWn97U5mm5iDOz9GFb9GfZ8ZC3rOwbOvs9K4HF9rFa8IwwqCqwW/pGNsOsD+Oz/4NP/hf5TYOqPYdB0CJI3faWatJvojTE3trPeAHe1sm4R1huBX2poNBwsqmT6sFS7Q/Ga09N7EBNh9dOHRKI/vA7e/y/IW2d1y0x/BM64xUrsLUlw/U6yL4Sp/wVVJbDhVVj5FLx8DaSNhSuegj6nd985KOVjIX1n7JGyamobGgOyDn1rrH76xOD/Qra6FBY9AM+dD2V5cPlTcN8mK3m3luRbEt0TJv0A7t1gvUZpLsw7Dz75PdTrMFUVHEI60X9dzCx4um4AzhqYxM6j5RwrD9K6Nwe/gLmT4MvnYML34e4vYdx3IKwL30mERVivcddqGHk1fPI7eP58KN7vtbCVskuIJ3praGVmgN8V29ykpnlkg230jTHwxdMw/1Irqd/+IVzyOET18N4xYhLhmudg9j+h5BA8Nx0Orvbe6ytlg5BO9PuLKgh3Cmk9o+0Oxaua6t4EVfdNXRW88T1Y8iBkXwRzPoEMH46UGXap9UYSlQB/vxw2veG7YynlYyGd6A8WVdKvVwxOR3CNsgi6fvqqEvjH1bDlbZj+KMx+xepb97XkwVayzxgPb94Ga17w/TGV8oGQTvS5xVVkBFm3TZNJA5PYVVAe+PXpTxyxumpy18A1z1tDILtz+GNMItzyNgyZCe//GNb9o/uOrZSXhHiir6Rfr+Dqtmly1iBrPP3ngdyqLzkEL14Ex/fBTa/B6dfaE0dYJFz3d2uM/cIfwfpX7YlDqU4K2URfXlNPcWUdGb2Cs0V/WloC8VFhfL7nmN2hdM6JI/DSFVBZDN9ZCIOn2xtPeJTVZZR1Drz7Q9jut7eHKPUNIZvoc4utETcZQdqiD3M6mJiVxMrdAdiirzgGL82y7mC9+U2rj9wfhEfDja9ad92+eTvkb7Q7IqU8ErKJ/tDxKgD6BWkfPcCUwUkcPF7JIVeFzoBQVQIvXQnFB6zumn5n2h3RqSJirWQf3RNenW198lDKz4Vsog/2Fj3A5EHJAHy+J0Ba9fW18PotULjd6iYZcLbdEbUsvg/cuMB6U3p1NtQG0BupCkkhnOiriA53khQbYXcoPjOkdxzJcRGsCoR+emPgvXtg3wqY9Rf7++Tb03eUNQoobz0sfsDuaJRqU8gm+kPHK8noFR005YlbIiKcNSiZlXuKsGrP+bHlj1nFxc77GYyebXc0nhl2CZzzAHz1so7EUX4tZBN9bnFVUHfbNJkyKInCEzXsLii3O5TWbXrDqi0z5ltw7k/sjqZjznsQBky1xtgX7rA7GqVaFMKJvjKov4ht0tRPv8pf++nzN8C7d0PmZLjs/wKvFrzDaXXhRMTC69+B2gq7I1LqG0Iy0ZdW1VFWXR8SLfrMpBgyekWzcrcf9tNXFMGCm627T6//u1VBMhDF94Gr51lfIi/7hd3RKPUNIZnom0bc9AvSm6WamzwoiS/2FtHQ6Ef99A318MZ3ofwo3PByx2rI+6NB58NZd0HOC7D7A7ujUeoUIZnom8bQB+tdsc1NGZxMWXU9W/JK7Q7lax/9yhphc/mfIH2s3dF4x/m/gOShVldUVbHd0Sh1Ukgm+lAYQ++uqe7Np7v8pPtmxxJY+ScY/z0Y0+ZMlYElPAquesa6o3dRgH2prIJaiCb6KuIiw+gZE253KN0iNT6K4X0T+HRXod2hWIXK3rkT+oyCi35ndzTelz7WGnK56XXYutDuaJQCPEz0IjJTRHaIyG4RebCF9Q+IyHrXY7OINIhIomvdfhHZ5FqX4+0T6Izc4uAfQ9/cOdnJrD1QTEVNvX1B1NfCG7da/fPXzbdawMHonPutN7JF91t3zypls3YTvYg4gbnAxcAI4EYRGeG+jTHmcWPMGGPMGOAhYLkxxn0eu2mu9X5RncoaQx8a/fNNpmanUNdgWL3PxmGWH/3aqis/6y+QNMi+OHzNGQ5XPAUVhfDBo3ZHo5RHLfoJwG5jzF5jTC2wAJjVxvY3An57m6AxJmRulnI3fkAvIsMcrNhpUz/9no9h1VMw7lYYeaU9MXSntDNg0g9h7Xw4sMruaFSI8yTRpwOH3H7OdS37BhGJAWYCb7otNsAyEVkrInNaO4iIzBGRHBHJKSz0XV9ySWUd5TWhMYbeXVS4k4kDk+zpp684Bm/fCSnD4KLfdv/x7TLtZ9AzE967F+oDfKYvFdA8SfQtdWS3NiD7cmBls26bKcaYsVhdP3eJyDkt7WiMmWeMGW+MGZ+SkuJBWJ2TWxz85Ylbc052MnsKKzhcUtV9BzXm6+GG17wAESH0e4+IhUv/CMd2wqdP2h2NCmGeJPpcoJ/bzxlAXivbzqZZt40xJs/1bwHwNlZXkG0OhdjQSnfnDLHeQD/rzlb9mudh52K48FfQ57TuO66/yL4ATrsGPvsjHN9rdzQqRHmS6NcA2SKSJSIRWMn8G+PGRKQHcC7wrtuyWBGJb3oOzAA2eyPwzvp6DH0ItSxdslPj6J0QyYruGk9/bJdVEmDwBTDxju45pj+a8T/WF7RLHrI7EhWi2k30xph64G5gKbANeN0Ys0VE7hSRO902vQpYZoxxr+rUG/hMRDYAXwLvG2OWeC/8jsstriIhKowe0aExht6diDA1O4WVu4/5vhxCQx28NccaQjlrbuAVK/OmhDSryuXOJbBjsd3RqBAU5slGxphFwKJmy55p9vN8YH6zZXuB0V2K0MusOvSh15pvMjU7mTfW5rLpcClj+vX03YFWPAF56+D6l6yiX6Fu4p1W3frFP4GB51nzzyrVTULuzti8kuqQ7J9vMjU7BRH4ZEeB7w6SmwMrHofRN8KItkbihhBnOFzyBJQctPrrlepGIZXojTEcLqkirWfoJvrE2AjG9OvJx9t9lOhrK+HtOyC+L1z8mG+OEaiyplpfzK78k5XwleomIZXoy6rrKa+pJz2EEz3AtKGpbMgtpfCED8Z2f/grKNoNV86FqB7ef/1Ad8EvAYH/PGJ3JCqEhFSiz3ONHw/lFj3A+cOs2u9e777ZtwJWPw0T5lj90OqbevaDs++DLW/D/pV2R6NCRIgm+iAtpuWhkWkJpMZH8rE3E33NCXjnLkgcCBf8t/deNxhNvgcSMmDJg9DYYHc0KgSEVKJvuiM0PYS/jAVrmOW0oal8uvMYdQ2N3nnRZT+H0kNw5dPWHaGqdRExMONXcGQjfPUPu6NRISDkEn2E00FybKTdodhu2rBUTtTUk7PfCzMh7f7QKt41+W7InNT11wsFI6+GzLPgo/+B6jK7o1FBLqQSfV5JNX17RuFwhPDNOy5nZycT7pSud99Ul8LCeyB5CEx72DvBhQIRuOg3ViljHW6pfCzEEn0VaT1Cu9umSVxkGBOzkvioq8Mslz4MJ/KsLhu9Cahj0sfBqBvg87k63FL5VEgl+sPFVSHfP+9u2rBUdheUc+h4ZedeYNcHVh/z5Hsgwy/mlAk80x8BccAHv7Q7EhXEQibR1zU0cvREdcgPrXTXNMzyg21HO75zdSm8dw8kD4XztFhXp/XIgMk/gs1vwKE1dkejglTIJPojpdUYA+khPrTSXVZyLINT41i2pROJfunDcCLf1WWjv9MumXIvxPWGpT+z6vcr5WUhk+j1ZqmWXTSyN6v3FXG8otbznZq6bKbcCxnjfBdcqIiMg/N/DrlfwtZ37I5GBaGQSfSHNdG3aObIvjSaDnTfNHXZpAzTLhtvGvMtSB0JH/y3TjuovC5kEn1Tiz7U69w0d1p6Auk9o1m25YhnOzR12cz6K4Tp/Qhe43DCjF9D8X748jm7o1FBJmQS/eGSapJiI4gKd9odil8RES4c0ZsVu45RUVPf9sbaZeNbg6fDoOmw4g9Qebz97ZXyUMgk+rwQL0/clpmn9aG2vpHlO9uYS7aqBBb+SLtsfG3G/1h1g1Y8bnckKoiETKK36tDr6JCWnDkgkcTYCJZsbqP7ZunDUH4UrtQuG5/qPQLOuNnqvinaY3c0KkiERKI3xpBXUkV6z9CdQrAtTodwwfBUPt5eQG19C0XOdi6D9S9bXTbp2mXjc9MeBmeE9cWsUl4QEom+tKqOytoGbdG34aKRfThRU8+qPcdOXVFV7BplM9ya4Fr5Xnwf601120I4+IXd0agg4FGiF5GZIrJDRHaLyDf+2kXkPBEpFZH1rscjnu7bHQ7riJt2TRmcTFxkGIs25Z+6YvGDUF6gXTbdbfLd1nSMSx/Wm6hUl7Wb6EXECcwFLgZGADeKyIgWNv3UGDPG9fhVB/f1qcPFOoa+PVHhTmaM7M3izUeoqXdNhrH9fdi4AKb+F6SPtTfAUBMRC+f/Ag7nwJa37I5GBThPWvQTgN3GmL3GmFpgATDLw9fvyr5eo3fFembWmHROVNfzyY5CqCiC9+6FPqfDOQ/YHVpoGj0bep9u9dXXVdsdjQpgniT6dOCQ28+5rmXNnSUiG0RksYiM7OC+iMgcEckRkZzCwjaG+XVCXmk1EWEOkuMivPq6wWbKoCSSYiNYuCEPFv2XNaTyqmchTH9vtnA4rZr1JQetuXiV6iRPEn1Ls3Q07zRcB/Q3xowG/gw0FezwZF9roTHzjDHjjTHjU1JSPAjLc4dLqkjvGY2ITjjSljCng0tO70vk9retyavPexB6j2x/R+U7A8+FIRfDiv+Fcu82gFTo8CTR5wL93H7OAPLcNzDGlBljyl3PFwHhIpLsyb7dIa+kir49dMSNJ67JdvCoPM/xXqNhyn12hyqadY0AABjDSURBVKPAKo1QXwWf/NbuSFSA8iTRrwGyRSRLRCKA2cBC9w1EpI+4mssiMsH1ukWe7Nsd9K5YDxnD6LUPEyEN/CH6PnCG2R2RAkjOhjNvt+blPbrV7mhUAGo30Rtj6oG7gaXANuB1Y8wWEblTRO50bXYtsFlENgBPAbONpcV9fXEiralraKTgRA1p2qJvX84LyN6PWD7gHt7YH0VxR0oXK98696cQmQDLdLil6jiPmmyu7phFzZY94/b8L8BfPN23Ox0tsyYc6ast+rYd2w3LfgGDptPv/Lup376S9zflc/Ok/nZHpgBiEq1kv/Qh2LUMhlxkd0QqgAT9nbH5pdawNO2jb0NDHbx1u3VD1Ky5jEjrwZDecbyxNtfuyJS7Cd+H5CGw5CGo109bynNBn+h1DL0HPvkd5H0Flz8FCX0REa4f34/1h0rYceSE3dGpJs5wuOh3cHwPfPms3dGoABL0if6ItujbdmAVfPqkVTFxxBUnF189NoNwp/DamkNt7Ky6XfYFkH0RLP+DVZpCKQ8EfaLPL60mPjKM+Khwu0PxP1Ul8NYd0GsAzPz9KasSYyOYMaIPb3+V+3VJBOUfLvot1FXCh7+yOxIVIII+0eeVVNFHW/PfZAz8+/9B2WG4+jmIjP/GJtef2Y/iyjr+s9XD+WRV90geDBPvhK9ehsNr7Y5GBYCgT/T5pdU64qYlX71sFcs6/2Hod2aLm5w9OJn0ntHafeOPzv0pxKXC+/dDYwtzCCjlJgQSfZWOoW+ucCcs/glkndPm3a9Oh3DtuAw+232M3OLKbgxQtSsqwZp2MG8dfPWS3dEoPxfUib6mvoFj5bX07aEt+pPqquGN70F4NFw1zyqc1YbrxmcA8HqODrX0O6dfB/2nWNUtdTJx1YagTvQnR9zozFJfW/ozOLoJrnwaEvq2u3lGrxjOG5LCP1cf1C9l/Y0IXPIEVJfBh7+0Oxrlx4I60eeVWIk+TVv0lk1vQM4LMPmeDt1ZeeuULI6V1/DvDfntb6y6V+8R1heza/8OuTl2R6P8VFAn+vxS62YpbdEDx3ZZE4n0mwjTH2l/ezdTs5MZnBrHiyv3YbTOiv+Z9pA17eDCe6y7nJVqJsgTvbboAaithNe/A84IuPZv1h2WHSAifG9KFlvyylizv9hHQapOi4yHS5+Agi3weYslp1SIC/JEX0XPmHCiI9r+wjGoNY2XL9hqjZfv0eIEX+266ox0esaE8+Jn+7wcoPKKYZfCsMvgk9/D8b12R6P8THAn+pJq+iSEeLfNmuetCb7Pe8i6fb6ToiOc3DQhk2Vbj3DouA619EuXPA6OcPj3j7WUsTpFUCf6vNLq0C5mdvALWPIgDJnplQm+bzmrPw4R5q/a3/XYlPclpMEFj8Lej2H9P+2ORvmRoE70+aUhPIXgiSPw+rehZ6Y1wbej65e6b49orhidxj9XH+RYeY0XglReN/42yJxslTIu6/ZZO5WfCtpEX1XbQEllXWi26OuqYMFNUFMON7wC0T299tJ3nT+YmvoGnluh/cB+yeGAWX+BhlprFI524SiCONHnNQ2tDLUWvTHWH/jhdXDNc9Y4ay8alBLHFaPTeOnzAxRpq94/JQ2CC/4bdv8H1r9idzTKDwRtos8vaapDH2It+s/+CJteh/N/bo3E8IG7z8+mpr6BeZ9qq95vTZhjlUdY8hCUavmKUBe0ib6pRZ8WSjdLbXvPqlF+2rUw9b98dpjBqXFcPjqNl1Zpq95vNXXhNDbA23da/6qQ5VGiF5GZIrJDRHaLyIMtrP+WiGx0PVaJyGi3dftFZJOIrBeRbrtHu6nOTcjUoj+0Bt68HdLHWX/gIj493I/Oz6ZaW/X+LXEgXPIH2P8prPyT3dEoG7Wb6EXECcwFLgZGADeKSPOO333AucaYUcCvgXnN1k8zxowxxoz3QsweyS+tIjkugsiwELhZqmgPvHqDdRv8Ta9ZlSl9bHBqHFeOSWf+yv1awtifjfkWjLwKPv4N5OokJaHKkxb9BGC3MWavMaYWWADMct/AGLPKGNN0b/wXQIZ3w+y4vJLq0OifryiCV661voS9+U2ITe62Qz9w0VBE4HeLt3fbMVUHicBlf4S4PvDmbVCjk72HIk8SfTrgPsVQrmtZa24DFrv9bIBlIrJWROa0tpOIzBGRHBHJKSws9CCstuWXhsAUgjUn4JVrrPHSN71mjbboRmk9o7njnEG8vzGfL/dpPXS/Fd3LGoFVckCHXIYoTxJ9S529Lf5PEZFpWIn+p26LpxhjxmJ1/dwlIue0tK8xZp4xZrwxZnxKSooHYbUtv7Q6uGeWqquGV2+E/I1w3XzoN8GWMO48dxB9e0Txq39vobFRE4jf6j/ZGom15S1Y/Yzd0ahu5kmizwX6uf2cAXzjljsRGQU8D8wyxhQ1LTfG5Ln+LQDexuoK8qnymnpOVNcH71yxDfXWLFH7P7UmEBl6sW2hREc4efDiYWw+XMYba3UYn1+b8v9g6KWw7Odw4HO7o1HdyJNEvwbIFpEsEYkAZgML3TcQkUzgLeAWY8xOt+WxIhLf9ByYAWz2VvCtORLMN0s1NsA7P4Ad78PFf4DRN9gdEVeMTmNsZk8eW7Jdh1v6M4cDrnraKovxr+9YZTJUSGg30Rtj6oG7gaXANuB1Y8wWEblTRO50bfYIkAT8tdkwyt7AZyKyAfgSeN8Ys8TrZ9FMUx36oKtc2dgA7/zQdUPUL2DiHXZHBFj16n979emUVdfxyMItdoej2hLVA2542fp+Z8FN1lwFKuiFebKRMWYRsKjZsmfcnt8O3N7CfnuB0c2X+1pTog+qUTeNDfDuXVbJ4Wk/h3PutzuiUwzrk8C907N5YtlOLj09n0tOb38+WmWT3iOtuQleuxneuROune+VonfKfwXl1W0qf9C7R6TNkXhJQ73Vkt/wKkx7GM7teslhX7jz3EGcnt6DX7yzWbtw/N3wy2DGr2Hru/DRr+yORvlYUCb6I2VBdLNUXbVVbnjjAmvUxLk/sTuiVoU5HTxx3WhOVNfz83c26/yy/u6su2Hcd636SDl/szsa5UNBmejzS6uDYwx9zQn453WuL14f98rkIb42tE88P54xhMWbj/C3lfvtDke1RQQueQIGX2hNN7npDbsjUj4SlIn+SGk1fRICvH/+xFH4++Wwf6U1ccjEVu818ztzpg7kwhG9+e2ibXojlb9zhsP1L1nj7N+aA9sXtb+PCjhBmejzSqoCu2plwTZ4fjoU7oDZr8Do2XZH1CEOh/C/148mMzGGH76yjqNl1XaHpNoSEQM3LoC+o+Ff34U9H9kdkfKyoEv0FTX1lFXXB27XzZ6P4YUZ1gxBty6y9WaorkiICufZW8ZRWVvPnS+vpbpOy+T6tagEq1ZScjb88wbY/r7dESkvCrpEf6SsaWhlgCV6Y2DlU/Dy1dAjA27/ENLOsDuqLsnuHc+T149m/aESfvjKOuoaGu0OSbUlJhG+8x70GQWv3QIbFtgdkfKS4Ev0J2+WCqA++toKq6TBf34Bwy6D25ZBz37t7xcAZp7Wl99ceTofbS/gx69voEHr4fi3mET49rswYAq8fQd8/lctghYEPLphKpDklQTYzFJHt8Abt8GxHXDBL2HKvT6fNKS73TQxk7LqOn6/eDtxkWH85srTcDiC6xyDSmQc3PQvq6zx0oegcLs1Oicswu7IVCcFXaJvatH39vfyB8bAl89ZBaaiesDNb8GgaXZH5TN3njuIsqo6/vrJHipr63n82tFEhAXdB8rgER4F1/8DPvo1fPYkFO22fo5Nsjsy1QlBl+jzy6pJjI0gKtyPb5Yqy4P37oNdSyH7Ipg1F+K6XprZ3z1w0VDiosL4w5IdHCuv4embx5EQFW53WKo1Dgdc8CikjrDKbzw7Fa6eBwPOtjsy1UFB16SyxtD7aWveGFj7d5g7EfatgJmPWROGhECSB6v42Q/PG8yT149m9d7jXP/M5xwoqrA7LNWeUdfBbUshLArmXwYf/Q801NkdleqAoEv0fjuGvmA7vHQFvHePNV75Byth0p1B1x/viavHZvC3W88kr6SKS5/6jHfXH7Y7JNWetDPgjhXWHLQrHocXLoS89XZHpTwUdIn+SJmflT+oKobFP4WnJ0P+Bmv+zm8v7PZp//zN1OwUFt07laF94rl3wXoe+NcGyqq1lejXIuPgyrnWjGalufDcNFjykM5DGwCCKtFX1TZQUlnnH+WJayth5Z/gqbHw5TwY9x340Vcw/ntaEtYlo1cMr82ZxN3TBvPGulzOf2I5b63L1WJo/m7kVXD3Ghh3K3zxNPx5nDWwoL7W7shUK4Iq4zTdLGVrH31tJXzxDPxpNPznEesj75zlVkteRyx8Q5jTwf0XDeXdu6aQ3iuaH7++gRue/YKc/Vojx69F94LLnrRu7EscBIvuh7+Mg69e0YTvh4Iq0ee7xtD3taOP/sRR60uqP46AJT+F5CFw6xK45S3oO6r74wkwozJ68vYPJvPYNaezp7Cca5/5nNnzPmfV7mPawvdnGeOsUh3fehOiesK7P4T/O93qx68oan9/1S2Canhlt88s1dgI+5bDupdg+7+tkQhDL4HJd0PmWSH5RWtXOBzCDWdmcvnoNF798hDPLt/DTc+vJjs1jhvO7MfVYzNIjNWbdvyOCGRfAIPOtwqifTHXavQs/4NVq2nUbBh8gd5wZSPxx9bS+PHjTU5OTvsbNjP34908vnQH2341k+gIH42jNwbyvrJm5tnyFpQctFoyo26ACXMgebBvjhuCqusaeOerwyxYc4j1h0oIdwpTBiczY0QfLhieSqq/DqNVVgXWtfOtGveVx6yunuwZMGQmDJ5u3SSovEpE1hpjxre4LpgS/c/f2cT7G/P56pEZ3g2ougz2f2q1VnYts5K7IwyyzoUxN1n1acI16fjSjiMneGPtIZZuOcrB49aE1sP6xDNpYBITshI5I7MnfRKiEP0U5V8a6qy/m81vWn87VcUgTkgbY33q7T8Z0sZCfB/9BNxFXU70IjIT+BPgBJ43xvy+2Xpxrb8EqAS+a4xZ58m+Lelsor9t/hrySqtZfO/UDu97UkMdFO2BIxshdw0c+hKObALTAOGxkDUVhl9hfSSNSez8cVSnGGPYebScD7Yd5fM9Raw9UEyVqwRyYmwEI/omkN07joHJsQxIjiUzMYY+PaKCY1rJQNfYYP1N7foPHFgFh3OsctwAsSlW1cyUoZA02CqX3DMTEtKtyVFUu9pK9O320YuIE5gLXAjkAmtEZKExZqvbZhcD2a7HROBpYKKH+3pNfml12+WJG+qh9oTVqqg4BhWFVjmC4v1QcgCK9sKxndDoGs8dHgvpY+Hs/2fVocmYoP2MNhMRhvaJZ2ifeO6aNpi6hkY2Hy5l0+FSthwuY0t+KQu+PHQy+TdJjougd0IUSXGRJMdG0Cs2gh7R4SREhREfFU5sZBixkU5iIsKICncQFe4kKtxJhNNBRJiDCKeDMKcQ5hD91NBZDidkTrIeYM2HnPeVdX/JkY3W48AqqK9y20kgvi/E94bYVOsu8pgkq7s0uhdExn/9CI+B8GjrERYFzggIi7T+dYT2G70nX8ZOAHYbY/YCiMgCYBbgnqxnAS8Z6+PBFyLSU0T6AgM82Ndr/rf4HpIqG+DP4VZfemOd1UJvqLWGPZ7yH8hNWJTVeuiVBdkXWrU9eo+AlOHgDKrvq4NOuNPBGZm9OCOz18llxhiOltWw71gFh4oryS+pJr+0ioITNRSV17CnoJziyloqazs3GUqYQ3A2ezik6WH1QDhEEKw3JnEtE1zPXa8jrm1OIS0+PXWTLrzR+Odb1HDX4wYkupFkU0RG42FSGwtJbSykd3UhvSqLSTQ76WXWEG/KiKC+Q0doRKjHSQNOGnFY/4oDg9CIkwbXAMRGHFh9HIJBrOeu37dxLWtyal+I+/KO/5abXqvS2YMRD6/s8P7t8SSLpQOH3H7OxWq1t7dNuof7AiAic4A5AJmZmR6EdarGRkNlwiDi45zQK8b11xZufexzRljv8pEJ1jt/VA+IS7U+Lsb3gbje2j8YRESEPj2i6NMjirNo/d6FuoZGyqrqOFFdT0VtPZW1DVTU1FNd10hNfQPVdQ3U1jdSU99IXYOhvqGRukbr3wZjaGgw1DcajDE0GmgwBmOsN5qGRoOBkz9bz83JP2hjmicKThlG2mqHahe+UjNd2blbJVBOFuXA3pZWG0O4qSGm8QTRjZVENlYS1VhFRGMVEaaGiMYawk0tYaaWMFOH09TjpN761zQgNOIwVnoXYxAacNCIuH7/Dhr5Ol270r4xuP/y3bOFnPJ77fjv2H3/+vCEDu/vCU8SfUsZsPnZtLaNJ/taC42ZB8wDq4/eg7hO4XAI436ss9grz4U7HSTFRZIUF2l3KEr5lCeJPhdwn+4oA8jzcJsID/ZVSinlQ57cGbsGyBaRLBGJAGYDC5ttsxD4tlgmAaXGmHwP91VKKeVD7bbojTH1InI3sBRriOSLxpgtInKna/0zwCKsoZW7sYZX3trWvj45E6WUUi0KqhumlFIqVLU1jj6oipoppZT6Jk30SikV5DTRK6VUkNNEr5RSQc4vv4wVkULgQCd3TwaOeTEcOwXLuQTLeYCeiz8KlvOArp1Lf2NMSksr/DLRd4WI5LT2zXOgCZZzCZbzAD0XfxQs5wG+OxftulFKqSCniV4ppYJcMCb6eXYH4EXBci7Bch6g5+KPguU8wEfnEnR99EoppU4VjC16pZRSbjTRK6VUkAvIRC8iM0Vkh4jsFpEHW1gvIvKUa/1GERlrR5ye8OBczhORUhFZ73o8Ykec7RGRF0WkQEQ2t7I+kK5Je+cSKNekn4h8LCLbRGSLiNzbwjYBcV08PJdAuS5RIvKliGxwncsvW9jGu9fFGBNQD6xyx3uAgVgTm2wARjTb5hJgMdYMV5OA1XbH3YVzOQ/4t92xenAu5wBjgc2trA+Ia+LhuQTKNekLjHU9jwd2BvDfiifnEijXRYA41/NwYDUwyZfXJRBb9CcnKzfG1AJNE467OzlZuTHmC6BpsnJ/48m5BARjzArgeBubBMo18eRcAoIxJt8Ys871/ASwDWseZ3cBcV08PJeA4Ppdl7t+DHc9mo+K8ep1CcRE39pE5B3dxh94GudZro95i0VkZPeE5nWBck08FVDXREQGAGdgtR7dBdx1aeNcIECui4g4RWQ9UAD8xxjj0+viyZyx/qYrk5X7G0/iXIdVw6JcRC4B3gGyfR6Z9wXKNfFEQF0TEYkD3gTuM8aUNV/dwi5+e13aOZeAuS7GmAZgjIj0BN4WkdOMMe7fCXn1ugRii74rk5X7m3bjNMaUNX3MM8YsAsJFJLn7QvSaQLkm7QqkayIi4ViJ8RVjzFstbBIw16W9cwmk69LEGFMCfALMbLbKq9clEBN9VyYr9zftnouI9BERcT2fgHXNiro90q4LlGvSrkC5Jq4YXwC2GWOebGWzgLgunpxLAF2XFFdLHhGJBi4AtjfbzKvXJeC6bkwXJiv3Nx6ey7XAD0SkHqgCZhvX1/L+RERexRr1kCwiucCjWF8yBdQ1AY/OJSCuCTAFuAXY5OoPBvgZkAkBd108OZdAuS59gb+LiBPrzeh1Y8y/fZnDtASCUkoFuUDsulFKKdUBmuiVUirIaaJXSqkgp4leKaWCnCZ6pZQKcprolXIRkf8Wkfs7s15EVrk9f9xVlfBxEfmuiKT5Il6lPBVw4+iV6gjXDTRijGn05XGMMZPdfrwDSDHG1IjIJ8Bm/PRuUxUatEWvgo6IDHDVLf8rVv2TfiLygIiscdX2/qXbtg+LNR/AB8BQt+X3iMhW1/YL3F5+hIh8IiJ7ReQet+3LXf8uBGKB1SJyAzAeeEWs+ujRbtunydd109eLSIOI9PfV70SFNm3Rq2A1FLjVGPNDEZmBVdxqAlaxqIUicg5QgVV24gysv4V1wFrX/g8CWa5WeU+31x0GTMOqib5DRJ42xtQ1rTTGXCEi5caYMQAi8gPgfmNMjntwxpg8oGmbu4BzjTEHvPsrUMqiiV4FqwOuOt4AM1yPr1w/x2El/njgbWNMJZxsjTfZiNUSfwerCmKT940xNUCNiBQAvbEKUHWKiEwBbgemdvY1lGqPdt2oYFXh9lyA3xljxrgeg40xL7jWtVYD5FJgLjAOWCsiTY2iGrdtGuhCY0msiSReAG5wm4hCKa/TRK9CwVLge65a5ohIuoikAiuAq0QkWkTigctd6x1AP2PMx8BPgJ5YnwI64wTWJ4dTuEruvg781Bizs5OvrZRHtOtGBT1jzDIRGQ587qpiWw7cbIxZJyKvAeuBA8Cnrl2cwMsi0gPr08AfjTElrn07aj7wjIhUAWcZY6pcyycDZwK/dPty+BJX371SXqXVK5VSKshp141SSgU5TfRKKRXkNNErpVSQ00SvlFJBThO9UkoFOU30SikV5DTRK6VUkPv/51iWHrv68igAAAAASUVORK5CYII=\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], - "source": [ - "z=np.linspace(0,3,100)\n", - "plot(z, nz1(z), label='nz1')\n", - "plot(z, nz2(z), label='nz2')\n", - "legend()\n", - "xlabel('redshift z')" - ] - }, - { - "cell_type": "code", - "execution_count": 5, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "DeviceArray(1., dtype=float32)" - ] - }, - "execution_count": 5, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "# Just for fun, let's check they are correctly normalized\n", - "from scipy.integrate import romberg\n", - "romberg(nz2,0,3)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## Now, build some lensing cls" - ] - }, - { - "cell_type": "code", - "execution_count": 6, - "metadata": {}, - "outputs": [], - "source": [ - "from jax_cosmo.core import Cosmology\n", - "from jax_cosmo.tracers import get_lensing_tracer_fn\n", - "from jax_cosmo.angular_cl import angular_cl" - ] - }, - { - "cell_type": "code", - "execution_count": 7, - "metadata": {}, - "outputs": [], - "source": [ - "# get a cosmology\n", - "cosmo = Cosmology(Omega_c=0.3, Omega_b=0.05, h=0.7, \n", - " sigma8 = 0.8, n_s=0.96, Omega_k=0.,\n", - " w0=-1., wa=0.)" - ] - }, - { - "cell_type": "code", - "execution_count": 8, - "metadata": {}, - "outputs": [], - "source": [ - "# Let's build a function to get some lensing cls\n", - "@jit\n", - "def lensing_cls(cosmo, ell, nz1, nz2):\n", - " \n", - " # We build a tracer function for each nz\n", - " tracer_fn1 = get_lensing_tracer_fn(nz1) \n", - " tracer_fn2 = get_lensing_tracer_fn(nz2)\n", - "\n", - " return angular_cl(cosmo, ell, tracer_fn1, tracer_fn2)" - ] - }, - { - "cell_type": "code", - "execution_count": 9, - "metadata": { - "scrolled": false - }, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in asarray is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n", - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in array is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n", - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype requested in astype is not available, and will be truncated to dtype int32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n" - ] - } - ], - "source": [ - "ell = np.logspace(1, 3, 100)\n", - "cl = lensing_cls(cosmo, ell, nz1, nz2);" - ] - }, - { - "cell_type": "code", - "execution_count": 10, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "[]" - ] - }, - "execution_count": 10, - "metadata": {}, - "output_type": "execute_result" - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAD8CAYAAACRkhiPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3dd3gVdd7+8ffnnFQgBQg1dAKhd+kidlCRFREp6ipNsa/6PKu/3Wf7rm3XgkpXUVQUsaPio6uIVAEBQWroASF0Ekogyff3B7jLwwIGUubMOffrunJdZs45M3dgzM3Md+Y75pxDREQiU8DrACIi4h2VgIhIBFMJiIhEMJWAiEgEUwmIiEQwlYCISASL8jrAuUhJSXF16tTxOoaIiK8sWrRol3Ou0ule81UJ1KlTh4ULF3odQ0TEV8xs05le0+kgEZEIphIQEYlgKgERkQimEhARiWAqARGRCKYSEBGJYL66RPR8HczN45u1O6mSGEeVxDgqJcQSHVT/iYj4ogTMrBfQKy0t7bw+v2HXQe547buT1geVE2KplhRPanI8NcrHU6NCGWqUj6d2hTLUKF+GmCiVhIiEP/PTQ2XatWvnzudmsSPH8snIyiEr+wjb9+ey/cARftx3mB/3H2HrvsNs3XuYo/kF/3p/MGBUT46jbko56qWUpV6lsqRVKkdalXJUKheLmRXnjyUiUqLMbJFzrt3pXvPFkUBRxUUHaZaaBCSd9vWCAkdWdi5b9h5i0+5DbNp9kI27D7FhVw6LNu7h4NH8f703uUw0DSqXI71qAulVE2lSLYFGVRMpGxsRf5QiEmb0mwsIBIyqSXFUTYrjgjoV/s9rzjl2HMhl3c4c1u7IZk1WDmu2Z/PB4m1k524Gjp9eqluxLE2qJ9KiRhLNU5NplppIQly0Fz+OiEihqQR+htm/C6JLWsq/ljvn2LrvMKt+zOaHbQf4Ydt+Fm/ex7TvfzzxOahfqRytaibTqmYybWuXp2GVBIIBnUoSkdChEjhPZkaN8scHkS9rUuVfy3fn5LJs636+z9zPki37+HJVFlMXZQJQLjaK1rWSaVe7Au3rVqB1rWTiooNe/QgiIpExMOwl5xyb9xziu817+W7TPhZs3MPqHdk4B9FBo1XNZDrWq0jHehVpW7u8SkFEit3ZBoZVAh7Yf/gYizbtYf6GPcxbv4dlmfsocBATFaBd7fJ0SUvhwgYpNKueRECnj0SkiFQCIS77yDEWbNzDnIzdzMrYxart2QCULxNN1waV6N6wEt0aVqJSQqzHSUXEjyL+EtFQlxAXzSWNqnBJo+NjC7tycpm1dhcz1+5k5ppdfLR0GwAtaiRxcXplLmtchWapibpfQUSKTEcCIa6gwLHixwN8tSqLr1ZnsXjLPpyDKomxXNq4Clc2rUqnehV1h7OInJFOB4WRXTm5zFi9ky9W7GDm2p0cOppPQmwUFzeqTM9mVbkovRJlYnSAJyL/phIIU0eO5TNn3S4+W76Dz1fuYM/Bo8RFB7g4vTJXt6jGJY0qqxBERCUQCfLyC/h24x6mL9/Op8u3szM7l7joAJc2qkKvltXpnl5Jl5+KRCiVQITJL3As2LiHad9v49Nl29l98CgJsVH0aFaV61qn0qFeRd25LBJBVAIRLC+/gDnrdvPh0m1MX76dnNw8qibG0btVda5vW4OGVRK8jigiJUwlIMDxMYQvVu7gve+28vWaneQVOJqnJnF9m1R6t0qlfNkYryOKSAlQCch/2JWTywdLtvHud5n8sO0AMcEAlzetQr92NemalqLTRSJhJGRLwMxqAc8Du4A1zrnHzvZ+lUDJWLHtAG8v2sL7i7ey99AxUpPj6deuJv0uqEG1pHiv44lIEZVICZjZS8A1QJZzrtlJy3sAzwJBYMLZfrGb2WVAfefcWDN71Tl3y9m2qRIoWbl5+XyxIos3F2zmm7W7CBhc0qgyAzvU4qKGlXV0IOJTJVUC3YAc4NWfSsDMgsAa4HIgE1gADOB4ITx6yioGA/nAVMABk5xzL59tmyqB0rNlzyHeXLCZtxZksisnl9TkeAZ2qEW/djU1h5GIz5TY6SAzqwNMO6kEOgF/cM5deeL7RwCcc6cWwE+ffwj41jk308ymOuf6nm17KoHSdyy/gM9X7OC1eZuYs2430UHjqubVuKVTbdrUKq/5i0R8oDQnkEsFtpz0fSbQ4Szvnw78wcwGAhtP9wYzGw4MB6hVq1bxpJRCiw4GuKp5Na5qXo2MrBxen7+JqQsz+WDJNpqlJnJr57pc06KabkQT8aniPhK4AbjSOTf0xPc3A+2dc/cUPaqOBELFwdw83l28lVfmbCQjK4eUcjEM6lCbQR1rUTkhzut4InKK0jwSyARqnvR9DWBbMW9DPFY2NoqbO9bmpg61mJ2xm5dnb+DZf65l9Ix19GpZnaEX1qVxtUSvY4pIIRR3CSwAGphZXWAr0B8YWNSVmlkvoFdaWlpRVyXFyMzo2iCFrg1S2LDrIBNnb2DKwkze+S6TrmkpDOtWj24NUjRuIBLCinJ10GSgO5AC7AB+75x70cyuAp7h+BVBLznn/lpMWXU6yAf2HTrKG99u5pU5G9lxIJdGVRMY3q0evVpWJzqoZx6IeCFkbxY7VyoB/ziaV8AHS7Yy/pv1rNmRQ/WkOIZeWI/+7WtqemuRUqYSEM845/hqdRZjZqzn2417KF8mmls71+WXnWuTXEZzFYmUBt+XwEljAsPWrl3rdRw5Tws37mH0jHX8c1UWZWOC3NSpNkO61tUVRSIlzPcl8BMdCYSHlT8eYNSMdXz8/TaigwEGtK/FHRfVp2qSykCkJKgEJCRt2HWQUV9l8O7irQTNuPGCmozoXp/qyZq0TqQ4qQQkpG3Zc4hRMzKYuigTQ2UgUtx8XwIaE4gMmXsPMWrGOt5euAXDGNC+JnddnEblRJ0mEikK35fAT3QkEBky9x7i+S8zeHtRJlEB45ZOtRnRPY0KevKZyHlRCYgvbdp9kGf/uZb3F28lPjrIkAvrMfTCuiTGRXsdTcRXVALiaxlZ2Tz1+Ro+Wbad5DLR3NU9jZs71dbMpSKF5PsS0JiAACzfup8nPlvNzDU7qZYUx68ua8j1bWvoiWciP8P3JfATHQkIwJx1u3h8+mqWbtlHwyrleKRnY7qnV9JEdSJncLYS0Ixe4jud66fw/p2dGTWoDUfzCrht4gIGjp/Pssz9XkcT8R2VgPiS2fHHXH7+wEX8qXdTVu/Iptfzs3jgrSVs23fY63givqHTQRIWDhw5xugZ63hx1gYMuL1bPW6/qD5lYzVjqYhOB0nYS4yL5tc9GvHlgxdxZdOqjPwyg4v/PoOpizIpKPDPP3RESpsvSsDMepnZuP37dc5Xzq5G+TKMHNCad0Z0plpyPA+9vZTrRs1m8ea9XkcTCUk6HSRhq6DA8f6SrTz26SqysnPp0yaVh3s20tTVEnF0OkgiUiBg9GlTgy8f6s6I7vWZtvRHLvn710z4Zj3H8gu8jicSElQCEvbKxUbx6x6N+OxX3WhXpzx/+XglVz37DXPX7fY6mojnVAISMeqmlOXlWy9g/C3tOHwsnwHj53H/m4vJyj7idTQRz6gEJKKYGZc3qcLnv7qIey5J45Nl27n071/z6tyN5OsqIolAKgGJSPExQR68Ip3p919Iy5rJ/O6DH+gzajbLt+oKNIksvigBXSIqJaVepXJMGtKeZ/u3Yuu+w1z7/Cz+Mm0Fh47meR1NpFToElGRE/YfOsbjn63ijfmbSU2O56/XNaN7emWvY4kUmS4RFSmEpDLR/O265rx9RyfiogPc+vIC7ntzMbtzcr2OJlJiVAIip7igTgU+ue9C7ru0AZ8s+5HLn57Jh0u34aejZpHCUgmInEZsVJBfXd6QafdcSM0KZbh38mKGvbqQHQd0OamEF5WAyFmkV03g3RGd+c1Vjflm7S4uf+prpi7K1FGBhA2VgMjPCAaMYd3qMf3+bqRXTeCht5cyeOICHRVIWFAJiBRS3ZSyvDW8E7/v1YS563dz+VNf895iHRWIv6kERM5BIGDc1qUun97XjQZVEvjVW0u5fdIidukKIvEpX5SAbhaTUFM3pSxTbu/Eb65qzIw1O7ny6Zl89sN2r2OJnDNflIBz7iPn3PCkpCSvo4j8y09jBdPu6UrVpDhun7SIB6cs5cCRY15HEyk0X5SASChrWCWB9+7swj2XpPHe4kx6PvMN327Y43UskUJRCYgUg5ioAA9ekc7bd3QmGDBuHDeXx6ev4mieHl4joU0lIFKM2tYuzyf3XUi/tjUZPWMd14+ew/qdOV7HEjkjlYBIMSsXG8XjfVsw5qY2bNl7iKtHzuLNbzfrUlIJSSoBkRLSo1k1pt/XjTa1k3n43WXc+fp37D+kQWMJLSoBkRJUNSmOSYM78EjPRny+Ygc9n53Jgo0aNJbQoRIQKWGBgHH7RfV5Z0RnoqMC3Dh2Ls9+sVaPs5SQoBIQKSUtaybz8b0X0rtVKk9/sYZBE+Zp/iHxnEpApBSVi43iqX4t+fsNLVm6ZT89n/2Gr1ZneR1LIphKQKSUmRl929Zg2r1dqZwQy20vL+DRT1dyLF/3FEjp80UJaO4gCUf1K5Xj/bu6MLBDLcZ+vZ7+4+axdd9hr2NJhPFFCWjuIAlXcdFB/nZdc0YOaM3q7dlcPfIbvlql00NSenxRAiLh7tqW1fnw7i5UTYzjtokLeHz6KvJ0ekhKgUpAJETUO3F6aED741NODJown6xsXT0kJUslIBJC4qKDPNqnBU/1a8nSzH1cPXIW89bv9jqWhDGVgEgI6tOmBu/f1YWE2CgGjp/HmK/Xae4hKREqAZEQ1ahqIh/c3YWezarx2KeruOO1RXpgjRQ7lYBICEuIi+b5ga357dWN+WJlFr2fn83q7dlex5IwohIQCXFmxtAL6zF5WEdycvP4xQuz+XDpNq9jSZhQCYj4RPu6Ffj4nq40S03k3smL+dNHK3SXsRSZSkDERyonxvHGsI7c2rkOL83ewKAJ89mZnet1LPExlYCIz0QHA/zh2qY8c2Mrvs/cR6/nZrFkyz6vY4lPqQREfOoXrVN5Z0RnooJGvzFzmbJgi9eRxIdUAiI+1rR6Eh/d3ZX2dSvw3+98z+8+WK5xAjknKgERnytfNoaJt13A8G71eHXuJgZNmM+uHI0TSOGoBETCQFQwwP+7qjHP9m/F0i37uPa5WSzfqqnX5eepBETCSO9Wx8cJAPqOmcNHup9AfoZKQCTMNEtN4oO7u9I8NYl7Ji/miemrKNBD7eUMPC0BM2tiZlPMbLSZ9fUyi0g4qZQQy+tDOzKgfS1GzVjH8EkLyda8Q3Ia510CZvaSmWWZ2fJTlvcws9VmlmFmD//ManoCzznnRgC3nG8WEflPMVEB/nZdM/7Uuylfrd7J9aPnsHn3Ia9jSYgpypHARKDHyQvMLAi8wPFf7k2AASf+td/czKad8lUZmAT0N7MngYpFyCIip2Fm3NKpDq8Obs+OA7lc+8Is5q7T8wnk3867BJxzM4E9pyxuD2Q459Y7544CbwK9nXPLnHPXnPKVdeLrLuBhYNd5/xQiclZd0lL44K4uVCwbw80vzmfyt5u9jiQhorjHBFKBk29bzDyx7LTMrI6ZjQNeBZ48w3uGm9lCM1u4c+fOYg0rEknqpJTlvbu60DkthUfeXcYfP/qBfA0YR7ziLgE7zbIz7mXOuY3OueHOuUHOuVlneM8451w751y7SpUqFVtQkUiUGBfNS79sx21d6vDy7I0MfWWBBowjXHGXQCZQ86TvawC6UFkkhEQFA/y+V1P+el0zZq7dRd/Rc9myRwPGkaq4S2AB0MDM6ppZDNAf+LCoKzWzXmY2bv9+3QEpUlwGdajNK7e158f9h7lu1GwWbdrrdSTxQFEuEZ0MzAXSzSzTzIY45/KAu4HPgJXAFOfcD0UN6Zz7yDk3PCkpqairEpGTdG2Qwnt3daFsbBQDxs/THcYRyJzzz8BQu3bt3MKFC72OIRJ29hw8yu2TFrJg414evLwhd1+ShtnphvjEj8xskXOu3ele07QRIkKFsjG8NrQD17VO5R+fr+Ght7/naJ6mpI4EUV4HKAwz6wX0SktL8zqKSNiKjQryVL+W1KlYlqe/WMO2fYcZc1NbkspEex1NSpAvjgQ0JiBSOsyM+y5rwNM3tmTRpr30GT1bVw6FOV+UgIiUruta12DSkPbsyjnKdaNm6xnGYUwlICKn1aFeRd69szPxMUH6j5vL//6w3etIUgJ8UQK6T0DEG/UrleO9O7uQXjWR219bxMTZG7yOJMXMFyWgMQER76SUi+XNYR25vHEV/vDRCv4ybYUeUhNGfFECIuKt+Jggo29qy62d6zBh1gbumbyYI8fyvY4lxcAXl4iKiPeCAeP3vZqQmhzPXz9Zyc7sXMbd0pbkMjFeR5Mi8MWRgMYEREKDmTGsWz2eG9CaJVv20XfMXDL36hJSP/NFCWhMQCS09GpZnVcGt2fHgSP0GTWHFdsOeB1JzpMvSkBEQk+n+hV5Z0RnggGj39i5zMnQwwH9SCUgIuetYZUE3r2zM6nJ8fzy5W/5ULOQ+o5KQESKpFpSPFPu6ETrWuW5d/JiJnyz3utIcg5UAiJSZEnx0bw6uD09m1XlLx+v5NFPVupeAp/wRQno6iCR0BcXHeT5gW24uWNtxs5cz0NvL+VYvqajDnW+KAFdHSTiD8GA8afeTXnw8oa8u3grQ19ZyKGjeV7HkrPwRQmIiH+YGfdc2oDH+jTnm7U7GTh+PnsPHvU6lpyBSkBESkT/9rUYfVNbVvx4gL5j5rB132GvI8lpqAREpMRc2bQqkwa3Jys7l+tHzWHtjmyvI8kpVAIiUqI61KvIW8M7ke8cN4ydy+LNe72OJCdRCYhIiWtSPZF37uhMYlw0gybMZ+aanV5HkhN8UQK6RFTE/2pVLMPUEZ2oXbEsQ15ZwLTvdXdxKPBFCegSUZHwUDkhjjeHd6RVzWTumbyY1+dv8jpSxPNFCYhI+Dh+d3EHLk6vzG/eW84LX2XgnO4u9opKQERKXXxMkLE3t+UXrarz5Ger+dsnK1UEHtGTxUTEE9HBAE/1a0VifDTjv9nAgcN5/K1Pc4IB8zpaRFEJiIhnAgHjj9c2JTk+mpFfZpCde4ynb2xFbFTQ62gRQyUgIp4yMx64Ip3E+Gj+8vFKcnIXMeamNpSJ0a+n0qAxAREJCUMvrMcT17dg1tqd3PLitxw4cszrSBFBJSAiIaPfBTV5bkAblmbuY8C4eezOyfU6UtjzRQnoZjGRyHF1i2qMv6UdGVk59Bs7l+37j3gdKaz5ogR0s5hIZOmeXplXB7dnx4Fcbhg7h827D3kdKWz5ogREJPJ0qFeR14d2IPtIHjeM1QykJUUlICIhq2XNZN4a3okCBzeOm8cP23RKuLipBEQkpKVXTWDK7Z2IiwowYNw8vtNU1MVKJSAiIa9uSlneHtGZCmVjuGnCfOas2+V1pLChEhARX0hNjmfK7Z2oUT6e215ewIzVWV5HCgsqARHxjcqJcbw5vBNplcsx7NWFTF++3etIvqcSEBFfqVA2hjeGdaRZahJ3vfEdHyzZ6nUkX1MJiIjvJMVHM2lIB9rVLs/9by3h7YVbvI7kWyoBEfGlcrFRTLytPV3TUvivqd/z2jw9pex8qARExLfiY4KMv6UdlzaqzG/fX85LszZ4Hcl3fFECmjtIRM4kLjrI6Jva0qNpVf40bQVjvl7ndSRf8UUJaO4gETmbmKgAzw9sTa+W1Xns01WM/OdaryP5hp7aICJhISoY4JkbWxETDPDU52s4mlfAg1c0xEyPqzwblYCIhI1gwHiybwtiooznv8rgWEEBD/dopCI4C5WAiISVQMD46y+OP7B+7NfrOZbn+J9rGqsIzkAlICJhJxAw/ty7GdHBAC/N3kBeQQF/vLapiuA0VAIiEpbMjN9d04ToYIBxM9eTX+D4c+9mBAIqgpOpBEQkbJkZj/RsRDBgjJ6xjvwCx9+ua64iOIlKQETCmpnx31emExUwnvsyg/wCx2PXtyCoIgBUAiISAcyMB69IJxgwnvliLQUOnuirIgCVgIhEkPsva4hhPP3FGpxzPHlDy4gvApWAiESU+y5rQDAAf//fNRQ4xz/6tYroIlAJiEjEufuSBpgZT362Ggf844aWRAV9MYtOsVMJiEhEuuviNMzgiemrcQ6e6heZRaASEJGIdWf3NAzj8emrcMDTEVgEKgERiWgjutfHDB77dBVG5B0RqAREJOLdcVF9nIPHp68CIqsIVAIiIhw/IoDIK4JSKwEzqwf8BkhyzvU9sawsMAo4Csxwzr1eWnlERE51chEEjIi4fLRQNWdmL5lZlpktP2V5DzNbbWYZZvbw2dbhnFvvnBtyyuI+wFTn3DDg2nNKLiJSAkZ0r89/XZnO+0u28V9vLyW/wHkdqUQV9khgIvA88OpPC8wsCLwAXA5kAgvM7EMgCDx6yucHO+eyTrPeGsCyE/+dX/jYIiIl566L0wB48rPVYPBk3/C9s7hQJeCcm2lmdU5Z3B7IcM6tBzCzN4HezrlHgWsKuf1MjhfBEnzyvGMRiQx3XZxGQYHjH5+vIWDGE9e3CMvZR4vyizcV2HLS95knlp2WmVU0szFAazN75MTid4HrzWw08NEZPjfczBaa2cKdO3cWIa6IyLm559IG3HdpA6YuyuSRd5dREIanhooyMHy6Sjzjn5BzbjdwxynLDgK3nW0jzrlxwDiAdu3ahd/fgIiEtPsva0CBczz3ZcaJR1eG14NpilICmUDNk76vAWwrWhwRkdBiZjxweUMKnOOFr9YRDMCfezcLm0dVFqUEFgANzKwusBXoDwwsllSnMLNeQK+0tLSSWL2IyFmZGQ9dkU5egWPs1+uJCgT4fa8mYVEEhb1EdDIwF0g3s0wzG+KcywPuBj4DVgJTnHM/lERI59xHzrnhSUlJJbF6EZGfZWY83KMRQ7rWZeKcjfz145U45/8z1IW9OmjAGZZ/AnxSrIlEREKUmfHbqxuTX+CYMGsDUcEAv+6R7usjAk0bISJyDsyM3/dqwrH8AsZ8vY6YoPHAFelexzpvvigBjQmISCgxM/7cuxnH8gsY+WUGUcEA917awOtY58UXN2hpTEBEQk0gYDzapwV9Wqfy1OdrGPP1Oq8jnRdfHAmIiISiYMB4om8LjuYX8Ninq4gOBhjSta7Xsc6JSkBEpAiiggGevrEVefmOP09bQUxUgJs71vY6VqH54nSQmfUys3H79+/3OoqIyH+IDgYYOaA1lzWuzP+8v5wpC7b8/IdChC9KQGMCIhLqYqICvDCoDd0aVuLX737P+4u3eh2pUHxRAiIifhAbFWTsTW3pWLciD0xZwifLfvQ60s9SCYiIFKP4mCATftmONrXKc+/kxXyxYofXkc5KJSAiUszKxkbx8m0X0LR6Ine+/h0z14TuNPi+KAENDIuI3yTERfPK4PbUr1yO4ZMWMm/9bq8jnZYvSkADwyLiR8llYnhtSHtqli/D4IkLWLRpr9eR/oMvSkBExK8qlovl9aEdqJwQy60vf8vyraF1RkMlICJSwionxvH6sI4kxkVz84vzWbMj2+tI/6ISEBEpBanJ8bwxrAPRwQADx89nw66DXkcCfFICGhgWkXBQu2JZXh/agQLnGDR+Hpl7D3kdyR8loIFhEQkXDaokMGlIe3Jy8xg0YT47DhzxNI8vSkBEJJw0rZ7ExMHt2ZWdy6AJ89mdk+tZFpWAiIgH2tQqz4RfXsCWPYe45aVv2X/4mCc5VAIiIh7pVL8iY25uy5od2QyeuIBDR/NKPYNKQETEQxenV2Zk/9Ys3ryXYa8u5Mix/FLdvkpARMRjPZtX4+83tGR2xm7ufuM7juUXlNq2fVECukRURMJdnzY1+MsvmvHFyiwemLKU/AJXKtv1RQnoElERiQQ3dazNwz0b8dHSbfzmvWU4V/JFoGcMi4iEkDsuqk/OkTye/yqDsrFR/PbqxphZiW1PJSAiEmIevKIhObl5vDhrAwlxUdx/WcMS25ZKQEQkxJgZv7umCTm5eTzzxVoS4qIZ0rVuiWxLJSAiEoICAeOxPs05mJvHn6etoFxskBsvqFX82yn2NYqISLGICgZ4pn8rujWsxJItJXN1pI4ERERCWGxUkHE3tyU2qmT+za4SEBEJcXHRwRJbty9OB+lmMRGRkuGLEtDNYiIiJcMXJSAiIiVDJSAiEsFUAiIiEUwlICISwVQCIiIRzEpjqtLiYmY7gX3A2a4VTTrL6ynAruLOVcLO9vOE8raKsq5z/Wxh31+Y9/3ce8Jt/4LS28e0f3m3f9V2zlU67SvOOV99AePO93Vgodf5i/vnDdVtFWVd5/rZwr6/MO+LtP2ruP/eS2s72r+K78uPp4M+KuLrflOaP09xbqso6zrXzxb2/YV5X6TtX1B6P5P2rxDcv3x1OqiozGyhc66d1zkkPGn/kpJUUvuXH48EimKc1wEkrGn/kpJUIvtXRB0JiIjI/xVpRwIiInISlYCISARTCYiIRLCILQEzq2dmL5rZVK+zSPgxs1+Y2Xgz+8DMrvA6j4QfM2tsZmPMbKqZjTjf9YRVCZjZS2aWZWbLT1new8xWm1mGmT0M4Jxb75wb4k1S8aNz3L/ed84NA24FbvQgrvjQOe5jK51zdwD9gPO+dDSsSgCYCPQ4eYGZBYEXgJ5AE2CAmTUp/WgSBiZy7vvXb0+8LlIYEzmHfczMrgVmAf883w2GVQk452YCe05Z3B7IOPEv/6PAm0DvUg8nvncu+5cd9zjwqXPuu9LOKv50rr/DnHMfOuc6A4POd5thVQJnkApsOen7TCDVzCqa2RigtZk94k00CQOn3b+Ae4DLgL5mdocXwSRsnOl3WHczG2lmY4FPznflUUVN5wN2mmXOObcb0P+cUlRn2r9GAiNLO4yEpTPtYzOAGUVdeSQcCWQCNU/6vgawzaMsEn60f0lJK9F9LBJKYAHQwMzqmlkM0B/40ONMEj60f0lJK9F9LKxKwMwmA3OBdDPLNLMhzrk84G7gM2AlMMU594OXOcWftH9JSfNiH9MEciIiESysjgREROTcqARERCKYSkBEJIKpBPEeB68AAAAlSURBVEREIphKQEQkgqkEREQimEpARCSCqQRERCKYSkBEJIL9f+6s7gXreH/kAAAAAElFTkSuQmCC\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], - "source": [ - "loglog(ell, cl)" - ] - }, - { - "cell_type": "code", - "execution_count": 11, - "metadata": {}, - "outputs": [], - "source": [ - "# Let's compare to CCL\n", - "import pyccl as ccl\n", - "cosmo_ccl = ccl.Cosmology(\n", - " Omega_c=0.3, Omega_b=0.05, h=0.7, sigma8 = 0.8, n_s=0.96, Neff=0,\n", - " transfer_function='eisenstein_hu', matter_power_spectrum='linear')\n", - "\n", - "tracer1 = ccl.WeakLensingTracer(cosmo_ccl, (linspace(0,6), nz1(linspace(0,6))), use_A_ia=False)\n", - "tracer2 = ccl.WeakLensingTracer(cosmo_ccl, (linspace(0,6), nz2(linspace(0,6))), use_A_ia=False)" - ] - }, - { - "cell_type": "code", - "execution_count": 18, - "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in asarray is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n", - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in array is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n" - ] - }, - { - "data": { - "text/plain": [ - "Text(0, 0.5, 'Lensing angular $C_\\\\ell$')" - ] - }, - "execution_count": 18, - "metadata": {}, - "output_type": "execute_result" - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZMAAAEMCAYAAAABLFv3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3dd3xUVf7/8ddnkkDoNSC9V0VBiqIgxVWRjgUQEAu7CBaKFMu6Lm6xIi4oXxEUUEREUZqALB1UBERpglQpAUGkhBpIOb8/iP5YDJhhZnInmffz8ZgHzsnMnTdZZz+ee+79HHPOISIiEgif1wFERCTrUzEREZGAqZiIiEjAVExERCRgKiYiIhIwFRMREQlYtNcBMpOZtQHa5MuX7y9Vq1b1Oo6ISJayevXqX5xzcen9zCLxPpN69eq5b775xusYIiJZipmtds7VS+9nOs0lIiIBUzEREZGAqZiIiEjAImoBXkSyj6SkJOLj40lMTPQ6SrYTGxtL6dKliYmJyfB7VExEJEuKj48nX758lC9fHjPzOk624Zzj0KFDxMfHU6FChQy/T6e5RCRLSkxMpEiRIiokQWZmFClSxO8ZX0QVEzNrY2ajExISLvsYqWdOwtlTQUwlIpdLhSQ0Luf3GlHFxDk30znXs0CBApf1/p+PJfLmK09y+pUaJM57Hk4dDnJCEclK9u/fT+fOnalUqRI1a9akZcuWbNmyhS1bttCyZUsqV65MjRo16NixIwcOHGDx4sW0bt3a69ghEVHFJFAnziTzS5G6fJlYgdgvX+Ls0BqcmDYAju72OpqIZDLnHB06dKBp06Zs376djRs38vzzz3PgwAFatWpF79692bZtG5s2baJ3794cPHjQ68ghpWLih4pxefl77/sp2XsGL1Z8l5lJDcj53Th+GNODzfuPex1PRDLRokWLiImJoVevXr+N1a5dm61bt9KwYUPatGnz23izZs246qqrvIiZaXQ112WoWTI/Nbu3Z+/R2xi54Gvmr93Ohv8s5a7K8FffuxS6ZRCUqe91TJGI8dzM79m471hQj1mzZH7+3ubKi/58w4YN1K1bN8Pj2Z1mJgEoVTAX/e5sxoQnutP/T1U5vXcDtusLeOdPHPm/23DbF0ME9j4TkcijmUkQFMqTg75/qsKpmwbwyfLWHFk6mk4HpmET2nG4yLUU6P1foqIzfvOPiPjnUjOIULnyyiuZMmVKuuNLlizJ9Dxe08wkiHLniObeJlfS66nXWNZyAcNy9ubDA6W55T9fMmV1PMk7lkFqitcxRSQImjdvzpkzZxgzZsxvY6tWraJy5cp89dVXzJo167fxzz//nPXr13sRM9OomIRAjmgfd11Xmb5PvEC5ji8TGxPF21NmEv1eaxJerUvSdx9CSrLXMUUkAGbG1KlTmTdvHpUqVeLKK69kyJAhlCxZks8++4zXX3+dKlWqULNmTcaPH0+xYsUAWLBgAaVLl/7tsXz5co//JsGh/UwygXOORZt+4tvPx9P66AdU9+0hIVcZct38BDnq3ANROtso4q9NmzZRo0YNr2NkW+n9frWficfMjOY1SzKg/1McunchQwv9jT0nozj72SDeW7yeU2c1SxGRrC2iikkw2qkE+PncWKUYA/sO5OR9CxhSfCTPzttHoxcXsmNkB86sHA8pSZ5kExEJREQVk0DbqQTTdZWKMrT3nXzSuyENS8DxAzvJObsvx165hjMr31VREZEsJaKKSTiqW64wI/9yKyk9FjAs7l/8eConOWf3IeGV2pz+aZPX8UREMkTFJExcW64wjz/yGKk9FvBa0X+w8mRxmo75kdFLt5O4f6suKRaRsKZiEmbqlCtM/0f7UvjPn1C1VBFenb2OhFG3cuTVepxdOwVSU72OKCLyOyomYapuuUJM6HEdEx9qzIeFe/PL8TPkmNqDw69dT9LGWWrTIiJhRcUkzNWrUJS+fQZz8N5FDC8wiGMJR4n5qAuLPv+E5BTNUkS8dsMNN3gdISzobrks4oYqxWnY768s2/wgH86awKglsZTfsIRhVb+ndu16+Mpd73VEkYj01VdfeR0hLKiYZCFmxk3VS9K42mDqbvqZ1+ZupPDqEfi+O8DBEk0p2vafWImrvY4p4o1xrX4/dmV7aPCXc1ttT7z79z+v3QXqdIWTh+Cj7v/7swdm/f716cibNy/79++nXbt2HDlyhKSkJP71r3/Rrl07Vq1aRY8ePVi5ciUpKSk0aNCAyZMnX3Rvk5dffpkJEybg8/m4/fbbefHFF1mzZg29evXi1KlTVKpUibFjx1KoUCFGjBjBqFGjiI6OpmbNmnz44YcMGTKEH3/8kZ9++oktW7YwbNgwvv76a+bMmUOpUqWYOXMmMTExLFiwgIEDB5KcnEz9+vV58803yZkzZ4b+vhej01xZkJlxS83ifNa3KRvbzWZ0jnvJsW8l9lZjDo7vBkd2eR1RJKLExsYydepUvv32WxYtWsSAAQNwzlG/fn3atm3LM888w+DBg+nWrdtFC8mcOXOYNm0aK1asYO3atQwePBiA7t2789JLL7Fu3Tpq1arFc889B8CLL77Id999x7p16xg1atRvx9m+fTuzZs1i+vTpdOvWjWbNmrF+/Xpy5crFrFmzSExM5P7772fy5MmsX7+e5ORk3nzzzYB/B5qZZGE+n9GybmWSa49gxoqHOL7gNe78cRZPT7qdjm0KULtMQa8jimSeS80kcuS+9M/zFMnwTCQ9zjmefvppli5dis/nY+/evRw4cIArrriCZ599lvr16xMbG8uIESMueoz58+fzwAMPkDt3bgAKFy5MQkICR48epUmTJgDcd9993H33uRnW1VdfTdeuXWnfvj3t27f/7Ti33347MTEx1KpVi5SUFFq0aAFArVq12LlzJ5s3b6ZChQpUrVr1t2OOHDmSfv36XfbfHzQzyRaio3zccUMtOj05hk+bzefzQ8VoP/JLlg7rwpHpT8HpI15HFMnWJk6cyMGDB1m9ejVr1qyhePHiJCYmAnD48GFOnDjB8ePHfxtLj3MOM8vwZ86aNYtHHnmE1atXU7duXZKTz/X4+/V0lc/nIyYm5rdj+nw+kpOTCVVzXxWTbCQ2JoruTWuxdHAz+t9cmaMJxynw7ZucGlqLY/NfgaTTXkcUyZYSEhIoVqwYMTExLFq0iF27/v+p5p49e/LPf/6Trl278sQTT1z0GLfeeitjx47l1KlTwLkiVKBAAQoVKsSyZcsAmDBhAk2aNCE1NZU9e/bQrFkzXn75ZY4ePcqJEycylLV69ers3LmTbdu2/c8xA6XTXNlQ3pzR9L2lGodv+IS3Z82h2oZhNPniXxxbMRruGkv+ao29jiiSbZgZXbt2pU2bNtSrV4/atWtTvXp1AN577z2io6Pp0qULKSkp3HDDDSxcuJDmzZv/7jgtWrRgzZo11KtXjxw5ctCyZUuef/553n333d8W4CtWrMi4ceNISUmhW7duJCQk4Jyjf//+FCyYsdPasbGxjBs3jrvvvvu3BfhevXoF/nvQfibZ396jp5k+7SNq7xjFX319uatJPR6sW5hc+QqBH9NqkXASDvuZHDp0iGuvvfZ/ZiLZhfYzkd8pVTAXD99/H0Ue/i+VKlTilbk/sPm1lvw8ojnJu1d6HU8kS9q3bx8NGzZk4MCBXkcJC9niNJeZ1QSGAIeABc65Kd4mCk/VrsjH2/fVZ+WOX/j60+bcefh9osfewoHSLSjW4XmsSCWvI4pkGSVLlmTLli1+v2/9+vXce++9/zOWM2dOVqxYEaxonvD8NJeZjQVaAz875646b7wFMByIAt52zr14iWMMAFY655aZ2QznXNtLfWakneZKj3OO+Wu3s3f2UO4+8yk5LZmdt42jcsN2XkcTyZBNmzZRvXp1v66AkoxxzvHDDz9kudNc44EW5w+YWRQwErgdqAncY2Y1zayWmX12waMYMAHobGavAEUyOX+WZGbcUrsyXZ/4P+Y0m80H1orW01N45INv2bd1ja78krAXGxvLoUOHQnapa6RyznHo0CFiY2P9ep/nMxMAMysPfPbrzMTMGgJDnHO3pT1/CsA598IfHCcK+NQ597v/vDaznkBPgLJly9bNjgtmgThxJpnRS3cwdulW5vj6kS+nj5hb/k6eup3BFw7/zSHyv5KSkoiPj7/kvRtyeWJjYyldujQxMTH/M36pmUm4FpO7gBbOuT+nPb8XuM459+gl3v80kAd40zn3xaU+T6e5Lu7AsUSmT51Mw+2vUcv3Iwfz1aRA+5fJUUmXE4tEunA/zZWe9E6CXrTqOed2Oud6Oue6/lEhkUsrnj+WnvfdR47eS3iryBMkHdtPjgmt+XrBVJ1OEJGLCtdiEg+UOe95aWBfoAc1szZmNjohISHQQ2V71UoU4KHHnmZHp6X8J9ej3DMvmrtGLWfrillqzyIivxOup7migS3AzcBeYBXQxTn3fTA+T6e5/JOS6vj4mz288d/1fJb0F2Kioki6aTAFG/eCqJg/PoCIZAthfZrLzCYBy4FqZhZvZj2cc8nAo8BcYBPwUbAKifgvymd0blCWuYNuY3rt0axNLkvBxc9waGg9Er+foy2ERSQ8ZiaZTTOTwOw7copZU8Zy857Xqejbz9ybPuWWps3x+XS9v0h2FvZXc2UWM2sDtKlcufJftm7d6nWcLO+7Hw8we9pExhyoxlWl8vPaVTupcl0ryF3Y62giEgIqJhfQzCR4nHPMWLuPt2Z/zdQzPUmKyk3STU9SqPFDEJUtuvWISJqwXjORrM3MaFe7FJ8MbM/H177P+pQyFFr8NL8MrU/i5vlexxORTBJRxUSXBodOrhxRdGt3O+X6zeftUv/k5MkT2KTOzFm+RveniEQAneaSkPh2+098PGMakw6U5dqyBRlefRNlbuwMOfN6HU1ELpNOc0mmu7ZSCf7dtxcv33U1OQ9tpMzSASS8UptjqybpUmKRbEjFRELG5zM61ivDW4PuZ3z1t9h9Ni/5Z/Vi//DmJO1d63U8EQmiiComWjPxRv7YGO7v3JlcDy/h7UL9iTmylWNvt+OrzXu9jiYiQaI1E8lUzjkWr93Kx3MWMDuhLK2ujOMfVbZSpEFn8EV5HU9ELuFSaya6EUAylZnRrHZVGl5ZiRpLd7B9yQSKbB/OgS/eoOBd/yFn+eu8jigilyGiTnNJ+IiNieKxm6swqP+TjC3+Vzi+j5zjb2Xvuz3g5C9exxMRP6mYiKdKFcrNg70Hs6PTEibHtKfYjqlsHdGWXYdOeh1NRPwQUWsm6s0V3s4mpzJt3gKmLN/CmtSK9LmxOH+pkaRTXyJhQr25LqAF+PC2PyGR52dvotL3w+kbPZW9Fe6i1F0vQ54iXkcTiWi6aVGylCsKxDLinjpc120Ik3N0oNiOqZx49RoOL30LUlO9jici6VAxkbB1ffXy3PHEWKZd9yEbU8pQeOFgvnv7YRKTUryOJiIX0GkuyRJ+OnqKuZNe593dRaBIZZ6/pRgNq5aEXAW9jiYSMXSaS7K8EgVzc3/vJ3jugXYAHJ3Sh2NDa3N0xUT1+hIJA34VEzNrZGaPmFnF88YqBD9WaKidStZ3U9U4Pu/XmKP1+vJjcmEKznmYvSNuJfnAD15HE4lo/s5M4oAGwN/NbJiZ1QOeD36s0HDOzXTO9SxQoIDXUSQAOaOjuKddGwo9uoR3C/ch7+ENuDdvZPuyyV5HE4lYfhUT59xU4EFgJLABaAxsCUEukT9UNi4f3R/7B6vb/JeZvpu5a7bjqU/Xc/ToYa+jiUScP+zNZWZ/A045514FcM6lACvTHiKeMjOa16vFiVofsHHeFt77ajvd13XnRMmqlOr8Gpa/pNcRRSJCRmYm9wJvXjhoZn82s6eCH0nEf3lzRvNM65pMe/gGvs1zI3F7F3D6tbocXPA6pOpSYpFQy0gxOe2cO5XO+ASgW5DziASkZuki3DPgdeY2mcYaV5m4Zc+wf1gjEg/t9jqaSLaWoWJiZiUuHHTOnQGSgx9JJDA+n9G2eSOqPD6P90r+jR0J0HbsZr7cpm7EIqGSkWLyKjDdzMqdP2hmxQD1tpCwFZc/lu49B5LafSZniabn24v58aUbObZ2htfRRLKdPywmzrmPOXf11moz+8zM/mVmzwNfAkNDHTCYdJ9JZGpUpSif97uJAQ3zk3TqKPmn3sueUXfhErRtsEiwZLidipnlAzoAVwIngdnOuSzZk0TtVCLXtp8OsXLiP7jj+ERSfTGcaPxXijV9GHxqBiHyR9SC/gIqJpEtNdUxZ+mXFF78JC7VsaLRWB5uXpmc0dqDXuRStAe8yHl8PqNV00YcvHY+r8z4ho8WbmPF2nW8WnkdpVr/FWJivY4okuVobi8RKy5/LC93a8T4B+rT4OwKSq0dwS+v1OXEDwu9jiaS5WSomNg5ZUIdRsQLTasVo9egF/mg2uucOJNE3g87sGfcg7hTassiklEZKibu3MLKtBBnEfFM7hzRdLmnOycfXMrHsXdRYudUFox6nJ8STnsdTSRL8Oc019dmVj9kSUTCwJXlrqDDoDFMb/ABfz3SiluGLWX6vEWkHo33OppIWPOnmDQDlpvZdjNbZ2brzWxdqIKJeCU6ysedrW7n436tqFO2ICWXDSZxeD1+XviG9qAXuQh/7jMpl964c25XUBNlAl0aLBnlnGPusq8puHAQ17OeffmvoWiX0eS4orrX0UQyXVC27U0rGseA4kC58x4i2ZaZ0eKmhlQeMJ/3r3iS3AnbYNSNbFn5udfRRMJKhouJmf0ZWArMBZ5L+3NIaGKFhtqpyOUqmi+Wbr2eYn37eXwc1YrWU8/w3MzvOXnimNfRRMKCP2smfYH6wC7nXDOgDnAwJKlCRNv2SqAa17mStoPeofP1lfj4y40ce/Va9kweAGfT26VBJHL4U0wSnXOJAGaW0zn3A1AtNLFEwle+2Bj+0e4q3nuwAauj61Bm09scHFqPEz8s8jqaiGf8KSbxZlaQc/ebzDez6cC+0MQSCX/XVi3HnwZ/yOSaIzl1Jom8H7Zn97s9IfmM19FEMt1lNXo0syZAfuBz51xS0FOFmK7mkmDbuHs/P3zwJEVObufDqq/xXPurKJZPPb4kewmoa7CZHQfSe5Fx7ub4/IFHzFwqJhIKSSmpvL10K68t2EH56CO8VXY+5Tu9guUp4nU0kaAI6NJg51w+51z+dB75smIhEQmVmCgfvZtVY3afxrTM/yOld03j2Kt1OfzNFK+jiYScugaLBFnlYnnp0+9pZjf8gPiUghT+rAe7Rt2NO37A62giIePPHfDPpjfunPtHUBNlAp3mksyy++Axvpzwd+5IeI/FeVtS48G3KFskt9exRC5LUO6A59xWvb8+UoDbgfIBpxPJxsrG5adz/2HMv2kKQ0504Lb/LOWTuQtJTdCFkJK9ZHinRefcq+c/N7OhwIygJxLJZsyMVjc3o07d0zz96ToqfPkwp77ez+nm/yKu0f1g5nVEkYAFsmaSG6gYrCAi2V3JgrkY90ADfr55OFtcaeIW9GP3G61JObrX62giAfOnN9f6tNbz68zse2AzMCJ00USyHzOjRZNGlOq/iEmFHybul5UkDq/Pno3LvY4mEpDLbUGfDBxwziWHJFWIaQFewoFzjvlfLCdhwVD+nvwAfW6pyZ8blScqKsrraCLputQCfIbXTIA70zlwArDaObfmcsOJRCoz45bGN/Bz7Q+4ceoG3pjzLbcv7UDOm/pQ/KYeWkuRLMWfNZN6QC+gVNqjJ9AUGGNmg4MfTSQyFMsXy1v31mVou0r8nJKH4osGsPv1lqRoq2DJQvwpJkWAa51zA5xzAzhXXOKAm4D7Q5AtXWZW0czeMbMp543lMbN3zWyMmXXNrCwiwWJm3NbwWso9vpBJRR+j6KHVJA5vwIEl78Bl9M8TyWz+FJOywNnznicB5Zxzp4EMtUk1s7Fm9rOZbbhgvIWZbTazbWb25KWO4Zzb4ZzrccHwHcAU59xfgLYZySISjuLy56LzI//ky1un84Mry9YF4xi9dDspqSooEt78WTP5APg6rfW8AW2ASWaWB9iYwWOMB94A3vt1wMyigJHALUA8sMrMZgBRwAsXvP9B59zP6Ry3NLA+7Z9TMphFJCyZGbfc2JCfay3kX5+sZMaczaxat4F/1zlKsRvv01qKhCV/blr8p5nNBhqlDT3knPv1kqgMnVpyzi01s/IXDDcAtjnndgCY2YdAO+fcC0DrDMaL51xBWcNFZltm1pNz6zyULVs2g4cV8U6x/LkZfn8Tmq/Zx5HpT1Js/gx2rZlKme5v4ct/hdfxRP6HP/eZ5OTczop5gAJAy4v16/JTKWDPec/j08YulqOImY0C6pjZU2nDnwJ3mtmbwMz03uecG+2cq+ecqxcXFxeE2CKhZ2a0r1OKlv1GMblQL644+CUnX6vHweUTtZYiYcWf01zTgQRgNRlcI8mg9ObsF/2WOOcOce6qsvPHTgIPBDGTSFgpXjAPHfu8yNwlLSm5+HGunvswX+/eS4O7B+Hz6bSXeM+fYlLaOdciBBnigTLnfw4h2g7YzNoAbSpXrhyKw4uElJnRomkTfrpmMRPf/TcvfVeeq46t4JW2FSlVvJjX8STC+XM111dmVisEGVYBVcysgpnlADoTogaSzrmZzrmeBQoUCMXhRTJFiUL56NL3BZ7scD0b9/zC8f/7EzvHdMWdOuJ1NIlg/hSTRsDqtEt41/3aq8ufDzOzScByoJqZxZtZj7SWLI8Cc4FNwEfOue/9Oa5IpDEzulxXlpmPNWZDvkaUjp/N0VfrcmTtbK+jSYS63N5cv3HO7Qpqokyg3lySnaSmOj6bO5saXw+misXzY/mOlO/yHyxHHq+jSTYTlM2xnHO70nsEL2bomVkbMxudkJDgdRSRoPH5jLa3tyK61xKm5r6TkztW0nfyOg6dCOZ1MiKXluGZCYCZFQKqALG/jjnnloYgV0hpZiLZVXJKKmMW/8CwhTspGXuWsdVWUanDMxCTy+tokg0EZWZiZn8GlnJubeO5tD+HBCOgiARHdJSP3jfXZMajjWiVcy2VNr7BgaHXcWLnKq+jSTbnzwJ8X6A+sMs51wyoAxwMSaoQ0WkuiRQ1SuSn3+N/4+Maw0lNPE7s+NvY/emzkJLkdTTJpvwpJonOuUQ4dze8c+4Hzt0Rn2Xo0mCJJDmifdzd6X5+7raIRdGNKLtuOCve6k1iktrXSfD5U0zizawgMA2Yl9bwMSQ3F4pI8FxTpTyNBk9lcoV/02/3TbQasYz123dDaqrX0SQb8WsB/rc3mTXhXH+uz51zZ//o9eFGC/ASqb7Y+guDP/6OlxOHULZQLCXuG0tM4XSv+hf5naAswJ/PObfEOTcjqxUSrZlIpGtUpShz+jVhb8nbKXx0A2dfv54Dy8araaQE7LJmJlmdZiYisHj5SgrMfZQ6bGZn3M2Uvf9tfHkKex1LwljQZyYikvU1bdiAUv0X8XGhP8OBDTw08Vv2HT3tdSzJovxpp/J4OsMJwGrn3JqgpgoxzUxE/j/nHB99vZ3n5mwjpy+FD6ovp3qHJ7Cc+byOJmEmWDOTepzbR6RU2qMn0BQYY2aDAw0pIt4wMzo1rMycvo25o+B2qm16g1+GNuD41i+9jiZZiD/FpAhwrXNugHNuAOeKSxxwE3B/CLIFnRbgRS6uXJE8PN3nMabXfoszZ8+Se2Jrdn38tG50lAzxp5iUBc6/eisJKOecO01wd14MGd20KHJpUT6jQ4dOHLt/CQtimlLu+5GseaMbp8/qRke5NH92WvwA+DrtZkUD2gCTzCwPsDEU4UTEGzUrlKbi4Cl8NOktxm4yzr6+jOF3VqdWueJg2iZYfs/frsF1ObdJFsCXzrksuYqtBXiRjPty2y8M/HgtA04Np37RJErd9w7RBUt6HUs8EKyuwTk514srD+fufm9pZs8GJ6KIhKsbKxfl8z6NSS1Rh2KHv+HU8Ov4ecXHXseSMOPPmsl0oB2QDJw87yEi2VyBPDno+PBzfPWnT9njilJszp/58e37cIm6mEXO8WfNpLRzrkXIkohI2Lu5cWP21VzC1HGDabpnJk998AWDO91C4Tw5vI4mHvNnZvKVmdUKWZJMoEuDRQJXskh+2j3+JjNums2n233c9toSNs0ZpUuII5w/d8BvBCoDP3LuUmADnHPu6tDFCw0twIsEx8Z9x3h34jheOvkse3PXoGj3d8l5RZba5kj8cKkFeH+KSbp9qp1zuwLI5gkVE5HgSUxKYeakN7ll+/PktGQONxpCqZt76xLibCgoV3M553al9wheTBHJimJjori7+6NsvmMu66wapb54iu/H9CAlNfI6kkeyPywmZvZF2p/HzexY2uP4r89DH1FEsoLrrqlF1QHz+LhIb/69sxpd3/6afUdOeR1LMskfFhPnXKO0P/M55/KnPfL9+jz0EUUkqyiUN5a7Hn2B9nd0YX18AjOGP8qOCY9CklrbZ3f+3LR4t5nlS/vnZ8zsUzOrE7poIpIVmRkd65Vhdp9GlMiVSsXtE/hpaENO7s5SO1WIn/y5NPhvzrnjZtYIuA14FxgVmlihoUuDRTJPuaJ5aTlwLJ/WHEFU4hFixt7MntmvQGqq19EkBPwpJr+2DW0FvOmcmw5kqTuV1DVYJHPFRPm4o+N9/NRlPit8dSi+4gXGz5xPcooKSnbjTzHZa2ZvAR2B2Wm9urTtr4j8oWuqVaH2oNkMrzSGIcuT6PjWcn76YZXXsSSI/CkGHYG5QAvn3FGgMDAoJKlEJNvJlysHg7rfyYh76hD38xeU+PBP7BjXE3dWLf6yA79a0GcXumlRxFvxB4+wevxA2p2cwv4cZcnTZTz5ytf1Opb8gaC1oDezLmb2tJk9++sjeDFFJFKUjitE6wFvM63W/8GZ48SOv5Vdc/7jdSwJgFrQi4gnonxG+zu78nPXhSyPqsvwL/Yz7L+btTifRakFvYh46uqqFTkxaBbTpm9gxMJtRK2fTJem1xBXt53X0cQPEdWCXkTCU97YGIZ1qsOIztdw0/HPiJvZnR3je+nO+SzEn2LSCPjWzDab2TozW29m60IVTEQiT9vapYl7ZC4zct9BxZ2T0u6cX+t1LMkAtaAXkbCTnNA/giEAAAzuSURBVJLKzE/fp9GGZ8hvp9nSaSm1atTwOlbEC8rVXMBuoDFwX1oBcUDxIOTLNGqnIpI1REf56HB3d/Z2ns+r0X+h/YQfeWPhVlKSzngdTS7Cn5nJm0Aq0Nw5V8PMCgH/dc7VD2XAUNDMRCTrOJaYxF+nbmDvusWMzjUSaz+SIlff5nWsiBSsmcl1zrlHgEQA59wRslhvLhHJevLHxjCic2163XoNR1NyUuTTjuz44HFIPut1NDmPP8UkycyiOHd6CzOL49xMRUQkpMyMW5s1J7rXYubE3k7FLe8Q/2pjEvdv8TqapPGnmIwApgLFzOzfwBfA8yFJJSKSjvIl4rh54Ad8UuUF8p3azfvjXmfTT9rwNRz41ZvLzKoDNwMGLABuc85luR4IWjMRyfpWrN1An5l7OZKYyquNHK2b3ojFanuJULrUmklAjR7NbLdzruxlH8AjKiYi2cMvJ87w9OSV/Gt3N4jJRc7O4ylQ+XqvY2VbwVqAT/fYAb5fROSyFc2bk7cebMTK+q+RlJREnvdbsnPGv7WbowcCLSaR179eRMKKmdG69R0k3LeQL6MaUP7bl9k5vAXJp497HS2i/GExMbPjZnYsncdxoGQmZBQR+UM1K5aj3qCZTCkxkNWHY+g0bh3xR055HSti/GExcc7lc87lT+eRzznnT9dhEZGQyhMbw10P/Y3oO0ez5cAJHhr+EdsnDdQ9KZlAe7iLSLbTrnYpZvVpzJ151lNp85i0e1K2eh0rW1MxEZFsqWyR3Nz7+Ct8WvUl8p3aTeqoxuxb+q7XsbItFRMRybZionzc0aUXmzt8zhbKUXJhH776ZASB3BIh6VMxEZFsr0HtayjVfyHvF+zNg6tK8/DEb0k4oY23gknFREQiQlyBPHTp8wKPt7yGrzbu5OCr17F7zmugWUpQZLliYmYVzewdM5tyqTERkQv5fEbPmyrx/oN1OWhFKbtiCDveaE/qycNeR8vyMrWYmNlYM/vZzDZcMN4ibTvgbWb25KWO4Zzb4Zzr8UdjIiIXU6tyea4c9DlT4h6m9C/LODysAUd+WOp1rCwts2cm44EW5w+ktbUfCdwO1ATuMbOaZlbLzD674FEsk/OKSDaVP1cO7nz4eRbdOJGTyT5+/HAwSzf/7HWsLCtTbzp0zi01s/IXDDcAtjnndgCY2YdAO+fcC0DrYH22mfUEegKULZvlelOKSAiYGbfdejvbalzDyx99zdfjVvH4jUXo3bgsMQXV4MMf4bBmUgrYc97z+LSxdJlZETMbBdQxs6cuNnYh59xo51w951y9uLi4IMYXkayucpmSjHusHfc0KEullX/j5PCGHFwzy+tYWUo4FJP0Og9f9PIK59wh51wv51yltNlLumMiIv7IlSOKF+6oRe5bn+Fgaj7ipnVhx6SBkJLkdbQsIRyKSTxQ5rznpYF9ofggM2tjZqMTEhJCcXgRyQaaNW5CbO/FfB7bgoqbx7B7WFMSD+32OlbYC4disgqoYmYVzCwH0BmYEYoPcs7NdM71LFBAu7GJyMWVuaIozQdO4pOK/+Ts8cPc/946dhw84XWssJbZlwZPApYD1cws3sx6OOeSgUeBucAm4CPn3PeZmUtE5EI5on3c2b0Puzot4IfjOWj3+hLWTXkBks94HS0sBbRtb1albXtFxB8/JZxm/PjRPHXkWeJjq1Lk/onkuqKq17EyXSi37c1StGYiIpejRIFcDHr0MaZVf4W8p/eSOuom9n0x0etYYSWiionWTETkckVH+WjfuSdbO3zONspScv7DrJ/4lDoQp4moYiIiEqj6ta+mZL+FzMjXmac2lKTf5DWcOJPsdSzPRVQx0WkuEQmGuIJ5adV/FC1uacHMtfv4bGgP4he943UsT2kBXkQkAKu2/UTUxDu51n3PtlLtqdR9JJYzr9exQkIL8CIiIVK/cgnK9Z/H9PxdqRg/nf2v3sCJPeu9jpXpVExERAJUJH8e2vQbyWfX/B8xZ45w9p1WbNi53+tYmSqiionWTEQkVHw+o+0dXYjvNJ8hMX25Y8x3jP9iBy4pMrYH1pqJiEiQHTl5loEfryVu64f0zT2PvPe+T76yV3sdK2BaMxERyUSF8uRgTPd6NKpfj+izR4kZezN7FryVrfebVzEREQkBn89o3f4e9nWex3qrTpllg9k2uhvuzHGvo4WEiomISAhdU6MaVQb8l2kF76P8vtkMHz+RhNPZb4+UiComWoAXES8UzJuLdn2H88mN03ljV1lav76MzetWZKvTXhFVTNSbS0S8YmZ0uvUmJj/UkIpJO6j0SYu0017ZY5+UiComIiJeq1uuEP/p05VZBbtScd8s9g/NHjc5qpiIiGSyQvly0abv68yu/X/EnD1C1Ds3Z/neXiomIiIe8PmM1h26sLfTPDZaFT5auIIJy3dm2Zb2KiYiIh66pmZ1Kjw+j3XlH+Rv079n+DvjOLl3o9ex/KZiIiLiscL5cjP2gesYfFsVWu4eim9MM+KXvut1LL9EVDHRpcEiEq58PuPhZlU5efdHbLYKlF7Yh23v/DnL9PaKqGKiS4NFJNzVqXUVpfvNZ2a+TlTe8zF7hzbmVMIvXsf6QxFVTEREsoKiBfLSsv9bzKj5GstOlabd2xvYeiC827ComIiIhKEon9G244OUve9tjpxO4tE3PmHL+/0h+azX0dKlYiIiEsZurFyUWX0a06Xg91TdNpbdw5qQ+MtOr2P9joqJiEiYK54/lq79XmJalecpdPJHzo5sxIFvpnsd63+omIiIZAHRUT7ad32EDa1nsM8Vpvhn3fn28/e8jvWbiComujRYRLK6hvUbkPeRRXyQuxv3LM7HczO/52xSitexIquY6NJgEckOSscV4a7HX+eeG6ry0Zeb2PJyE35ZP9/TTBFVTEREsosc0T6GtL2S19uVJffZwxT65C52fPoPSE31JI+KiYhIFta8YQN8Dy1mWUxjKq57lR0jWpNy4lCm51AxERHJ4sqXLMb1g6YytUR/Sh1Zwao37uPg8TOZmkHFREQkG4jNEU2Hh4awtPEEnjzRiVYjlrFyy95M2xpYxUREJBu55U8tefORduTN4ePk+13Y9mZHXOKxkH+uiomISDZTo0R+pj96I8fi6lHhwDz2v3ojx3eHdmtgFRMRkWwoX66ctH10KHPrjibm7FGix97MniXjQ/Z5KiYiItmUmdGybUf2dZ7HD1aJHAuHMOWrTSH5LBUTEZFs7uoa1Snbfz4vlxjGuoOhuVs+OiRHFRGRsFIkfx5e7tmBlNTQXN0VUTMT9eYSkUgW5TNyRIfm//YjqpioN5eISGhEVDEREZHQUDEREZGAqZiIiEjAVExERCRgKiYiIhIwFRMREQmYuUxqTxxOzOwgsCudHxUAMnITSlHgl6CGCn8Z/d1khszKEszPCfRYl/t+f9+X0ddn5HWR+D2B7P1dKeeci0v3J845PdIewOgMvu4br7OG6+8mO2UJ5ucEeqzLfb+/7/PjO/CHr4vE70kw/rfOqll0mut/zfQ6QBgLp99NZmUJ5ucEeqzLfb+/78vo68Pp34dwE06/m0zLEpGnuQJlZt845+p5nUMknOl7Elk0M7k8o70OIJIF6HsSQTQzERGRgGlmIiIiAVMxERGRgKmYiIhIwFRMAmRm7c1sjJlNN7Nbvc4jEq7MrIaZjTKzKWbW2+s8ElwqJukws7Fm9rOZbbhgvIWZbTazbWb2JIBzbppz7i/A/UAnD+KKeMbP78om51wvoCOgS4azGRWT9I0HWpw/YGZRwEjgdqAmcI+Z1TzvJc+k/VwkkozHj++KmbUFvgAWZG5MCTUVk3Q455YChy8YbgBsc87tcM6dBT4E2tk5LwFznHPfZnZWES/5811Je/0M59wNQNfMTSqhFu11gCykFLDnvOfxwHXAY8CfgAJmVtk5N8qLcCJhJN3vipk1Be4AcgKzPcglIaRiknGWzphzzo0ARmR2GJEwdrHvymJgceZGkcyi01wZFw+UOe95aWCfR1lEwpm+KxFIxSTjVgFVzKyCmeUAOgMzPM4kEo70XYlAKibpMLNJwHKgmpnFm1kP51wy8CgwF9gEfOSc+97LnCJe03dFfqVGjyIiEjDNTEREJGAqJiIiEjAVExERCZiKiYiIBEzFREREAqZiIiIiAVMxERGRgKmYiIhIwFRMRMKImeUxszfM7Hqvs4j4Q8VEJLz04lyL9kZeBxHxh4qJSHhpAWwB1ngdRMQfKiYiYcLMYoEo4FpgicdxRPyiYiISPqpwrpj84JxL8jqMiD+006JI+IgDqpK2X7pIVqKZiUj4KAl8AvjMrJDXYUT8oWIiEgbMLJpzayVXAKOAFG8TifhHm2OJiEjANDMREZGAqZiIiEjAVExERCRgKiYiIhIwFRMREQmYiomIiARMxURERAKmYiIiIgH7f2LdDgyNQujMAAAAAElFTkSuQmCC\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], - "source": [ - "ell = np.logspace(2,3.5,100)\n", - "\n", - "loglog(ell, ccl.angular_cl(cosmo_ccl, tracer1, tracer2, ell),label='CCL')\n", - "loglog(ell, lensing_cls(cosmo, ell, nz1, nz2), '--', label='jax_cosmo')\n", - "\n", - "legend()\n", - "xlabel(r'$\\ell$')\n", - "ylabel(r'Lensing angular $C_\\ell$')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Ok, so we get a very good match to CCL \\o/" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## But wait! There is more! Let's differentiate" - ] - }, - { - "cell_type": "code", - "execution_count": 22, - "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in asarray is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n", - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype float64 requested in array is not available, and will be truncated to dtype float32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n", - "/home/francois/.local/lib/python3.8/site-packages/jax/lax/lax.py:5190: UserWarning: Explicitly requested dtype requested in astype is not available, and will be truncated to dtype int32. To enable more dtypes, set the jax_enable_x64 configuration option or the JAX_ENABLE_X64 shell environment variable. See https://github.com/google/jax#current-gotchas for more.\n", - " warnings.warn(msg.format(dtype, fun_name , truncated_dtype))\n" - ] - } - ], - "source": [ - "# We define a fake likelihood where the parameters are an nz\n", - "ell = np.logspace(2,3.5,100)\n", - "# Our fake data is the auto-spectrum of one redshift bin\n", - "data = lensing_cls(cosmo, ell, nz1, nz1)\n", - "\n", - "# Here is our likelihood\n", - "@jit\n", - "def likelihood(nz):\n", - " return np.sum((lensing_cls(cosmo, ell, nz, nz) - data)**2)\n", - "\n", - "# And here are our gradients :-D\n", - "grad_likelihood = jit(grad(likelihood))" - ] - }, - { - "cell_type": "code", - "execution_count": 23, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "a : 0.000000 ; b : -0.000000 ; z0 : 0.000000 ; " - ] - }, - "execution_count": 23, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "# So for instance, we get the derivatives with respect to the nz\n", - "# But there are very small ^^'\n", - "grad_likelihood(nz2)" - ] - }, - { - "cell_type": "code", - "execution_count": 45, - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "2.357229 4.0668464 0.5074219\n", - "a : 1.000000 ; b : 2.000000 ; z0 : 0.500000 ; \n", - "a : 2.500000 ; b : 4.000000 ; z0 : 1.500000 ; \n" - ] - } - ], - "source": [ - "# Let's do some optimization\n", - "# Fitting the redshift distribution\n", - "a = 2.5\n", - "b = 4.\n", - "z0 = 1.5\n", - "\n", - "nzs=[]\n", - "z = np.linspace(0, 3, 100)\n", - "\n", - "for i in range(30):\n", - " gr = grad_likelihood(smail_nz(a=a, b=b, z0=z0))\n", - " a -= 0.2e15* gr._params['a']\n", - " b -= 0.2e15* gr._params['b']\n", - " z0 -= 0.2e15* gr._params['z0']\n", - " nzs.append(smail_nz(a=a, b=b, z0=z0)(z))\n", - "\n", - "print(a, b, z0)\n", - "print(nz1)\n", - "print(nz2)" - ] - }, - { - "cell_type": "code", - "execution_count": 54, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "Text(0.5, 0, 'redshift z')" - ] - }, - "execution_count": 54, - "metadata": {}, - "output_type": "execute_result" - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEJCAYAAACE39xMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOy9eXhc5Xn3/3nONqs2y5ItW/KKzWLjjX0J4IQSQtjSFl4SQkqahtImadP0JeRtmy5vSUt/TUnCCwklISENKSEJJCELXIQGsy82iwHbgDfZliVr1+xz1uf3x9EcS9bIlheQl+fjSxeamXPOPCOke+65n+/9vYWUEoVCoVAc+WiTvQCFQqFQHBpUQFcoFIqjBBXQFQqF4ihBBXSFQqE4SlABXaFQKI4SVEBXKBSKo4R9BnQhRFwI8ZIQYq0QYp0Q4p+qHCOEELcLITYJIV4XQqx4d5arUCgUivEwJnCMDbxfSpkXQpjAM0KIR6SUL4w45kPAguGvM4BvDf9XoVAoFO8R+wzoMuw8yg/fNIe/9uxGugL4r+FjXxBC1AshWqSUXeNdd+rUqXLOnDkHtmqFQqE4Rnn55Zf7pJRN1R6bSIaOEEIHXgaOA+6UUr64xyEzgR0jbncM3zcqoAshbgBuAJg1axZr1qyZ0AtQKBQKRYgQYtt4j01oU1RK6UsplwGtwOlCiMV7Pke106pc524p5alSylObmqq+wSgUCoXiANkvlYuUcghYBVy8x0MdQNuI261A50GtTKFQKBT7xURULk1CiPrh7xPAhcBbexz2MPCJYbXLmUBmb/VzhUKhUBx6JlJDbwG+P1xH14AfSyl/JYS4EUBKeRfwG+ASYBNQBD75Lq1XoTiqcV2Xjo4OyuXyZC9FMcnE43FaW1sxTXPC50xE5fI6sLzK/XeN+F4Cn5nwsyoUiqp0dHRQU1PDnDlzEKLa1pTiWEBKSX9/Px0dHcydO3fC56lOUYXiMKJcLtPY2KiC+TGOEILGxsb9/qSmArpCcZihgrkCDuz3QAX094iH33mYi390Mb2F3sleikKhOEpRAf094q9X/TU77Z2c/8PzJ3spCsUh46mnnmLFihUYhsFPf/rTyV7OMY8K6O8itz53K+d9/zzy5Tz/57T/Q9Eu4goX3/cne2kKxSFh1qxZ3HvvvXzsYx+b7KUoUAH9XeU/1/4nPX4PbV9rI2/nGcwNErfifPznH5/spSkU49Le3s6JJ57Ipz/9aRYtWsRFF13Ezp07WbZsWfSl6zrbtm1jzpw5LFmyBE1ToeRwYEJeLor9p7fQSzKeZDA/SMbJcPPjN0MA/hSfZ3Y9M9nLUxwB/NMv17G+M3tIr3nSjFr+4bJF+zxu48aN3H///Xz729/m6quv5oknnuC1114D4M477+TJJ59k9uzZh3RtioNHBfR3ic899jk0TSNbzJI0kliGRcpI0dHbQdkpM1AYYEpqymQvU6Goyty5c1m2bBkAp5xyCu3t7QA8++yzfOc73+Hpp5+exNUpxkMF9HeJF3e9SDwWp1AuMD01nbgRJ6bH6BzoJCDgC499gXs/cu9kL1NxGDORTPrdIhaLRd/ruk6pVKKrq4tPfepTPPzww6TT6Ulbm2J8VEB/F3A8h1gsRrYYflxOW2kSZoIZ6Rn4vk9nuZPHdz6O53kYhvpfoDj8cV2Xq6++mn/7t39j4cKFk70cxTionYx3Acuw+Pr7vo7u6ixrXkZzqpmaWA1JK8nJ008mZsZoqGlg5X0rJ3upCsWEeO6551i9ejX/8A//EG2MdnZ2snr1alpbW/nJT37Cn/7pn7Jo0eR9qlCACG1Y3ntOPfVUeawMuBgqDbGqfRXr+9bTlevikY2PkKhP0DPYQ/fnuyd7eYrDiA0bNnDiiSdO9jIUhwnVfh+EEC9LKU+tdrzK0A8xeTvPvP+cR92/1/Hnv/xzym6Z+kQ9V554JZ87/XO0pFs4u+1sbNcmZsX2fUGFQqGYICqgH2L+9sm/JRVPIaXkW698i6V3LeWOF+/g1a5XiRkxPnryR7n8+MspO2XiVpwgCCZ7yQqF4ihBBfRDzKNbHkVKSb6cxxQmpm6ypnMNd79yN/e8cg8za2dy0XEXUXJKSCn5y1//5WQvWaFQHCWogH6IsaWN7dpIKTF1k4SRIGbESBgJklaSFzte5J3+d3BLLps6N3Hfm/dN9pIVCsVRggroh5iYFaPshB7GMSNG2kqTNJOkzBSZcoaufBd+4HPL+bcAkHUObSegQqE4dlEB/RCSt/OUnTKFcgEdnYSRwNIsLN1C0zRm1MzgA3M/wGB5ENMyaW1qpbmhmclSGikUiqMLFdAPIelYmm1/to0Nf7KBL537Jc6bfR5zGuZQE6shZaQQCJ5ofwI/8Fk2fRmGZpCKp3h4w8OTvXSF4oC47bbbOOmkk1iyZAkf+MAH2LZtW9XjSqUS559//l6dRv/lX/7lgNbw9a9/nWKxGN2+5JJLGBoaOqBrHSwTee57772Xzs7O6PY111zDxo0bD8nzKx36IaTklEhYiTH3bx3cyo7sDgZLg0xPT6e1tpXXdr3Gnzz6J0ypnYKf9XnrL96ahBUrDjeONB36E088wRlnnEEymeRb3/oWq1at4oEHHhhz3J133onnefzlX44vAkin0+Tz+f16ft/3mT9/PmvWrGHq1Kn7vf7J4IILLuCrX/0qp54aSsmffPJJ7rvvPr797W+POVbp0CeRE75zAnPvmkv9v9Rz32v34fkeAHMb5nLe7PO44oQrOKP1DHoKPWhCQ5MamtBoz7dP7sIVimGqWeeWSiU6Ozur2ueuXLmSZDIJwJlnnklHR0fV6/7whz/kiiuuAKCrq4vzzjuPZcuWsXjxYp5++mm+9KUvUSqVWLZsGddeey0AV155JaeccgqLFi3i7rvvjq6VTqf5+7//e8444wy+8pWv0NnZycqVK1m5Muy8njNnDn19feO+FoDVq1ezZMkSzjrrLG666SYWL148Zs2rVq3ivPPO4yMf+QgnnXQSN954YyQzvv/++zn55JNZvHgxN998c3TOvp77pz/9KWvWrOHaa69l2bJllEol3ve+9/H446EVyEEjpZyUr1NOOUUebSz49gI5685Zkn9Ezrptlrzw+xfKu166S24d2Cpdz42OC4JA2p4tP/rgR+XiexfL+n+vl7/b8rtJXLnicGH9+vW7b/zmZim/e8mh/frNzXt9/q1bt0pd1+Wrr74qpZTyqquukj/4wQ9GHXPHHXfIq666asy5n/nMZ+Q///M/j7nftm05bdq06PZXv/pVecstt0gppfQ8T2azWSmllKlUatR5/f39Ukopi8WiXLRokezr65NSSgnIBx54IDpu9uzZsre3d8ztvb2WRYsWyWeffVZKKeXNN98sFy1aNGbdTzzxhIzFYnLz5s3S8zx54YUXyp/85Cdy586dsq2tTfb09EjXdeXKlSvlz372swk/9/nnny9Xr1496rkuvPBCuWbNmjFrGPX7MAywRo4TV1WGfojYmd2JZViUnBI6OjEjhiEMNvRv4DebfsNDGx7iqW1PUfbKCCGwdIu7L7mbXDGHH/h8/7XvT/ZLUCiA8a1zYbd97ne/+91R59x3332sWbOGm266acz1+vr6qK+vj26fdtppfO973+Mf//EfeeONN6ipqam6jttvv52lS5dy5plnsmPHjqjOrOs6f/AHf3DAr2VoaIhcLsfZZ58NsNdpS6effjrz5s1D13U++tGP8swzz7B69WouuOACmpqaMAyDa6+9lqeeempCzz0ezc3No+rqB4qy+jtE3LHmDoQQlJ1yuNlpptA0LVS66BZlv0ymnKHgFNg0sIm6WB1tdW1s6wk3kZ7f9vwkvwLFYceHbp2Up61mnQuMa5/7+OOP85WvfIUnn3xy1LkVEokE5XI5un3eeefx1FNP8etf/5rrrruOm266iU984hOjzlm1ahWPP/44zz//PMlkkgsuuCC6RjweR9f1A34tcj/2DYUQY25P9Pzxfo7VKJfLJBJj99/2F5WhHyIe2fIIEG6MWnooVUxbaVJmKpIupqwUvvTZld+FL32klJw09SSEEGzMbKTslvfxLArF5DCefe6rr77Kn/7pn/Lwww/T3Nxc9dyGhgZ8348C8rZt22hububTn/40n/rUp3jllVcAME0T13UByGQyNDQ0kEwmeeutt3jhhRfGXVtNTQ25XG7Cr6WhoYGampromj/60Y/GPfall15i69atBEHAAw88wLnnnssZZ5zBk08+SV9fH77vc//993P++RMf/l5tve+8884hcapUAf0QMT05nf5sP57vETNi0UALXdcxdRNDM2hONfN239ukrTSz62YjhOC4puM4adZJGLrBCzvG/6VVKCaT8exzb7rpJvL5PFdddRXLli3j8ssvr3r+RRddxDPPhKMXV61axbJly1i+fDkPPvhgpHy54YYbWLJkCddeey0XX3wxnuexZMkSvvzlL3PmmWeOu7YbbriBD33oQ9Gm6ES45557uOGGGzjrrLOQUlJXV1f1uLPOOosvfelLLF68mLlz5/KRj3yElpYW/vVf/5WVK1eydOlSVqxYEW34ToTrr7+eG2+8MdoU7e7uJpFI0NLSMuFrjMc+ZYtCiDbgv4DpQADcLaX8xh7HXAD8Atg6fNdDUsr/u7frHo2yRYBMKcMvN/6SHZkdDJQGaEo1MbtuNkkzSUu6hSF7iNNnnk7eyTNYGuQfn/pH3iq8xbaebdx8+s387fl/O9kvQTGJHGmyxYny6quvctttt/GDH/xgspcCQD6fj8pGt956K11dXXzjG6PCGqtWreKrX/0qv/rVr97VtXzta1+jtraWT33qU2Me21/Z4kRq6B7w11LKV4QQNcDLQojfSinX73Hc01LKSyf2Eo4uCk6BW56+hS+f+2XqEnV8fMnHo8d6C73syu+iv9TPkD1EW20bju+wZXALCSPBn5/y5/zFU39B3IyrDF1x1LJ8+XJWrlyJ7/sTrn+/m/z617/mX//1X/E8j9mzZ3PvvfdO2lrq6+u57rrrDsm19ruxSAjxC+AOKeVvR9x3AfC/9yegH00Z+u2rb+fb67/N9p7tNFlN3PL+W7hq0VXo2uhf3IJTwNItNg9upuAUOGXGKQCccM8J5Eo5DMfgtT97jYZEw2S8DMVhwNGaoSsOjHe1sUgIMQdYDrxY5eGzhBBrhRCPCCGqVveFEDcIIdYIIdb09vbuz1Mf1jywPuyMKztlOrId3PrsrXz2kc/y36//N+1D7QQybEZIWSlM3eT4xuNZ3Ly7kcHxHCzDoiffw4sd1X60CoVCsW8mHNCFEGngQeDzUso9LQJfAWZLKZcC/w/4ebVrSCnvllKeKqU8tamp6UDXfNixs7CTQAY4nhNJFj3fY2P/Rp7Y+gS/2/o71vWso+AU8AMfIQQxI5Q0SSkZzA+SKWRwpcuq9lXKrEuhUBwQEwroQgiTMJj/UEr50J6PSymzUsr88Pe/AUwhxJFhrHAoEOC4DkAkWUwYCWJmDD/w0YVOV76LgdIAqztXk3fySCmj4F6v1zOYH8THZ0PvBt7pf2eSX5BCoTgS2WdAF6Gy/h5gg5TytnGOmT58HEKI04ev238oF3o4Y5kWtmsDoQd6TI+FlrlCI2WlMDSDmB5jqDyEJjSSZpL2oXY2DWxCSsnfve/vMA0TIQSvdb/G/2z+n0l+RQqF4khkIhn6OcB1wPuFEK8Nf10ihLhRCHHj8DF/CLwphFgL3A5cI4+huoFdtBnIDWAJi5geI2Em0DU9Kr84vkPSTJKxM8ysmUnOzjFQGqA2VosQgm25bRzfejypWArXd1nXt26yX5LiGOb222/nxBNPpKGhgVtvDbtVf/7zn7N+/Z7CNsXhxj5li1LKZwCxj2PuAO44VIs60tj2F2H7vuM5/HbLb9me2U7JDTtGdU1HIjG08EfdnGpme2Y7hmYwPT0dgGXTlvGzbT/DNEzqY/U0JhrDUo02+fIuxbHHN7/5TR555BHmzp0b3ffzn/+cSy+9lJNOOmkSV6bYF8rL5SD53qvf47aXbuM/Vv4HFy28iA8v/HD0WHe+m4HSAEPlIaSUpK00QggydobmVHPkE3HerPMIXgywTIuB8gA5J4cXeCqgK95zbrzxRrZs2cLll1/OH//xH7N582Y+9rGP8fDDD/Pkk09yyy238OCDDzJ//vzJXqqiCiqgHyS3vngrVsLikgcu4bSW0/irs/6KKxZeQcyMMS09jWnpaQBk7SyBDBgoDSClZGoy3DPO2Tm6C924notlWAxkB3i953Ue2/wYlx1/2WS+NMVhwCcf/eSY+z4454Ncc8I1lLwSf/74n495/IrjruDK465ksDzIF1Z9YdRj37v4e3t9vrvuuotHH32UJ554IuqQPPvss7n88su59NJL+cM//MODeDWKdxsV0A+Skl9C8zX8wGfzwGa+8fw3WL1zNe+f+36WTVvG9JrpCCGojdUCoUwxaSaJG3ECGdCV7yJmxHDcUIvu+A62a/Nq16usnLuStJXexwoUCoUiRAX0g0TX9UiyGDNipKwUbuCyeXAzA+UBWmtaWdi4kKSZpC5ehxAiCtICQWOikbgRZyg/RECAi4ulW7zd/zYbejdw2szTJvPlKSaZvWXUCSOx18cb4g37zMgVRxfKbfEgiZkxbC+ULFq6FTot6nE0oWEIAyEEnblOduZ2sjO7k65cV3SuEILGZCMpK8VFsy4iVwwtNd3ARUrJSztfikZeKRSTyf5a1ComBxXQD4KuXBemYeK4DoYwiBkxdE0nYSYwNROBIKbHKHklaqwaego9OL4TnZspZ6JrffGcL5KMJdE1ne2Z7QyWB+nKdbEtU32KukLxXnLNNdfw7//+7yxfvpzNmzdP9nIU46BKLgdBjVVDMVskV8hRa9WG3aHDQT1phsG5omQxNINABkxJTKHslRkoDURSRoCXul9iXss8tvdsx/VcGmINxPQYWXtPlwWF4t2lMirt+uuv5/rrrwfgnHPOUTr0IwAV0A+CdCzN5s+F2UogA17ufJn1vespuAUEYa3cD3xiegwv8NCERtpK01sMjclGuipeOPtCvvXmt7AMCy3QsAyLhmRDZBGgJIwKhWJfqIB+EFz2o8t4qesl/u38f+P6Fddz2szTOG3maUgp6Sv2kbEz9BX7aEg0kHNykQ49a2cjS4AK9Yl6PN/DMi3y2Txlr8yM9AxiRkyZdSkUigmhaugHwdqBtTTUNvDJX36Sc+85l1+89QvKXhkhBE2pJo6bchxnzDyD2XWzSRgJ6uP12J6N7dmRjBFgoDRA3slHNrq2b7M9u50d2R08tuUxOvMHPw1cceSg3sAVcGC/ByqgHwSmYUamXNsz27nn1Xv4p1X/xGObH6Mz20kgA4QQmLrJ3Ia5NKWakIQdo5WA7vgOQ+UhamO1uK4bXjOwqbFqEAgGigO8sOMFXN+dzJeqeI+Ix+P09/eroH6MI6Wkv7+feDy+X+epkssBEgQBMSPGUGEIgLgZx9BCmeKOzA7KXpmpmanMrJnJzNqZUXklbsSZXT87uk7eySMQTElMwbZt+nJ9SCQCgS99AIpekYJToD5R/96/UMV7SmtrKx0dHRxNA2AUB0Y8Hqe1tXW/zlEB/QB5o+eNsKnIc9DRw6YiI4WGhiY0EkYCx3coukXe6nuLtJWmrbZtlFEXwJTEFNJWGk1ofPmcL/O5Rz8HQM7JUXSK4SAMGd5WAf3oxzTNUaZYCsX+oEouB8j33/g+ALZrY2gGpjCJGTEM3cAyLCQSTWjEjLDxKGEkGCgN8E7/O/iBP+palm4BcMnCS6hL1WEaJn35PlZ3rWagOACEJZ3B0uB7+yIVCsURhQroB8gt599Cg2xA93WaU82kzFC1kjSSCASaCDN1hkuhtbFasnaWhJGIJIiduU5y9u7uu7f63qKtqY10PE3ey9OcbMbQDOZPmY/jOaO6TBUKhWJPVMnlAElaSZ66/ikg3MDYOriVN3veZKg8RCADdBEG7YAAQzPQNZ2yV4480EtuibJXHqV2WTlnJfJZiWmYGL5B2kqTMBPUxevIlrMU3SK2Z0fzSBUKhWIkKqAfIPPunEemnOHvzvw7Pn/255k3ZR7zpsxDSknGztBb6MXxHUpuiZpYDQWnABAZc2XsDLrQSZmp6JoJM4Hrh0qXfClPV76LRU2LaEm3gISiWyRrZ2kyjp4B2wqF4tChSi4HiDAFNakavvD4FzjnnnP4ybqfUHAK4dDneD0LGhewqHkRrbWtNCWbKLrFcLaoEXaNFt0iNbGayBoAwm5T13MxdRPbt9k6uJVt2W08s/0Zeou96EIn7+Qn8VUrFIrDGZWhHyCWYZEpZBAI+op9/HT9T1nfu55z2s5h6fSlNKXCLLrS3h8zYtTHQ5VKyS0BjPE6LzgFPN8jbsVD2WNyKo3xxqiM4wc+DYnQDmDkG4FCoVCACugHxKb+TRi6geOFLou6FsoWDc1gR3YHgQyoT9Qzt34uKStF0kxi6VakZrF0i/p4fXQbwjp8wS1QKBTYNbALHx83cNF0jbJXpsaqQdM0auPh5mpNrCbcdFUoFIphVEQ4AO5fdz8ArudiaAYJM4GlhV7oQggQ4PouvcVedmR2UHSLDJYGo+6/mBFjSmLKqGsW3AKBDPizFX8WdYXm3Tztg+34gY8XePiBT97Ok3fyo9QxCoVCASqgHxDPdjyLlDLM0DWDlJHC1M1Qjz7sg540k9ieTcpKMVQeorvQjRAC13ejgRgjMTWTtJXmgrkXMLV2KqZuki1leaf/HXqLvZT9MggYsofoK/ZRcAuqPVyhUIxCBfQD4NGPPcon5n2C2cnZHNd4HKZhEjfi6CL0Qa9IFQMZkDSTFJwCSTMJhMOiO3OdY4JxzIhRG6vlqY6nmD5lOnErTsbORMZei5oXsXzacgA0oYUSRn/sG4NCoTh2UTX0A+SL532RL573RSCcPrS+dz09hR50oWNpVlTfjukxHN+JNkeLbpG4EY82NaWU5J08KSuFJjQuaLuAn239GaZhYuomCSMRjqpLNFKfqCdmxOgr9lHySpTcEnFj/8x7FArF0YsK6AfAzP83E8dzuPb4a/mPi/+DlpoWWmpakFIyVB5isDyI67vomh4ZbKXMFK7v4gbuqGYi27fJOTlMPczyz2k9Byklpm6Sc3JszWylpaaF+lg9G/s20pRuwtRMCrKgNkUVCsUo9hkRhBBtQognhBAbhBDrhBB/WeUYIYS4XQixSQjxuhBixbuz3MODVCJFPBbnG6u/wfu+9z7ue/0+MuUMQggaEg3Ma5jHcVOOY2HjQmzPRhMacSNO0S0CROUXCCWMmtCI6WH3Z8yM4fouhmHgei7tg+0Mlgd5o+cNXu1+lfbBdrJ2Fj/wR6lkFAqFYiIZugf8tZTyFSFEDfCyEOK3UsqRAwY/BCwY/joD+Nbwf486Sm4Jy7DIFrLo6AyVhnhk4yN0ZjtZMWMFi5oWMT09PZwrqiVJmkkak40IISh5JUwtLKVAWG4pe2USZmKUrjxqLvJspqWnMT09ne58NwIRZfNTk1MpOAUkctQbhEKhOHbZZ0CXUnYBXcPf54QQG4CZwMiAfgXwXzLc6XtBCFEvhGgZPveo4tHNjyKE2K1wsVIYuoGu6/QV+3i161Xa6tpoTDbSlGzC1M2oNDItNQ0v8KJr2b6NRJIwEqOeo7u/O8zcCTXotmfjBi6BDNCEhqmZNKeaydpZhkpDxI24Kr8oFIr9U7kIIeYAy4EX93hoJrBjxO2O4fv2PP8GIcQaIcSaI9XA/6cbfgoQBfSYHiNhJKKySaVu3pXrIu/k2Z7ZHunKK9OLKniBF26i7lE6ufqkq5FIfHwy5QzbMtso2KFOPZABbuBScAvhG4EIPzUEMnjvfggKheKwZMIBXQiRBh4EPi+lzO75cJVTxoikpZR3SylPlVKe2tR0ZBpM7cjtwPEcHNfB1MywS1SPoQ3/KONGHEEYuN3AJWfn0DWdnJ0b42eettI0p5rHtPEvaVlCy5QWNE3D9m36in3sKuxCEm6WBjKgK9cV+bv0Fnuj+rxCoTh2mZDKRQhhEgbzH0opH6pySAfQNuJ2K3BUTjZ+6hNP8U7vO9z825vJubmw5DJsuqVpGrqmI6UkZaUouaXwfqGRc3JIKSP5YsWPpZony6ptq2isbWQwP0i2nKWlpoWUlWJO3RzOaj2Lnbmd5OxcNESj0uTkGd6oaUgKheLYYp9//SKMOPcAG6SUt41z2MPAZ4UQPyLcDM0cjfXzCgubFvKzj/0MCBuF3ux+k/5SP0KEgy186ZM0k6HnilWDlBLbs0fJFYfKQ0jkGAsAgHn189hS2oKphxuoMS1GTI9RE6vB1E1m1sxE1kg6sh1IJJZu4QYuZa88xvBLoVAcO0wknTsHuA54Qwjx2vB9fwPMApBS3gX8BrgE2AQUgU8e+qUeHsz65ixs1+Z9097HfX9wH7WxWs6edTYA2XKWwfIgA6UBLM3CD3wSZiLa/Kw0AY1Ut1TjI8d/hMe7Hg990d08Gwc3kjATTElM4Y3uN4gZMRY2LkQX4aeBxmRjdE0vUFm6QnGsMhGVyzNUr5GPPEYCnzlUizqcScVTeL7HQ28/xPu//34+fcqnufz4y2lMNlIbr6U2XsusulmUvTJJJ5Qtlr0yQBTQHd+pqm6BMNif3nJ61FzkOOHouTn1c9ie2U7GztCSbmFXbldYavEd4kY8+t72bAxLBXSF4lhE/eXvB3va5jqBw9Pbnqa/1M+KlhUsblpMY7IRXdNJmAnmNoTT20tuKVLAQChXFIiqjUElr4QTOLiei6ZpuL7LtNQ05jXMo7/Ujxd4OL7DztxOpqWm0ZRqwg/8cGDGsMWuQqE4NlEBfT/Y0zbXFCZxM04gAwaKA7zd/zZTilNoSjXRlGyKNjwbEg3RZihA2Stj6daYDdFKrT1mxNjcuRlf+ljCCuWNmo4f+ARBEI6p000QoamX7YX2AZZuEdPUvFGF4lhFpXP7wW/bfwuEGvRow9KIEdfjkb685JVwfId3+t+hr9hX9TopM0XKSo253/EdIDT0+uC8DwLgSpeB8gDre0Lzr4AAOwhdFsteOZxyNNys5PgOZa+svNIVimMUFdD3A5kOxpgAACAASURBVNdzKZQLkQY9ZaXQRTitSBAqXAQCXejR5mTOztGR7RjV+JOyUlVdEh3fQdd0dE2nua6ZtqY2JJKyWyZTzjBQGsD1XBJGgpSZwvEdBkuDZO0slm5FpZzKMAyFQnFsoUou+8Hzf/w8RafI157/Gq91v0bciGNpVqg9R0Y18gpxIx4ZaVVa8x3fQRf6mGMDGeAFXqR8WduzlrpUHZ39nQyVh2hJt1AXq6OltoVzZ51Lzs7Rle+Kmo9SVopsebdM0vEdElp1FY1CoTg6UQF9P0laSf72/L8FwPZs1veup7/Uj+3Z6EKP6uIVK4CyVx6VjQ+UBojpsVE19crxtbHaKPDXWrX0B/2YhomlWcTNsKxTF6tDExo1Vg11U+sYKg9R9sroIrQcsD0bX/rknXzU1KRQKI4NVEDfD467+zgK5QLH1xzPg1c/SGOqkeUty5FSUnSL9BZ68aUfBXFfhrNA62J1QDhnNJABMaP6xuXIrP30GafzSMcjmHqoRV/Xs475DfN5/9z3s3bXWtzA5aSmk5BSRuWcWXWzwuu4Olk7G01K2vPTgEKhODpR6dsEcX2XmBVDInlm+zN8+P4P8+2Xv82u/C6EEGFrfsMc5k+ZT22sloZEwxj9eWVkXMX7vIIXeBScwqg6+5ULrwTANExKXokhe4iyX2awNEh3PnRj3DSwiaHyEAKBG7iRlUDSTFIfr8fUzarWAgqF4uhEZegT5Hftv0MTGo4buiz6vs9ru14j7+Q5peUUTmg6geZUMwCNyUaAqBW/oje3PRtDM8ZkzBV73CS7fc1Pn3k6jheqXrzAozHRyKLmRfQUenB8J9p0jetx5k+Zj67plL0yeSdPY6Jx3C5UhUJx9KIC+gT50fofAbuHT1iGhaGHply9pV60fo3ufDdz6ueQMBNhsDXio9r9Hd8ZM4xCSokbuGN06YZmsLFjIxJJjVkDhHV2BJT8EjEzFnm4jJxPWgn2pm7iBV7oLRP4ey31KBSKowNVcpkgL3e9DOzWoKesFJawMISBhoYudBzfYag8xMb+jUgpRw2zEELQnGoeY57lBV7Y5q+Z7MmSaUuAMNPvyffw0s6X6Mx2ggRDGIhhR4ZMOUNPoSe6huM7BDIgZ+dwfAdNaJGKRqFQHL2ogD5B0maaofwQrudiaRaWZoUKFN2KSiiVrDhuxHEDl+2Z7RScQnSNisZ8JJUMe+TgiwrxeJy2ptCVuOgV8X2fITusmdfEamhOhyWeglMg7+TxZSiPdAMXTYRWvpUuU01oKqArFEc5quQyQZ775HMEQcBjmx/j1xt/HU4M0hMITURllYoRV328PtoQrdTPc3YOQzPG1LZHDojek95SLzWJGnbIHZTcEs3pZmqsGhoTjZw761wgzMbjZhzHdrB9G0u3oo5TQzNCI7Bhoy/bs5Ubo0JxFKP+sifIxt6NLGhawMULLubiBRfz/JYevvvCczyzYSeZ4lvowuQPli7m/BPqmJ6OR7p0UzeRUpJ38iTN5JiAXq1jtIKQAk3T0DQt9Fof/jc1NTUq6cyum42u6eSdfNhMNMLB0dTMSJduaEaUpauArlAcnai/7Aly2cOXkS1maTam8YHmf+OhNS41sQRtU5ZwXD3syvVw22+3cPczNp88U3D16XUkzDDzdgM3GkQxksrQ5/GYEp9CiRKmblJ0i2zo24AXeKyct5LXu18n7+Q5YeoJNCYbo8x8anIqKUKfmErgdv3QTKyau6NCoTh6UAF9AvQX+zF0A8/3WDvwBus6P8v7513NHVd+mgVNs4FQYfLC1h7uenI9dz6xlTXbA27/X7vLIsAYlUneyaMJbdwpQ8ual/F8//OYuknWzlLySniBR87OkXfyBAT0FHroK/YxJTEFN3BHnS+EoDZWG9XtlSZdoTi6UZuiE+DBDQ8ihMD1XEDn+OZ6jp+Z5fH23/A/W/6HnkIPAGfNm8b3P7mSf7lyCS9t8fjcf68nV3ZxfCcqeVQIZIAf+FXVLRUuW3gZRbuIRBIEAQ3xBpa3LKe/2E9vsTeUKQahTDFtpZlRMwMI34AqA6n33IStWPSqDVKF4uhDBfQJ8It3fgGA4zqkDIuWurrQD10zyTk53u57m7f732aoPITru3z8zDnc/r/O4tXtea79zotky86Yckcla6+mbqlw6YJL2dK1hUK5gCe9aJNTExpBECAQyEACu7tQIdxorVxfSknJLeH6YfZeydKVG6NCcfShAvoEeL13AxAG4YQZD9UqVgJNaCSMBEIIDM1gZ3YnWTtL0S3yoZOb+M+Pn8K6ziy3PdoV+blUcH03khOOhxCC6enpQNhNujO3k+d2PMfmwc0ElX8yQAiB44ej6jLlDKZu4ks/esz27VHlGF3TCWQwympAoVAc+aiAvg9Wtw+Qy9XRm+nF8z1iRgxTCx0QK4FcF7sDc8JM0FfsY6g8xIUnTeOvLlzAw2s7+cnLHdE1KwqVvZVbKsel02lap7biBV6olLGSlP1yVIKZXT87UrN4gRdNQ4IRnwI0M8rQAXQRlmFUlq5QHF2oTdG94AeSv3noDZaYt/HrG85lZ24zP173Y4bsITShRZJDS7cgrHxgamFzUUVb/tEzp7Jq41b+4RfrWDGrnuOaaxAibAzal7WtF3h40iNmxaIMe2bNTGrMGmpiNaxoWRF1gQoh8AMf27ejNwrXd4kb8ahUU1HVCCHCkXbSx2TvbyoKheLIQWXoe+FXr3eysSdPS9tDbOp/kxOaTuDvL/h7vnrRV7lu6XXMqps1StNd8VaB3YoWL3D4ykdOJmnpfPa/X6Xshlnxnpuke1LJ4n3fDztQZbiJWXJLSCkjq9yCU6AmVsPU5FQs3Yo2O2usmihTr6xx5EaooRlKj65QHGWogD4Onh/w9cc3csL0Gt4o/Q8ffujDLLx9IWu71qIJjTn1czhn1jmcPvN0Wmtb8aVP3AgbiiC0yPUCD1/6tNbX8tWrl/LWrhx3PbmZslfeZ7nDl+HjQgoM3UAIQdkrs6FvA9sz24kZMd7ofoOXO1+mt9A7yg/G8R1qYjXRm0rFckBKGV1fE5oK6ArFUYYK6OPws1d3srWvwB+eqaHrOo7nsD2znc/85jN8/YWv0z7UTiADEmaCaelpzG2Yy7T0tKj9vrJRCWHmvvL4Zj68pIVvrnqHTb39+5QNVsbW1VnhZqqhG2TKGTzpIYSg5JToK/YhkWTsDOt7149xc6wYfwHUxmqrui2OPEahUBzZqIBeBccL+Mb/bOTkmXW0lx4FRtvmbhvaxmObH+PZ7c+yK78raqe3dIvmVDPTUtPCc3wXwW7jrb+55EQQAV/77ca9yhWBqM69tHkpmUIGCINvnVXHkulLGLQH6Sv2Abs3NwMCpqenEzNiOL5DT6FnlJxxT6SUuL4bfRpQKBRHNvsM6EKI7woheoQQb47z+AVCiIwQ4rXhr78/9Mt8b/nJyzvoGCzxhd9byO+2/Q4IbXMt3cLUzMgvpeyVaR9sZ0d2R5gtS4kmtChYG5pBykpF151Zn+BPzp3Fb9f38OKWwb2uwdItDM3ga7/3NXb07sD13LCFXxgUnSIaoTd6xRa3Yt8LjCq/VNQtfuCTKWdGqV2EEJFfukKhOPKZSIZ+L3DxPo55Wkq5bPjr/x78siYPKSXfeXorS9vqueD4JnbmdwK7A3rMiGHoYTZeyb4936O/2I/t2wyWBiN9d8pKURurHXX9T5zdxoy6NP/0y3V4fnUd+Eh9+NTUVFJm+KbgBR5bhrbwfMfzvNHzBsjdxwohsD2b/mI/HdmOSFJZ2aQdzxO9oklXZReF4shnnwFdSvkUMPAerOWw4KWtA2ztK3DdmbMRQrB06lI6+zsJgoCYHiNpJjFEOI1I08Ifn5SSuBGn6BajGZ/VgqQf+MRNnb+5ZDFv7crx0xHa9JHH2J49KqjPnDaT6Q3TKXklSl6JpmRTlFk3p5pZ2LiQlJXC9m10TccNwmHUI/XnlQA/JqBXNOmq7KJQHPEcKpnDWUKItUAn8L+llOsO0XXfc368poN0zOCSk8MOzYeveRiAwdIgD61/iI5sR7TxKBCh7lx61Bl12J4dDWbO23nyTp7p6elRu72u6dTH67lsSR33PLODO57YxO+vaMUydr+v+tKPSiEVpJSYholEoqMzLT2NhJkgYSY4YeoJmLoZda0aYrfDoqmHA6Yr+nNDMyh7ZaSU0Zoqg6VVhq5QHPkcioD+CjBbSpkXQlwC/BxYUO1AIcQNwA0As2bNOgRPfWjJlV1+80YXVy6fSdIKfzRn3HMGZ7Scwe2X3M6nTvkUUkp6i71sH9pO2S9jaEaUoQ+WB6PySMWQq5rDoRCCz1+4gE9+bzUPvtLBR08PfxZSSvzAHyMnrGzI2q6N0ARD5SEs3WLBlAXomk6mnEEIQWOyMcrIQ5uCcLZpZVRd5bq+9KPAD3v3ZFcoFEcOB61ykVJmpZT54e9/A5hCiKnjHHu3lPJUKeWpTU1NB/vUh5xfru2i5PpcfWorALZrk9fz/HjTj5l922yebn8agOZUM6fOPJVTWk6hMdkYlTNGDmJ2/NGGXFJKsnY2CrgXLGxiaVs9d/xuE44XllcqZY89A7rnexiGgY+P4zls6N3AjswONE1jTecaXu9+ne589+7yCiJ6Q4kb8ehNxdDCodaVAK9QKI4uDjqgCyGmi+GIIYQ4ffia/Qd73cnggTU7WDgtzbK2egAe3fxo6FzoOvQV+/jyqi9z+4u3s3VwK1LKSIN+4tQTw5Z6BDE9huuPHWjhBi5+4I8qdXz+wgXsHCrx4CthLb2iPR+T1cvQlTEgYKA8gGVY4bg7t0ymnEFKSdkrs3FgI/2lfurj9VHW7fhO1OwkhCBpJsdY6kJYoqmoZBQKxZHJRGSL9wPPA8cLITqEEJ8SQtwohLhx+JA/BN4crqHfDlwjj8CC7Nu7cqzdMcTVp7ZFAfWH634IhAqXSqt8d76bx7c8zgs7X6C/2B/5qCTMBHMb5kYacBhtjesFXpTJV7hgYRPLRmTpMSNWdapQvVlPf64/9GQPXFJmiuMajyPrZMnaWSQy/ESgWdieTUOiIZJL5uwcWTs76nrjyRSVfFGhOLKZiMrlo1LKFimlKaVslVLeI6W8S0p51/Djd0gpF0kpl0opz5RSPvfuL/vQ8+M1OzB1we+vaI3ue73ndWBYsmhYpIwUuqYTN+IUnALbMtvYkdlBppwZda2YEaMuVjcqeFfGwI1kZJb+0HCWXq3mftsHbmPXwK6wEchzQULezqOhIYf/hSfvdlisdICaujmqG9TxHbJ2dkzwrmTtKqgrFEcuqlMUCALJw2s7+cAJ05iS2p0h57xcFETjepyUlSKux9E1PezkFOFwZl/6dOW6yDt5YGxDUcV7vJpd7vkLm1g8s5a7nnob16tuB3Dh/AsxRLjBKpG0Z9pZ07mGlzpfwg/8KAhLGWbqBafA9sz2yHlRIiO54ngyxYqqRnmkKxRHLiqgA691DNGbs/nQsFSxwgdnf5Bt3dsQiNCGVjciP3QhRZTV6kKn5IUuiIEMKLmlMYGx0vm5J0II/uTcuWzpy/G7t7qrrm9t91pOmH0CDekGck6OnJOjtbaVtJVGImmrbePEphNpSoYbzZVsvCJdBKIGo8qbUTUvmYqlrkKhODJRAR14bF03hia44PjmUfd/9/Lvkrs5R+lvStz2odtYPHUxCTOBpVtVNxYt3cLxHQbLg6MCpiY0Ulaq6jkAH1zUxIy6BN95elvVx0+aehJSSgzdwA98aqwamlJNxI04lm7RVtdG0kySslK01raSMBOjlC4CMWbARdWALvQog1coFEceKqADj63fxVnzG6lLjC6JtN3exmnfPg1TN7ls4WXcdO5NXHXSVbTVtYVyQETkQV753vGdqOGowr7q0kKTXH/OXFZvG+TV7WM9XizDCrNtw6TslQlkQHehG9uzWdS0iJgeo6/Qx2B5kLp4HbqmR2sBaEo1jbIgMLTwjWHPvWtd0/dpGqZQKA5fjvmAvqknz5beAhedNG3U/et61lFfV8/W/FZmf202j2x8BIDmdDMrWlawomUF9fF6UlYq0pxXLHMr3aIQBvOsnR1XElixCLjmtDnUxA2+8/TWqsd5nhcGWwlFp8hbfW+xM7sTKSUv7HyBt/reojvfTcktUXAKowL6ng1Olm6RttJVN2Ara1IoFEcex3xAf2z9LgAu3COg37P2HiBUuGTsDP/x3H9wx0t30D7YHqlHZtbOZEbNDAzNIGkmIzvakdl5pbSxt2ESuqZTF4/x8TNn88ibXWzvL445ppKhO9JhyB6iPlZPQ6KBvJun6BYj6WJXvovuQnc0xaiyhpHKlr1l4q7vRrp1hUJxZKEC+rpulrTW0VKXGHX/qm2rAHDcMPtOmkkGSgM8s/0Znu94noHSbr+yplQTDYmGUB5YpaFIE9q44+Y0oUXZ/fVnz0HXBN99dmyW7rs+/dl+JBI3CN80WmtbKXtlsnZ2d1YtQ2li3IhHSptABuSdfLQxCmGQH1lXr6DkiwrFkcsxHdC7s2Ve2zE0ptwC0FvqBUIflbgRj3xRLMOi4BTYmd3Jxv6NeP7uzUVTN5mWmjbKG8ULvHGzYSnlqDr2tNo4ly6ZwU9f7iBvj960/NLpX2IgF76JOF5YSsnZoaxypI0uhMHYCzzKXnnUJ4aRAdz2bIru2E8CSr6oUBy5HNMB/fENoUzwokXTxzwmhcTzvcifRdd0kkZonVvJwH3pk7Ez7MjsiM6ryAKBaONxvHKLL/3I/bDCH509h7zt8eAe1rp/tPyPiBkxNBGOxNs6tJWXu17muY7n8KSHJrRwc1YLn9v2bLpyXWTt7Bhv9Mo6K/r4PVHyRYXiyOSYDuiPretmTmOSBc3pMY99dsln2dG9g5SeImEmSJmh7FDTtMjwKm7EsX07ymoz5cyo+nNFrlitoQiqe7csa6tnWVs933+unSDYHejve/M+FrQuIJVIMVAcIO/kWdi4kKmJqQgEs+pmsbxlOa01YadrZaM2siHQzFFSxch5sUppRRPamE8PCoXi8OeYDehl1+f5Lf184MRpVdUeXzrvS+T/T57M32S469K7WDF9BaZuYggDXeiRZa7t2cSMGIEMKLiFUVmwECKqj+9JpQmpmjb9+rPnsKWvwFMbe6P7zpx5JrDbH2Zqaiq1sVosI7z+1OTUcESebjK3YS61sVpM3RyldBk5dKOiNx9Pjx4zYuOqYBQKxeHJMRvQ17QP4ngB5y4Y6/Q7WBxk+u3Tmfm1mQQy4NxZ5/Jnp/8Z1yy+hra6tjBTH97olEhi+m5DrpEbohXNeDUqJY1qjTyXnNxCU02M7z/XHt23ZNqS0D7AMCk4BcpemY5sB7lyjuXTl2PpFrtyu2gfao8cFS3dwpc+gQxIW+lRwzYqQ6irlVb2HLChUCiODI7Zv9pnNvVh6oLT50wZ89h9b95HU10TRb/IgtsX8Ku3f4Uf+DQkGlg6fSmLmhfRkGiIgnHFYXFkQ5EXeJTc0l6dData5QKWoXHtGbN44u1etvYVgDDDrgy6EAjyTp6N/Rvpznfj+A5ru9eyaXAT/cV+8k6e/mJ/9Obi+E7V50lb6Wggx54EMqiqglEoFIcvx2xAf3ZTH8tnNZCKjd2wfOjth4BQTZJ38nznle9w15q7aB9qxw98LN1iRs0M0laaulhdVKse2VC0L/15pTwyHh87YxamLkZl6a4XatFt36bgFJiWnsa0mmkU3ELo+Dhc8u4r9rErv4uYHmNGzQxiejh0I1POkLNz0fXGe0OBsCTkBZ5SuygURxDHZEAfLDi82ZnhnPlVByuxNRPqwB0v1HObeljmWLNzDWt3rY3scmNGjMZkIxAGwFH6c98dpXjZk32VNZpr4nxocQsPvtJB0QnfHHKFHL2ZXjzpUXSKCARNySb8wCfv5AkIg29lIpEXeKMmFrmBi+3v3rSVUlJyS1UzcSVfVCiOPI7JgP78ln6khHMXNFZ93JUuQRDg+R4JM0HMiGHoYfu8L33ah9rpzHbi+E60yTjSL6WS3Y6nbnF9d0KNO584aza5sscvXusE4I8X/zH5Uh5JuKHq+A79pf6wtDNcC/cDP/re9m1KbinKyk3NHBW8hRDYvj1qI3fkY5rQVIORQnEEcUwG9Gc39ZGOGSxpra/6uK7rUfNOwkiEm4xCJ6bHkFKGroSaTke2g4JbGHN+JautVm7Zn1LGKbMbOGF6DT94fhtSSq5fcT2peApNaBi6wc7cTtZ1r+OVXa+MauuvlFgc3yHn5Bgsh4Zfpj7aG72yxmpKFwizdJWhKxRHDsdsQD9z3hRMvfrLf+jShyhlSzTHm0maSRJGYtRgCyFEVNaI6TGGykMMlYei83VNpz5eX7VGXgmQ41npjkQIwXVnzWZ9V5ZXtg9x5yt3Mnf6XGJWjF35XRTcAidPO5mpyanoQmdW7SyWTl/KzNqZGJoRadErbyDVOkZ1oVd1XqysUQihgrpCcYRwzAX0HQNF2vuLnD1O/Rzg1Fmn0vHXHez64i7uvuxuTplxCkkriUDs3kgUYQZr6ia2Z48JiOPVzn3p75cs8MplM6mJGfzg+XaOn3I8EGbauqYzJT6FuBl6oiOgLr577N38KfOZUTMjqutXRuDtWQaKGoyqyBc1ETZRKQmjQnFkcMz9pT63uQ+gqv4c4MH1D9J6ZytT/r8pFNwCS6Yv4bol13HFCVfQVtcW+rnooT95TI9FNetK4JRSkrNz45Yx/MDfryESqZjBH5zSym/e2MXZLb8HgGmY5O08Ja/E1sGt9BZ6WTptKaZusjO7kw29G6LzKwG8Il1sSjWRMHcbkVWycNUVqlAc+RxzAf2ZTf001cSqtvsDfPOVb9KQbsAObJZ/czmPvPMIgQyojdWyqHkR86fMpynZhOM7kf4cdjcUeYE3bjCXUoZeLBMot4zk42fOxvEDnn87LI+YRhiks3aWt/vfZqA0QNEtsnlgMzuyOyi6RYbKQ3TmOqPyULWNTwiz8PHKQxC+Ae3pN6NQKA5PjqmALqXk+c19nDO/cdySyKbBTQDYrk05KPOD13/AD17/AR2Z0CwrbaWpjdUyLT2NGqtmd0ORvruhCKp3gAohIoOt/eG45jRnz2/kR6s7cD0XS7coe2XKXpm59XNpq23DCzy68l3R5mjBKTBYGsT2bGbVzWJKImygKrklduV3TbguXsneVR1doTj8OaYC+pa+An15hzPnVZcrAnh4YRkl8Inroad4f7GfV3e9yrqedRScUNWSNJOht4tmRL7jEGq995wQdCi47szZ7Bwqsau/m+6hbvzAD+1vJdTEagAousXQr4Xd2bTt26M+EVSUKyM3Rl3fJWtnq2bhlT0DFdAVisOfYyqgr2kP/cRPrdLuX8EyLWw3bL5JWAkSZvglhKDklWgfaidn5yh7ZQBSVmqU/twP/HHVLSW3dMCB8cKTpjGtNsZU7Sxs18aRDp7vkbEz7MrtYrA8uFszLsNPCprQsD0b27PpLfSGSpfhte1pIlbxUK+GJjRlp6tQHAEcYwF9kIakyfym6v4lABJJ2SmjoYVZuGaGkkWhR6WVjJ0hU86Mci8EooBZTX9eCbYHqhgxdY2Pnj6LUv406lP10Vp6ij105DpY37MeKWVU26982b6NL31yTg7Hd8IavtBHBe9KeWi8oK3sdBWKI4NjK6BvG+SU2VP2Wg5Z98l1LKtdxsLGhdTH6qPW/7gRRyCI6THcwCVmxCg4BXbld+22pNV00la6ekCX/n5vhu7JNafNgvTjtDa1YhpmqEV3CpzYdCJTU1ORSOY3zOfEphNprW0lbsRHWRJE3uj62I5RXdPHzdB1oe91JqpCoTg82OdfqRDiu8ClQI+UcnGVxwXwDeASoAhcL6V85VAv9GDpzdls7StwzWltez0uFUvx60/8Gikl7/S/w9rutZEhV8Uy15c+MT1G0S2OqpdLKau+WVQy+YMN6NPr4jQlmoGtmEaoRW9INGDpVtQdWhOriYLvjJoZ0bma0KKAHjfiY4K3LvRxlTBCiL0aiSkUisODiWTo9wIX7+XxDwELhr9uAL518Ms69Ly8LWx/31v9/Nx7z2XuXXNpuLWBvJPn+KnHc/Wiq7lo/kU0Jhsj61rYbZlbyX4DGTBUHoqC5khGtuUfLL9/4oeAMMvOOTnKbpmNAxvZOriV46cej6EZtA+2s3bXWrJ2Njpv5PSipJmM6v4VTN3E1My9llXUxqhCcXizz4AupXwKGNjLIVcA/yVDXgDqhRAth2qBh4o17QNYhsbimbXjHrM1u5V0Ik3WyXL+985n1dZV+IFPY7KR4xuPZ17DPHRNx9TMsF0eOaoTE6rLFSvDJg4FnzvzKiBsLpJSMlAaYGPfRrJ2loJToLfYS0euAy/wKDgFtg1tI1POYOnWqGC9pxTR0i1SVmr8DtfAx/ZsFdQVisOYQ1EYnQnsGHG7Y/i+rj0PFELcQJjFM2vWrEPw1BNn9bZBlrXWEzPGz5KFJnA9N1SkeCUeWPcA/cV+3jf7fTSnm0mYCeJGHF/60ezQkQ1F4zUNaUKDQ6RinJqeiud7mLpJ3i3gBA4nNp3IjJoZaJrG5oHNpKxUKE0MXIpukZgRY3p69CDs7kI3CSNBXbxu1P3jlY1G2uketlYAvgdOHowYGHFQI/QUxxiHIqBX+6up+rldSnk3cDfAqaee+p5JJkqOz7qdGW44b95ejxspWawoXPJunte7X2dWeRZtdW0kzASGMOD/b+/N4+Oq6v//55m7zJ59b9o0LV0p0EIpoOwCIlAEFdkUwQUV9339un31p378uoEooiig7MriB0EFBWQtpS3QjbZ0TZp9ksxk9jtzz++Pk5lmT2iTZuE+H495NMm9Mzmnd/Kac9/n/X69dSikMC/gGTszbHbLeLd029e6j5SVwpYuEmnlZx4wA7hwkcwk8yvtVCaFqZn9Glfn0F36oJh5NB1FSpnPa+9LUeCBpgAAIABJREFUzpRsSqzQO3fBnmegdQu0b4WO1yHRBX2dL106uINQVAcVS6FiMcw+AWatBM3Z4HWYmYzHO7sR6LvTWAs0jcPrjhsvN3STsSUr5xaPeJ7bcBOJRzAwCJgBfKbqzWnoKl4diofwGl6KvcXoLh3d7DW2slXfzqEEPbf56Nbd4zafixZexJ2b7gSZJZW16Ex0srNrZ96jJRfbT2aSlHiVJ42UkpZoSz5v3nAZqjCpD5rQSGaTw/5eTQzdg3TCsW3Y+yxsvh92/ge69qifGz4oXwRzTwZ/GbgLwB2ATApSPZAMK/Hf+R945U71HE8RzD8DllwIiy8AfXxCYQ4OU4HxEPS/AZ8UQtwNnACEpZSDwi2TSa6g6Lg5w2+IdsQ6SKQTxFNxPLoHt+bGb/hV/rlU6YoSSSQVochTRDKbxK2586vvXOVoX3Jx6vHOEDm1/lQea3iM9u52wnGbcCqML+mjMdxIma+MVCZFqbeUYm8xbl3Z+0oklm2RzCSVoGsG0pL97ixy/w53t+ESrrwV72EJu4T3w7pb4ZW7IbwPzADUnwonfRLmnQ4l88E1xnHEQrDnv7DjcXj9Mdj8APjKYPnlsPJDUFI/gRNxcDg8jCVt8S7gdKBMCNEIfBswAKSUNwGPoFIWX0elLV4zUYM9WNbu7WJRZZBC3/DCWuYvY/uHt/O9J7/Hsw3P5lexXsOLy+XKhxxyxTrdyW4q/BXoQs//fCC51ewbcVccC3dtuYvK4kq6o920xJvpToU5ec7JlHpLSWQSLC9dzvzi+Qihmkn7TX8+9TKfi97HGz0n3rnwUdbODivoufTNCaVtKzx3A7x6L8isEu+3fQsWnw+m7+Be018KR16sHnYWdj4B62+FF34Dz/8ajrkcTv2iI+wO05pRBV1KefkoxyXwiXEb0TiTtSUb9nZx4fKaUc81dZPvn/V9bGmzp3sPm1o34XKpysqcZW6BuyBfcZkTvXQ2jeEyBm0mZu3siI2YD5a6gjpCXSEM3SBuQSLpxtRMdJeeT0nMxbu9upe5RXPV/DSThJVASonu0il0F/a7e8jl2WfsDG4Gf0AJIcb9w6kfnbvg399Tq2fdCys/CCddB8Vzx/f3uDRYcJZ6RJrhueth7S3w6t1w7AfgzG+Cb/i7OQeHqcqM3x3a2R6lJ5VhxZyR4+c119fg9/qJdkfZ9dldzCueR31RPeFUmIZwgyp/R+LRPcSsWL/sllg6RsAM9BNHKSUSOSEVlhctuoj1L6xXm7ipNNtaQ+zs2olAsHrRajy6h9c7XyeSjFARqKC2oBZQgp4LveTSFAfiNbz5XPuhyPdLHc8wUrwTnvqxElXNgFO+CCdep1bVE01BNZz7Q3jLp+Dpn8JLf4QtD8LZ/xeWX+FkyjhMK6Zo/tn48XKDag23fPbQ/UNz6LoKnbTEW3j3ve9mXdM6AIo8RRxVeRQV/go0oXLQM3ZmUP75QOEWQuDRPRMi6OfUn6NK+nUTHUlLtJOXm7Zh2RY9qZ68JYEQgmQmyb7wPhrCDZiaiVf35vPRbWnnTcZymJo5oljb0h5zT9RRkRI23AE3HAcv/g5WXAmf3gBv+z+HR8z7UlAD5/8UPvqUis0/dB3cegF07zu843BwOATeFIIe9OjMKxvekAsOpCxqaCTSCf665a88tvOxvF1u0B2krqguHxfPldpPlF3uSPhNv/JF102EsBEiQ0e4gKXlSyn2FLO2aa1Kl0SQsBIIlLCbmkl1sDof709YCToTnQdcGnsZSbD7xtkPiY7XlWA+dB2ULYCPPQ2rfwnBqtGfO5FUHQUf/Cesvh6aX4HfnAyv3je5Y3JwGCMzXtBfaejmmNoiXK7hBTeaiuI23KSsFKbLVFkuupvOZCeb2jbR1HMgC9NreKnwV2BoBra0h7TLlVKSzCQnNGe7ub2Zxo5GLJmmyGfzWkuYdEagazqmyySVTYFQ8f1cw+iB4xnKShegJ9UzZO56jpyn+kFh27Dmt3DTydC6UYn4Nf+AyiMP7vUmApcLjvsAfPwZqFgC938Y7r8W0rHRn+vgMInMaEFPWllea+nhmNmFI553/drrlXd4OoXH8FDoKcStufFoHmxpE0lGaAg3kLASQJ/Gyr2r1IGNl7Myqyoux6s8dAi+cNIX1O8iS3HAhWXHuHvDU+zu2k0qm8LKWoN8WVKZFKF4iIZwQ79x93VeBDW/4ZwXQa3SB1oHj4lwI/zpInj0yyp3/Lo1cNzVY089PNwUz4Wr/w5nfAM23ge/P1tt3Do4TFGm6F/S+LC5KUzWliyfPfKGaCgRorOnk3gqjs/w4TN8eA1vPqyiuZQToRCC7mR3XgANzaDQUzio3H+islv6UhmspKa0RuXA6walAcmm/V20x0JEUpG8Q2JtQW2+8jOZSebnYktbneMyBhmKaa6RC4hyaYty6ILgoXntEfjNW6HxJbUqv/I+tSE51dF0OO3LcOVfILIfbj5d5bI7OExBZrSgb9inNkSPqR15hf7Ts3/Kvuv28fOzfs5bZr8Ft+7Ga3jRNA2P7kFzaXkRi1vxfkI2MCc7V0w0Hs6KI/H4nscpCZZgGibt8XbKCuMkU+Xs7dDJyiy1wVqWVSyjMlBJwAxQ5CnC1MyhvdHtwSv0XPeloXAJl8rPH0s+eiYN//g63H05FNepWPlxV0+/7JEj3gbXPgmFs+HOS9QmroPDFGNGC/orjWFqCj1UFHhGPO/+LfdjZSyuO+E67rnkHj534ueoClTl889zjS1SmZTqFNTrthhNRweJ3kQVEw3khJoTAJU7b2NTX1pBmd/LvzY3UuAuoNxfDijhjqVjzCqYhd/05+86cjHygBmgzFfW77X7VoweEuH98Md3wAs3wqqPwoceg9L5h/aak0lJPXzoX7Dg7fDIF+Ff31R7Ag4OU4SZLegN3RwzSroiwFef+yoLfreAd9z+DizbYmHZQi5YeAFLy5eqrj+9+eeprDK7EkJg2RZW1hoUVskVHE101st7l7wXUIKetJJkbIu5VV1s6dhAQ8iiyFPE5rbNbO/Yzq4uFfe1shaaS0MTmto0RYn3wNRKl3ARdAdHtPy1pU0qkxo+jr7nGbj5NGh/Dd57O5z3P8oFcbpj+uGyO+D4D6tq1r9co7xjHBymADNW0EPRFPs646PmnzdFmjB1lbL41N6n+Myjn2FP1x5cwkV1sJq5RXMJmsF8l59cyl9OHAeGHVzCdVi6+9QV1ZHJZjANE8u2iKQiBHwRvKaHh17ZSSwdoyPeQVZmSWfThOIhtoe2k86mKfIU4dW9+deKW/H8hm+O0T6UBKoSdVCsXUp44Sa47UJlhPWR/8DSd47r3Ccdlwbn/T9VfLTlQbjzUicDxmFKMGMF/dXGMMCoK/TvPfM9ZTVrpfDoHrqT3Ty07SE2NG3Ii3a5vzy/knVr7gPVkgOyW2xpH1Z72ZSVUtk5mRTRdJTawhredsRRbG0s4NFtz/dLncytpJOZJIWewn4WuXErTszqL0i2tElYiWHnM6SdbiYN//tp+MdXYOG5SszLF43zrKcIQsBbPw3v/DXsfgpuvwgS3ZM9Koc3OTNW0F9u6MYl4KhZI2+IPrb7MQCS6SQBdwCP7sGjeWiJt7C1YyvRVBRQm4e5/PPcJuLAlbiVtQalAE4kqWiKxo5GsnaWhJXAq3s5YV4hhiZ5ZFMLcCC10kYJb64yNGNn+qVdDhx3Lpd+pPloQjuwhxALqZTE9bfDqV+CS/8MnuG7Q80YVlwJl9wKTRtUoVSsY7JH5PAmZsYK+iuN3SysDOJ3j1x6n8iqVWg6k857heu6Ks7RhEZbrI3OxOAOfANjz4cru6Uvvzj3FwBkyJDKpmiPtbOj61VWzIvxxI6dhBMWGZnJN4h2626SmSRZO8u+8D6iafVhlfd46SPemktDCDHixmgu3JRt3QK/O0OlJL77FmVuNVVzyyeCpe+EK+6G0A64bTVE2yd7RA5vUmbkX52UMl8hOhpe4aWlU61mA0YAv+nHrysf9FzrON2l0xJt6dd2bmBXn5zwTXR2S186U53UVdRhaAYu4SIUD2G6TM5eOouMneDxLW34DT91RXWU+8rx6J58LnrfjdHhKkbHUmCk7XkWccvbwYrDNY/AUe+ZuAlPZY44C664Bzp3O6LuMGnMSEFv6EzQFbc4epQKUYCdn9pJ22fbuO3C21hRvQKP5sGje5QHem+lZ24DMJefPVRmR1ZOfDHRQLa2byXoC6pOS6kIkXSEmoIaFpdXc2RNkH9v7WRu4WJKvCW4dTfFnmIq/ZWA6qCUy0XXXXreNrcvuksfeV9gwx2Yd12Kq7AGPvxvqF05ofOd8sw7Ha68V3VUuu0CR9QdDjszUtA37lcbokfPGnmFviu0i+U3L+fZvc9y1Yqr+M0Fv+FdS99F0B3Eq3tV4wrNjWVbaEJDc2kkM0nCqXC/18mJ/OEMtwC8c6HKHjENk7ilqlyllHQnuzlnaS2pZBm3P78XK2vRHmvHo3vyTaHdmhL03IdTpb+SAnf/mHfO432QoEsJT/xQGWvNPRl5zT+QRYe36feUpf5UVQXbvQ9uf6eyBnZwOEzMSEHf1BTG0AQLqwIjnvfZxz9L1p3l/PvO58Y1N+ISLpZVLOOUulOYXTA7XxGZiz/DAXfFvgghVPPoCbDKHYlT6k5RLe50I78x2hhuZF3zOqoLTU5ZUMuvnn6StftfZl94HzErRsJKkMqk8vPJhV2GurPQXBpFnqL+88pa8NAn4akfwfIrkVfcR9LwTE6v0alK/Slw+V0Qeh3+dLGT/eJw2JiZgr4/zMLKIG595BXz+tb1AKTSKW5ceyO/XvNrIskIpmZSGahkVsEs3JobW9q4NTdZO6vcFV0Tn2c+FnSXTtpKqzz6TIpQIkQoGaLAXYAtbS5bVUJXooe/v7ofUOmJe8N7CSVCuDW3ytrpnUvGztCZ6Bw5SyfVo3KuX/4znPZVeOeNCN0cehX/Zmfe6SrTp3Uz3HGJ+r9zcJhgZpygSynZuD88aroigC2U/W2utdzuyG4e3/143i7XZ/hw624CZgC37h4yXTFrZ0eumJxgkmmVtWJJi6SVxG/4mV04myXlSygtjLO0RvDXDXuwbRdxK95vYzRgBvJhIpdwkcwk8yv2HFZWFS3JSAvcej7sehIuvAHO+Frej6Vf+qLDARaeo1Ia96+Duy4HKznqUxwcDoUZJ+j7uxN0xy2OHIOguw03qbTyZwl6gvh1PwJBQ7iBxkhj3ku8wF2AS7iGrA7NyiwSeVg3Q/tyauWpNIWakEgyMoNEksqk6Ep0ETSDrD6miu5Emn9saiWWjuUFXUqVpphLXcxZFgz0QRdCkG3fQeaWs6FjB1x+Nxx7Vb9zxq3pxUxkyQVw8U3KCuG+D6iQlYPDBDHjBH1T74boaCv0bDaLx/SQtJJ4NS+FZiEBdwCXcOHW3MSteL7RQ2717TN8/Urmc46EhzNVcSC/Xf3b/O9PZVPErTjPNTzH+pb1vBZ6jQVVAZZUBbj3pUZ6UgkMl4GUklRWVZe2xdry88ttlPZF278B8eeLsdJRuPphteocQO4Dzgm7DMPR74ULfgbb/wEPfBScDz6HCWIGCnoEzSVYXBUc8byszNIWaiMUCRF0B/N+LV7Di8ul/FgEglAilA+1aC6tf7gl56x4mLNb+rKlYwvzqudR4CsgYSXoSfdgS5v5xfNJZVLoQuf9JyykO+pj4+7ifP58wkoM2hgdVGC07R+I21ajewrIXPUAzDpu2HG4dfdh3xSeVqz8IJz1Xdj0V3j4cypTyMFhnJlxf4Eb94dZUBHAY4wssqZuEvpSiISV4J5N97C3ey9e3YsudDSh5XPQXcKFqZkkM8lB1aG5RhZj8gWfIJaWLcU0TLyml2g0Sleii+Oqj6PAXUA6myaeiXPOklU8ubmA3zy1m8tPmMvcorl4jQN3GgkrgUf35P3SJRLW/wn+9zNQdRT6pX8i4Q5iS3vYuU7m/8G04eTPQjIMz/wMvMVw9ncne0QOM4wZ9VcopWTTGDdE5984n4qfVRBNRrl6xdV88a1fZEXVCnyGLx92SWWVYVfOqGqoVm2TvSoNuAOkrBRu0006kyaVSWFLm9ZYK3ErTn1RPdWBaj539jxC8Q5++vgL+E1//oMo92EF6k6jzFuK+cwv4G+fVJkaV/8do2CWEvpRVpUZO3PoHuoznbd9C1Z+CJ79BTzz88kejcMMY0YJekskSSiWZtkYBB0dAv4A9b+oZ1vHNvymn2NrjmVZ5TIMl6GqQnt90HNCPtAfXHNpkxpuyZG20nhMDxZqnI2RRjY0b6A93o4tbdpibcTka5y22Mvtazaxv7uHtlgbWTub93kHIJuBhz+HfOL7cMzlqpTdrTJh/KZ/1Lna0nYEfTSEUNa7y94Dj38HXvrDZI/IYQYxJkEXQpwrhNgmhHhdCPHVIY6fLoQICyFe7n18a/yHOjobey1zRxP0SCKCx/SQSCWwsfnek99jQ9MGhBAEzADzS+Zj6gc6FaWz6UHinbEzk5aqOBBTmJi6icvloivZRWuslapAFV7di+bSaI21IqXkA2+Zg2Wnuf4/m2iPtRO34hR7iplTOAfScbj3KpLr/kDLqmuxVl8P2mB74JFwCVfepMxhBFwulfmy4Bx4+PMqru7gMA6MKuhCCA24EXgHsBS4XAixdIhTn5ZSLu99fG+cxzkmNjVFcAlYWj2ybesXHv8CLpeLeCpOwB3Ab/hZ37KedU3rVJs5ISh0F+Zbsw30Ps/lrr+hJskTyEnVJ9Ed7Vbe6NkUmtBwa27KfGWUeEtoi7WRzCSpKfJw9tIq7l23j/3dCeJWXKVhxrvg9gth2yMY5/4P8pTPkx5g1GVlLcLJ8MhmXcJJXxwzmgGX3AZzToL7r3UaTzuMC2NZoa8CXpdS7pJSpoG7gSnZgmbT/jBHVATwmiOHBh7d/SgAiXSCIncRPrePgBEgaSXZG95L3IojhMDQjLzpVt9wS8bOIISYMhuBN59/M40djWSyGRLpRL6RdTQdxcpYSCkRCNLZNFecMBtTz/K7pxqJW3Ho3EXP78+kq/kVuPRPaCd+DE1og9IXc3sFI1WS5v5PnBX6GDF9yna38ki4532w9/nJHpHDNGcsijQLaOjzfWPvzwZykhDiFSHEo0KII4d6ISHEtUKIl4QQL7W3j78T3ab9YZbVjK1CNJPNYGUsCr2FFLoL0TUdXVOug3ErTiQVAZSQFXoK8+GWnPvgZG+G9sVjeFhRtQKXcCl/d2y2hbaxvWM7j+16jNZYK0IIEpkEtYXFXPPWOl7YGWXti/9G/u4skolOIu/+HSxZDai9goEVo0KIUe10YXJTOKclnkJ43/1QWAt3vheaXp7sETlMY8Yi6EOVQA6MNawH6qSUxwA3AA8O9UJSypullCullCvLy8vf2EhHoS2SpK0nNaYN0b9e8FcaWxop0Aso8ZZQYBbgN/x5n3Ap5bDl/LlwwmQWEw1F3BWnrrKOVCZFZ6yTUDxEbWEtC8sWYksbUzM5ouQIFpUt4vNnnsQVhTtY+PyXSJk+PFf9jWz10flVuVt3D7nBmRP0kfYOdJeez293GCP+MrjqIdWD9U8XQ9trkz0ih2nKWAS9EZjd5/taoKnvCVLKiJQy2vv1I4AhhCgbt1GOgc1NakV9ZM3obc9Oqj+Jnq/3sO1T27h48cUUeApw6240oWFqqulyLvsjnAz3EzCJHLWB8mSQyqTwmB7SdhoEBMwAATOAV/eStbMIIajwV6jw0dqb+GHqF9h2Pfce9Uc8lcuAA+3p3Jo7b3fQl+EaYQzFVNkwnjYUzoKrHlSx9T9dpBplODi8QcYi6GuBBUKIeiGECVwG/K3vCUKIKtGrcEKIVb2vGxrvwY7ElmYl6EtHEfQL776Q2htrWXTDIqoKqrh25bW844h3EDBVP9GcUHsNL6lMalDTClMzB/USnQoE9EC+krUt2kbGzhBKhNjcvhmf4WNl9UriyQjtD32Crf/8Mt0Lzubmub/kx/8N0RlT9gV989EDZmCQoOsuHb/pH9VtMmtn+zWodhgjpfPh/Q9CJqk2qbsbRn+Og0MfRhV0KWUG+CTwT2ArcK+UcrMQ4mNCiI/1nvYeYJMQ4hXgeuAyeZiXaJubwtSV+gh6RhabF1tepDhQzM7unfzo6R9h2za1hbUcXXk0Jd4SBCIvWLa0+4UPpvKq88w5ZwKofHTbIp1Ns7V9q/JBT8cIRxpouPM97N5wG+mjLyPyjh9z3dvnkbBb+OaDm/p51ICaa8JKDJqzqZmj3p3ke4062S5vnMql8P4HlIf67RdCpHmyR+QwjRhTmoaU8hEp5UIp5Xwp5Q96f3aTlPKm3q9/JaU8Ukp5jJTyRCnlcxM56KHY0hQZNV0RwDRN0laarJ3l6X1Pc8/me/INHyr8FXh0Dz7DRzqbVpkuveI+1Vednz7+0wC4TTepTIrmnmZMzWRhyUJK7Aztf34X7Q0vkH7rZzBP+QIRK8asYg/vf0sl/9qynxd3ZajwV+RfL51N05XsGrQ5KqUkmUmOuDkqhEBzaU7Ti4OlZgW8768QbVOi7rSycxgjUyPv7hDpSVrsCcXHJOhet5dEOkFAC1DqLaU72c2G5g35VMVib3Fe0PuuRqdaquJA5pfMp7WrlXgyTiqbImbFsLIWnq7dBB7/DiLWQfTMbyIXnosQIm8F/K5ja1lcbfLthzbTFTvQks7UzLxHel+EEMo3fYDN7kBym8vOKv0gmb0KrrgXwo1O02mHMTM11ekN8lqL6gZz5KyRBf36F6/H1E3iyTjF3mKK3EUEjAC2tGnuac6LlxACv+HHralwy1RMVRyKZUXLiKfiJLNJetI9uEM72fXUD9lruLEuvhmtZgWp7IHsHStr4TXcfPm8OYQTFl978BlaY60A+X6qAwUdyG8cjxSC0lya8lJ3VukHz9y3Kv/5rj2OqDuMiRkh6Jt7PdCXVo+csnjHxjtIWSmiySilvlJK/aVqk08z0IRGZ6KThJUAVEZHLqc6tzqf6oJ+17vvwu/2q7FGO4ntfZrtwUpaVryPV9Nd+cYWqWyKMl8ZpmYSNIPUFLn4+Gnz+PvGVv65ZV8+rJQzJhtYZJTzVB8tJ91wGVOmXd+0Zd5pcOW9vaJ+gQrDODgMw4wQ9C3NEUr9JpUFI+c/r/nQGr6x/BuUukqpCdYoH3R3EFMz8+Ktu/R+m4G5sMFUF3OAe7beQ311PQFvgOZMD68HqzAWnM3Rtaso9hSTyCQocBewrGIZdUV1FHoKKXAXEHQH+dgZ9SyrKeeHj25lV4fqVJ/bEB4k6JqBEGLQzweSW6U7HCL1p8KV90H3PtUGMNI0+nMc3pTMGEFfWlMwonikrBTheJgPrPwADV9p4Csnf4UqfxWGy1BOgr3+J+lsut9GoBACj+6ZcoVEQ/HBoiVk7Sx+jx9b6Gil86kKziKWjuE31Mq9yFPUr4zfrbupLajFb7r59RUnkc3C5+9bQyarvM8r/ZUEzMCg35X3TR8FW9ojN552GBv1p6iN0kgz/OFcJ0/dYUimvaBbWZvtLdFRN0Qve+AyTrz7RAI/CpC20pw691TOX3g+Zf4yJdZC5Z4P3AwFJepTeqVp2/DMzwncdSmJVBy/x0+SDKF4iOaeZvZ17+Pllpfx6T6OKDmC5p5m9nTv4dXWV+mIdwBqFV5f5ucb5y1nfUMrv/z3DmD4Un6f4RtS6AeSC804m6PjQN1b4AMPQSoCf3wHtG+b7BE5TDGmvaC/3hYlnbVHLSh6vul5dE0naSX54EMfJJaKUegp5MjyI/EaXjRxoPlzbjM0nU2PGlaYdOKdcNdlylt7yWrcuPGYHjXXTJJoOsqG1g2k7TSRVIRQIkRXoouWaAse3UMkFaE72c2O0A5SmRSXrVzERUcv5IYndvDMjg6klHQnu4mlY0P++tFy83NhF8cnfZyYdRxc/QhIW63UG9ZO9ogcphDTXtDHWvLv8/pIpBJk7SzhVJiHdzxMwkqguTQq/ZWUeEtIZVL5zVBb2vkWc1OWPc/ATSfDridU04RLbuXKZe8DwO/xE0lF6LF6qA3WMrtgNnOL5rK3ey9t8bb8Rm80Hc0XFUXTUTy6hx+/60QWVRRw3R3r2N0RI2NniFmDBT2dTdOd7B41N1936flMIYdxoHIpfPAfytjrttWw7R+TPSKHKcIUVquxsaUpgsdwUV82/O3/q82v4nV7iSajlHnKqAnWkLWzbA9tJ2El8hksukvHo3sAFV+espktWQv+83249QLQPfChf8Gqj4AQfPmEL7O7ZTeReISklaQj3kFWZtFcGpFUhIAZUOEPeaBJRyKjeor2pFX6p9d08cvLF2NoLj5820tYlkHGzgxrqTtaTnruPGeVPo6UzIMPPQYVi+Huy2HdrZM9IocpwLQX9M1NYRZXFaC5ho9xX/LAJbiEi2giSk2whvJAOUWeIjSh0R5vzze18Jt+dJdO1s5O3bzzjh3qVvu/P4EVV8JH/6sqC3txG27OmnMWUkp6rB46451E01Febn6Z1ztfZ1fXLhKZhGqG0XtHEoqHCJgB4lacrJ0llo5huqNcf8UyGrrifOHeLWSyclDYxSVcGJrRL7d9OAzNmNp3O9ORQDl84GGYf6Zq6P3Pb4CzV/GmZlr/hUkp8xkuI1FkFNEebieeijO7cDYVvgo8hgdTN5WAy2y/TbspmXdu2/DCTXDTKdC5E97zR3jnjeAefGfy8VUfp7KoEpfmwhAGnfFOupJdeHQPsXRMNb6wLbyGlzmFc6grqiPoDiKlJJKK4DN8CASLqg3+v4uP4rmdnfz4kV3E04lBYROP7kFKOaoD41RoqD0jcQfg8ntg1bXw/K/UfkoyMtmjcpgkpvVfWGNXgp5kZtT4+ZqPrCGSivC1f30Nl+ZSHujuAnShRMbKWmTsDAVu9TpjTck7bIR2qhXYnqdVH8oLb4Bg1bCn7w3vpbyonJSVYl+O3hUIAAAXgElEQVRkH6ZuMqdoDkF3kFJvab559JKyJf3y72sLagm6g7iEi4AZoCfVw7uPq6M1kuQn/9qCaWj8v3dX9nPI1106mksjlUkNaqI9FLl9iSmdNTTd0HQ47ydQvhge+RL8/iy49M9QvnCyR+ZwmJnWK/TchuiSEVIWv/Pkd5h7/VzIwI2rb+QjKz5Cmb8Mr+7F1E0Ml4Etbby6Nx82mDKeLdkMPPNz+M1boPkVWH298vcYQcwBrjnmGjLZDH6vn6SVJJVVNsAJK0E4GQYJS8qW4BIuMnaGcDLMvvA+Cj2F+XkXuAuQSHpSPXzyzAV88oxFPLCui//78GuDwit+w4/f9I86HSkl6Wx6TH7qDgfB8R9STo3xDvjdGbDp/skekcNhZlqv0Lc2RxACFlcFhz3nplduoqywjNKfl9LzlR6Orj6aumQd7fF2DJeK6+ZiwemsMqeaEh13Gl6Ev38eWjbC4gtUFktB9ZieamiGaoDtCbBf7ieUCNGT7GF7aDu2tKkvqseyLdrj7ewI7aC2oJaOeAcV/goSVgLLtvLOk7kiqy+es4hE2uaWZ7djywzfXn10ft9irG3n+raxs1321PjQnGnMOw0++jTcdzX85RrY9wKc/T0wPJM9MofDwLT+i9raHKG+zI/PHPpzKZqMUhQooifRA1n448t/xJY2hZ5CZgVn5b1KvIbq6pO1s5PfEzPWAQ99Am45G2IheO/tcNkdYxbzHIuLFmPoBgFPgHg6zp7wHhoiDZT7y8nIDIlMgo54B+lsmmQmiRCCUDxEMtObGWNnqQpU5S11hRB88/zFXHpCMX98fiufums9SevAvoMtbXpSPaP7u/TaBjjVoxNI4Sy45hE48Tp48bfwuzOhZdNkj8rhMDCtBX1Lc2TEcMvKP6zE0A26o93UF9UTToXZ0rYFKSVew4tbd+eNuSzbwiVck7dxl0nBczfADcfCK3fDWz8Dn1wLS995UC9378X3krbSGLpBJB0hnAwrO13dQ4W/gnJfOU09yhOkM9GJz/DRmeikyFOElJKuZFd+BZ3bCHW5XPyf847l02+r4+8bG7jqDy8SjithFghsaRO34qOOLRfmcqpHJxDNgHN/CFf+5UAI5tlfOlkwM5xpK+iRpEVjV2LEkv8YMTLZDD3xHhaVLWJWcBaJTILORCdSSjy6h4AZyIdaJqW1nJSw+QG4cRX865tQuwo+/py6TR4ig2WsFHoKWeZfRle0C4mkI95BKpsimUnSFmtjX3hf3p8mK7NYWSvvY+MzfHQlugCIpWPs7d6bX1EHzABXnTSP7188nw37urjwxmfYtD+MECJ/pzOqV7pLm/w7oTcLC86Gjz+vNtMf+5baMG1+dbJH5TBBTFtBf61ZFcEMJ+gt4RZ8bh/hWJhKbyULShdQ7ivHZ/hIZVL50MDAhg6HDSlh+7/gt6eqeKfhh/fdD+/7C5QvGpdf8ZfL/kKlrxJd0+mIdxBJRtjavpXNbZt5oeEF9nTvyeekJyzlxCgQlHhLSGeVVUCu0KozoRwYhRAE3UHetrSUP15zDCnL5l2/eY471+zDcBnKrTIzuHXdQPo6XDpMMP5SlfXynj9AuAFuPl0tHlI9kz0yh3Fm2gr6liblgT5cyKWqsIqzS8+mrbuNo6qPor64nkJPIUEziK4p0YHeRg66+/CJi5Sw/Z9wyzlw5yWQDMNFN8HHnoYj3jauv0oIwayyWcyvnk/UitLU08T20HYkksXliynxltAZ78TKWhxTdQwLShfkLXUL3AX5lMRibzExK5YvLPIZPtyam+PqSvj7p0/mhPoSvv7ARj5x53qiSVfec30sWFnLiacfDoSAZe9WYbwV71PhvetXqApTJwwzY5i2gr61uYdinzGsB3oileBXF/2Kji92cMERF1BbUEuRpwhDU00XfLrv8BpvZS3Y+BflvXLne6GnGS74BXxqHSy/HCboA2VxsdocDfqCtMfayWazlHhL2Bfel7fULXQX5j/QbGnTGmulJliDz/ABUOguxNRMQolQPp5e6ivFa3gpDbi57ZpVfPncRTy+tY23//wZHnk1lDc4GwuOG+NhxFsMF14PH/4PlMxX9Q03nQxb/1cVrzlMa8RkdbJfuXKlfOmllw76+Rf+6hmCHp07PnzioGNH3XQUaSNNW0cbXV/pIp6Os7V9K2X+MoJmEK/hzRe2uDX3xBa5xEKw/lZYewtE9kPZIjj5c3DUe9TG1QSTslIcfdvRpDNp9rbupTpQzdvmvY2jKo6iM95JTUENZ9SfwaLSRbzW8VregbEyUElNsIa2WBtFniIEgqaeJioDlXmhBxVj1106bt3NrvYoX7t/I2t2d3JcXTFffvsiVs4tGvXuJ5VJIZGHP+z1ZkdK2Po35dTZuQsql8GpX4IlqydsgeFw6Agh1kkpVw51bFr+9WSyNq+19LCkanC4pb2nnaSm+mCGU2Haom34TB9HVhypion6eJ0P9D0fN2wbdv8X/voR+NkS+Pf3oGyB6g953QtqRX6YNmDdhhstqxHwBgj6gnTFu9jTvYfnGp4jkUmQsTPsj+wnbsXpSnbRFmvDrblpi7XlC5GaepowNZM5hXP6ibmUMv+8jJ1hXnmAuz5yIv/z7qNp7Ipzyc3/4ao/PsXGxu4Rx5irME1lUo4j4+FECJVF9Ym1cPHNKtPqvg+oTKvnb4TEyNfNYeoxLQV9d0eMdMYeMn6+9PdL8ZgeWjpbWF6+nHVN61RGi+Gh0KNCC1JO0GqwYwc88UP1B3HbahUrP/b9cN0auOohWPQOcB3+//K116wlnoxTVlBGwk6ws2NnPnZd6FF9WPeG95LMJLFsi0gqgm3b7OraRYm3hISVoDXWml9pR9PRvKFZibcEIJ855HIJ3nv8bJ784hl85e1H83JDJ+f/6jHee9Oz/HNzC1l78B1hriG1S7hG3Ux1mAA0HY65FD6xBi65DYI18M+vw8+WwoOfgN1PO+GYacK0rBTd0jx0yf8PnvwB5YXlRBNRkvEk5xxzDkFPkGg6mu8dmsqmxk/MpYS2LbDtEdjykKrqRKh2Yad/DZZeCIb30H/PIeI3/dxw5g1ccNcFALTGW3m19VV0l86C5AKaIk2cUHsCXYkuCtwFWC4LUzfJyiydiU5KfaWE4iGklFQFquhMdGJLm+pANW7dTbGnmFAiRCgRothTjObS8JoanzhjEZceP5s71mznzhf3ce2fQtQWBbhweQ0Xr5jFwsoDFb65zekctnQqSQ87Lg2OvEg9ml+BF38Hmx+El/8MhXPUz5eshlkrJ2Vh4jA60zKG/sNHt/KHZ3az+bvnYuoH3lhlPy+jqriKHU07uGLRFaxespqjKo6izFdGsa94fAaeDKsVy64nYMdj0L1X/XzWSpVFcOTFb7iq83CxqXUTx958LD6vj55YD/WF9RxZcSSmbrK4fDFLy5cST8cpcBdQ4a9gedVyTM3Ea3hpjbbSlexiXvE8BILmaDO2tKkKVOHRPSQzSboSXRR7i/OpjjmydpbuRIQntnXw8CtdPL2jg6wtWVQZ5MwlFZyxqIJj5xSha+paSilJZpJ5SwZH2CeRdBxe+zu8chfsfgrsDPgrVH77vNNVA+tRvIUcxpeRYujTUtCv+sOLtPekePQzp/T7+cZ9Gzn+z8dzUsVJfHTVRzmq6igqfBX4Tf/BpSZKqTqsN66FhjXKF6P5FZBZlTc+92QVRln0jmnzpj7lT6fQbXcTjoVpCjVR7aumrqiOk+eczOKyxUgpCSVDnFh7IifPOTlf/dmd7MZwGVQGKpXXeqqH7lR33qWyzFfWb1WdsBL9cs3z1abCRWskwd9eaeDfW0O8tKeLjC0JunVW1BWzsq6Y4+qKWVTpJ+AVSCnRXBq6S3eEfbJJdKtFzGsPw64nIdkbYy89AmqPh9qVypu/fAmYvhFfyuHgOWRBF0KcC/wS0IDfSyl/NOC46D1+HhAHrpZSrh/pNQ9F0I//weOcsqCMn713OY2hRlbesZJMKkPHlzrY27WX5/c9zzE1x1DuKyfgDmC4jNHFPBmB0A5o3w7tr0HrJmh6WZVNg+oMNOs4mHMSzD9DVXTqo9vFTjVi6Rirbl0FJti2TXNnM+FYmNmB2ayctZKMzFDtr2Z5zXJWL1xNd6Kb3d270VwafsOfX5FHUhEK3AX5nxd7i/NuioZm0BptBVTjbZ/h62etG7fipDIpVYRkwdpdEZ7d2cm6PV1sb+sh95asKnCzqMpHXZmXulIfC8qLmVXso7rAg889LaOFMwc7Cy2vwq6n1GKncS3E2nsPCiidr+x8S+apr4vrobAWCmY5RmGHyCEJuhBCA7YDZwONwFrgcinllj7nnAd8CiXoJwC/lFKeMNLrHqygt/ekOP4Hj1M15xa2Rv9NcbCYoC9Ia1crbZ9tU+fE2vFobrxCoFtJ1SU90Q3xkHpEW6GnBXqaoLsBuvZAbyUkAC5DVWtWL4fqY9Sqo/qYaSngw/HNJ77JXTvuwuf20d7dTmt3Ky7hwmN6kBlJqbuUhaULCXgCzA7OxjRMvLoXj+FhVnAWcwrnUOgpzHvfBNwBKv2VtPf+Ubs1t8ovl1k8ugeXcOXDMba0iaVVr9Lcqt7QDAJmgHDcYt2+Dl5rjvBaSw/bWmLs7YyTsDKI3j18iUWhV6M84KE04KEs4KbY56bI56bQa+A3dYJeg4Bbw2voeE0Nr6Hh1l24DRem5sLQXRguF7om0F3C8Wc/VKRU4ceWjdC6Wf3bsQO6dsPAeg9fqQrbBMrBXw7eEpUf7y0Cd1A9zKBa5RtedTesu9WiSjdB6324dJWp8ybjUAX9JOA7Usq3937/NQAp5Q/7nPNb4Ekp5V29328DTpdSNg/3ugcr6Buf/Atv3/hpSgtKAXUr39LZwvt7evhFYB5kUkgricimVGhkOMyACpMUzobiueqRW1UUzz1saYWTSSQZ4bQ/ncb29u3ErTg+t4951fPyx21bNXZubG8kmozi9/ipKa3JH8+9d/Z37CeRThD0Bqksrhz0exraG0hZKQp8BVQUVQw6vrdVecUUBYooKygbdHx3y26ydpaSYAklwZJBx3c270RKSVlBGUWBokHHX9//OgAVRRUU+PtvpNvSZlfTLgAqiysJ+vpbMWeyGfa07gGguqQav6e/77uVsdjbpvZRZpXNwmv23wRPWSka2hsAmF0+G7fRv+AqkUqwP7QfgLqKOgy9//sulozR3Kn+jOqr6gfdafYkemjtUndD86vnD/pgisQjtHWrhc4RNUcwkO5oNx2RDlzC1e/a5+js6aSzpxPNpVFfVT/oeEe4g+5YN4ZuUFdRN+h4W3cbkXgEt+FmdvnsQcdbulqIJqJ4TS+zymYNOt4caiaWjOH3+KkuHbw3tb+9z3uvZAq+95r6vPd86r2XTCe5r+jDnHjltwedPxZGEvSx3LfOAhr6fN+IWoWPds4soJ+gCyGuBa4FmDNnzhh+9WB0XxHRRBQrY5FIJ1iQNth1xPkE5wZU+pXuQeQ+zU1/76d9AHwlaiXgK4FAhfr5m5wCTwEbPrIh//2axjV8/cmvs61zG5a08o0+cqX5WTtLIpUAodwVc+Ti41k7S8oaXPI/2vFcd6hsduTjmWxmyOM5hj3eO1TLtgYdt6WdP56xBz+/bwWrlR38/L62BVbGGhTnT2fSQ3491PPTmfSgPHwrc+B4ykoNEvRMJtPv+EBBH/j8gWSyvZ5GyBGPD/f83P+PlEM/f7Tjdm86pC3toZ8vsyCGP26jrl9WTtH3Xi993ztWxkIPDv7wGQ/GIuhD3dMMXNaP5RyklDcDN4NaoY/hdw9iyaqziKxyeiZOBCfUnsC/3/fvyR6Gg4PDQTKWtIFGoO+9Ui3QdBDnODg4ODhMIGMR9LXAAiFEvRDCBC4D/jbgnL8BVwnFiUB4pPi5g4ODg8P4M2rIRUqZEUJ8EvgnKm3xD1LKzUKIj/Uevwl4BJXh8joqbfGaiRuyg4ODg8NQjCmZV0r5CEq0+/7spj5fS+AT4zs0BwcHB4c3glN65+Dg4DBDcATdwcHBYYbgCLqDg4PDDMERdAcHB4cZwqS5LQoh2oG9B/n0MqBjHIczmThzmZrMlLnMlHmAM5ccdVLK8qEOTJqgHwpCiJeG8zKYbjhzmZrMlLnMlHmAM5ex4IRcHBwcHGYIjqA7ODg4zBCmq6DfPNkDGEecuUxNZspcZso8wJnLqEzLGLqDg4ODw2Cm6wrdwcHBwWEAjqA7ODg4zBCmtKALIc4VQmwTQrwuhPjqEMeFEOL63uOvCiGOnYxxjoUxzOV0IURYCPFy7+NbkzHO0RBC/EEI0SaE2DTM8el0TUaby3S5JrOFEE8IIbYKITYLIT4zxDnT4rqMcS7T5bp4hBAvCiFe6Z3Ld4c4Z3yvi5RySj5QVr07gXmACbwCLB1wznnAo6iOSScCayZ73Icwl9OBhyd7rGOYy6nAscCmYY5Pi2syxrlMl2tSDRzb+3UQ1dR9uv6tjGUu0+W6CCDQ+7UBrAFOnMjrMpVX6KuA16WUu6SUaeBu4J0DznkncLtUvAAUCSEGd5KdfMYyl2mBlPK/QOcIp0yXazKWuUwLpJTNUsr1vV/3AFtRPX37Mi2uyxjnMi3o/b+O9n5r9D4GZqGM63WZyoI+XOPpN3rOVGCs4zyp9/bsUSHEkYdnaOPOdLkmY2VaXRMhxFxgBWo12Jdpd11GmAtMk+sihNCEEC8DbcBjUsoJvS5janAxSYxbc+opwFjGuR7l0RAVQpwHPAgsmPCRjT/T5ZqMhWl1TYQQAeCvwGellAM7qU+r6zLKXKbNdZFSZoHlQogi4AEhxDIpZd89m3G9LlN5hT6TmlOPOk4pZSR3eyZVhyhDCFF2+IY4bkyXazIq0+maCCEMlADeIaW8f4hTps11GW0u0+m65JBSdgNPAucOODSu12UqC/pMak496lyEEFVCCNH79SrUtQkd9pEeOtPlmozKdLkmvWO8BdgqpfzZMKdNi+sylrlMo+tS3rsyRwjhBc4CXhtw2rhelykbcpEzqDn1GOfyHuDjQogMkAAuk73b4FMJIcRdqCyDMiFEI/Bt1GbPtLomMKa5TItrArwVeD+wsTdeC/B1YA5Mu+sylrlMl+tSDdwmhNBQHzr3SikfnkgNc0r/HRwcHGYIUznk4uDg4ODwBnAE3cHBwWGG4Ai6g4ODwwzBEXQHBweHGYIj6A4ODg4zBEfQHd5UCCG+I4T44sEcF0I81+frn/Q66P1ECHG1EKJmIsbr4PBGmLJ56A4OY6W3yERIKe2J/D1Syrf0+fajQLmUMiWEeBLYxBStvHR48+Cs0B2mJUKIub2e2b9GeXvMFkJ8SQixttdX+rt9zv2GUF70jwOL+vz800KILb3n393n5ZcKIZ4UQuwSQny6z/nR3n//BviBNUKIS4GVwB1CeXN7+5xfIw54dr8shMgKIeom6v/EwcFZoTtMZxYB10gprxNCnIMyaFqFMjz6mxDiVCCGslpYgXq/rwfW9T7/q0B97yq7qM/rLgbOQPlxbxNC/EZKaeUOSikvFEJEpZTLAYQQHwe+KKV8qe/gpJRNQO6cTwCnSSn3ju9/gYPDARxBd5jO7O31kAY4p/exoff7AErgg8ADUso45FfXOV5FrawfRDn25fi7lDIFpIQQbUAlykTpoBBCvBX4MHDKwb6Gg8NYcEIuDtOZWJ+vBfBDKeXy3scRUspbeo8N529xPnAjcBywTgiRW+Ck+pyT5RAWPkI1K7gFuLRPswMHhwnBEXSHmcI/gQ/2+mgjhJglhKgA/gtcLITwCiGCwOre4y5gtpTyCeDLQBFqVX8w9KDuBPrRawN7L/AVKeX2g3xtB4cx44RcHGYEUsp/CSGWAM/3OqtGgfdJKdcLIe4BXgb2Ak/3PkUD/iyEKESt7n8upezufe4b5VbgJiFEAjhJSpno/flbgOOB7/bZpD2vN7bu4DDuOG6LDg4ODjMEJ+Ti4ODgMENwBN3BwcFhhuAIuoODg8MMwRF0BwcHhxmCI+gODg4OMwRH0B0cHBxmCI6gOzg4OMwQ/n96gJDpS3zVDgAAAABJRU5ErkJggg==\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], - "source": [ - "plot(z, nz1(z), label='nz1')\n", - "plot(z, nz2(z), label='nz2 (starting point)')\n", - "for i,n in enumerate(nzs):\n", - " plot(z, n, '--',color='g', alpha=(i)/30.)\n", - "\n", - "plot(z, smail_nz(a=a,b=b,z0=z0)(z), '--', label='fit')\n", - "legend()\n", - "xlabel('redshift z')" - ] - }, - { - "cell_type": "code", - "execution_count": 53, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "Text(0.5, 1.0, 'Lensing angular $C_\\\\ell$')" - ] - }, - "execution_count": 53, - "metadata": {}, - "output_type": "execute_result" - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEeCAYAAABsaamyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3dd1yW9f748deHIYgIorhRUEEQFyjuLThTK9OyzKZ5rJPtYeOc6ndOZaddamblyMxKM0tzlHsP3LgXKg5UFETZ8Pn9caFfM0D2dV/3/X4+Htcj7uu+xhvsvt73ZyutNUIIIRyTk9kBCCGEMI8kASGEcGCSBIQQwoFJEhBCCAcmSUAIIRyYJAEhhHBgkgSEEMKBSRIQQggHJklA2C2l1B6lVHez4yhNSqlYpVSU2XEI+yFJQJSr8nyIaa2baq1Xlse9rEopVVkp9Y5S6rBSKlkpdUwpNV4pVd3s2ET5kCQghANQSrnksa8KsAYIAfpprSsDXQBXwL98IxRmkSQgbIJSqo5S6mel1Pncb6NP3fR+rFLqBaXULqVUklLqR6WUe+57LyulTuV+kz2glIq84ZyoQl6jlVJqe+41Zue+998C4h2rlDqSe/xepdSdhbnPre6llNJKqcAbjp2WXxwFxXBDHC8rpXYBV/NIBB8DF4EhWutDAFrrOK31P7TW0fn97sK+SBIQplNKOQHzgZ1AXSASeEYp1eemQ+8G+gINgBbAQ0qpYOBJoE3uN9k+QGwBt8vrGhWAX4BpQFVgFnBnfhfIdQTjW7M38BbwnVKqdkH3yf1di3Ov4sYAcC9wG1BFa511badSqh4wAnhNa51TzPsLOyBJQNiCNkB1rfX/01pnaK2PAl8Bw2467jOt9Wmt9UWMpBEGZANuQKhSylVrHau1PlLAvfK6RnvAJfe9TK31XGBzQQFrrWfnXidHa/0jcAhoe4v7UJx7lSCGa3Gc1Fqn3rQ/Cjivtd5Q0D2UUqFKqf9XnPiENUgSELbAH6ijlEq8tgGvAjVvOu7sDT+nAJ5a68PAM8CbwDml1A9KqToF3Otv1wDqAKf0X+dVP1lQwEqpB5RSO26Itxnge4v7UJx7lSCGgq5dEzhRiNu8CjQsTnzCGiQJCFtwEjimta5yw1ZZa92/MCdrrb/XWnfGSCYaeK+I9z8D1FVKqRv21cvvYKWUP0ZJ5Umgmta6ChADqPzOKcK9UgCPG17XKmEM+S0YciI3jnyfAUopNyATSMrvGGF9kgSEGVyVUu7XNmAbcDm3EbOiUspZKdVMKdXmVhdSSgUrpXrmPrDSgFSMKqKi2JB7zpNKKRel1O38vVrlRpUwHq7nc2N4GONbeGncawdwX+7foC/QrQxiAFiQ+99xSikvpZSrUqq5UuqbG7qHhgGVgRVFuK6wGEkCwgwLMR7W17Z/AQMxHjrHgAvA1xgNnrfiBozLPecsUAOjCqPQtNYZwGDgUSARuB/jIZmez/F7gQ8xHujxQHNgXSnd62mMv0UiMByYV9ox5J5/GegJNMZoS0gAfgDitdbncw+rD/gBcwt7XWE9SpaXFOLvlFKbgEla66n2dK+iUEqNBHKAH7XWV82OR5QNKQkIASiluimlauVW0TyI0a1zsdXvVVy5bQX1MUotD5kbjShLfxtFKISDCgZ+wujFcwRjANUZO7hXseSOHfi32XGIsifVQUII4cCkOkgIIRyYJAEhhHBglmoT8PX11QEBAWaHIYQQlrJ169YLWus8pwe3VBIICAggOlomNxRCiKJQSh3P7z2pDhJCCAcmSUAIIRyYJAEhhHBglmoTEELYtszMTOLi4khLSzM7FIfk7u6On58frq6uhT5HkoAQotTExcVRuXJlAgIC+Ots2aKsaa1JSEggLi6OBg0aFPo8qQ4SQpSatLQ0qlWrJgnABEopqlWrVuRSmKlJQClVXyn1m1JqilJqbJndKDsLds2GK+dvfawQokQkAZinOH/7YieB3Af3OaVUzE37+yqlDiilDhfiwd4Y+F1r/QgQWtxYbunUVpg7Ej4IhEld4M834OgqyMpzunghhB158803+eCDD/J9f968eezdu7ccI7ItJSkJTAP63rhDKeUMTAD6YTzU781dqLq5UmrBTVsNYDswTCm1nLJcvcivDTy2Anr+C9y8YMN4+HYQvNcAvr8HNn8FF4+V2e2FELbL0ZNAiWYRVUoFAAu01s1yX3cA3tRa98l9/QqA1vrdfM5/AdistV6tlJqjtR5S0P0iIiJ0qYwYTk+GY2vgyDI49Cck5g6mqxYIQb0hqBf4dwIXt5LfSwgHsm/fPpo0aWJ2GLz99tt8++231KtXj+rVq9O6dWu8vb2ZPHkyGRkZBAYGMmPGDHbs2MGAAQPw9vbG29ubn3/+meXLl//tOA8Pj1vf1Ebk9W+glNqqtY7I6/jS7h1UF2PR8GvigHYFHL8YeFMpdR8Qm9cBSqlRwCiA+vXrl06UbpUhpL+xaQ0JR+DwUjj8J2z5BjZOBNdK0LAbNO5jJAavOqVzbyEcxFvz97D39OVSvWZoHS/eGNi0wGO2bt3KDz/8wPbt28nKyqJVq1a0bt2awYMH89hjjwHw+uuv88033zBmzBgGDRrEgAEDGDLE+A5apUqVPI+zV6WdBPJqlci3qKG1jgEK/PavtZ4MTAajJFCi6PKiFPgGGlv70ZCRArFr4NAfcPAPOLDQOK5WCwjuZySF2uHgJB2rhLBFa9as4c4777z+7X3QoEEAxMTE8Prrr5OYmMiVK1fo06dPnucX9jh7UdpJIA6od8NrP+B0SS+qlBoIDAwMDCzppW6tgofxoG/cB/prOL8fDi6Gg0tg9fuw6j3wrGW8H9zfKC24Viz7uISwmFt9Yy9LefWSeeihh5g3bx4tW7Zk2rRprFy5Ms9zC3ucvSjtr7NbgCClVAOlVAVgGPBbSS+qtZ6vtR7l7e1d4gCLRCmo0QQ6PwuPLIYXDsOdX0L99hAzF2bdA/9rCD8Mhx2zIOVi+cYnhPibrl278ssvv5CamkpycjLz588HIDk5mdq1a5OZmcnMmTOvH1+5cmWSk5Ovv87vOHtV7JKAUmoW0B3wVUrFAW9orb9RSj0JLAGcgSla6z0lDbJcSwIFqVQNWg4ztqx0o9rowCLY/zvsXwDKGQI6QchACLkNvOuaG68QDqhVq1bcc889hIWF4e/vT5cuXQD4z3/+Q7t27fD396d58+bXH/zDhg3jscce47PPPmPOnDn5HmevLLXGcKn1DiptOTlwZruRDPYtgAsHjP11I6DJQGOr1sjcGIUoB7bSO8iRmd07yDE5OUHd1sYW+W+4cAj2/Qb75sPSN4ytZnMIHQSht0P1YLMjFkIIwCJJwGaqgwrLNwi6PG9siSeMZLD3N1jxtrFVbwJN74DQO6BGiNnRCiEcmFQHlafLZ3ITwjw4vh7QUCMUmt5pbL5BZkcoRIlIdZD5pDrIlnnVhnajjC35rFE62PMLrHjHKCHUbA7NBhubT4DZ0QohHIAlkoDlqoMKo3Kt/0sIl0/D3l+NbqfL3jK2uhHQfIhRQqhcy+xohRB2yhLDXk0bJ1BevOpA+8dh5J/w9C6IehOy02HxWPioCUwfBNtmQGqi2ZEKIeyMJZKAQ/HxNwanjV4L/9wMXV+EpJPw25PwQZAxMG3vr5Apy/cJkZeOHTve8piRI0denzn0nXfeKfL5np6eAJw+ffr6nENFkZiYyMSJE6+/Lu51SoM0DFuB1nB6m7Ewzp65cCUe3LyNLqct7jFmPJW5jIQNsGLDsKenJ1euXCnzc24UGxvLgAEDiImJufXBRVTUhmFLPDmUUgOVUpOTkpLMDsUcShljEPqNg2f3wohfjBHJe36B6QPgk2bGQjnn9pkdqRCmu/YtfeXKlXTv3p0hQ4YQEhLC8OHDufalt3v37kRHRzN27FhSU1MJCwtj+PDhfzn/ypUrREZG0qpVK5o3b86vv/76t3vFxsbSrFkzwChdhIWFERYWRvXq1XnrrbfyvcbYsWM5cuQIYWFhvPjii3+5TlpaGg8//DDNmzcnPDycFSuMpVamTZvG4MGD6du3L0FBQbz00kul8veyRMOw1no+MD8iIuIxs2MxnbMLNOppbLd9aMxyuvMHWP85rPvEmO205b1Go7JnDbOjFY5s0Vg4u7t0r1mrufFlqJC2b9/Onj17qFOnDp06dWLdunV07tz5+vvjxo1j/Pjx7Nix42/nuru788svv+Dl5cWFCxdo3749gwYNyncJx6+//hqA48eP06dPHx566KF8rzFu3DhiYmKu3zc2Nvb6dSZMmADA7t272b9/P7179+bgwYMA7Nixg+3bt+Pm5kZwcDBjxoyhXr16lIQlSgIiHxU8jIf9/XPg+f3Qd5xRaljyCnwYAjPvNnocSfuBcFBt27bFz88PJycnwsLC/vKwvRWtNa+++iotWrQgKiqKU6dOER8fX+A5aWlpDB06lPHjx+Pv71+sa6xdu5YRI0YAEBISgr+///UkEBkZibe3N+7u7oSGhnL8+PFC/z75sURJQBSCZw2jh1H7x41qoZ0/wK6f4NAScPeGpoMhbDj4RRiJQoiyVoRv7GXFze3/Vgd0dnYmKyur0OfOnDmT8+fPs3XrVlxdXQkICCAtreAvVKNHj2bw4MFERUUV+xoFtdOW5PfJjyVKAg7fJlBUNZpAr7fg2Rij/SCoj5EUvomC8W1gzYeQdMrsKIWwCa6urmRmZv5tf1JSEjVq1MDV1ZUVK1bc8lv3hAkTSE5OZuzYsbe8xs3TV9+oa9eu16ewPnjwICdOnCA4uOzmG7NEErD7cQJlxcnZaDu46yt44SDcPgEqVYdl/w8+bgoz7oTdc6S6SDi0UaNG0aJFi+sNw9cMHz6c6OhoIiIimDlzJiEhBc/z9cEHH7B79+7rjcOTJk3K9xrVqlWjU6dONGvWjBdffPEv13niiSfIzs6mefPm3HPPPUybNu0vJYDSJl1EHVHCEaNksHOWMQbB3RuaDzWqi+qES3WRKDYrdhG1NzJ3kLi1ao2g52vQ/RWIXQ3bvzO2LV9DjaYQfr8x/qBSNbMjFUKUMUkCjszJCRp2N7bURIj52UgGS16BP/8NIf0h/AFo1MOoWhJC2B1JAsJQsQq0edTY4vfA9plGddHeX8HLD8LuM0oIPv5mRyqEKEWWaBiW3kHlrGZT6PsOPH8Ahk43VkJb/T582tJoTN7zC2RlmB2lEKIUWCIJSO8gk7hUMFZAGzEXntkN3cfC+YMw+yH4KAT+eB0uHDY7SiFECVgiCQgbUKWekQSe2QXDfwb/jrDxCxjfGqbeZkxuJ11NhbAcSQKiaJycISgK7vnOmMwu8g24HAdzRxqlg8WvwoVDZkcpHNhnn31GkyZN8PHxYdw4Y9TyvHnzrk8dLf5KGoZF8VWuCV2eg07PGF1No6fC5i9h4wQI6AKtH4Img4xqJSHKycSJE1m0aBENGjS4vm/evHkMGDCA0NBQEyOzTZIERMnd2NX0yjmjm+nWafDzo+Dha/QqinhY1k0WZW706NEcPXqUQYMG8cgjj3DkyBHuu+8+fvvtN1atWsV///tffv75Zxo1amR2qDZDkoAoXZ41/q90cHS5UTpY/xms+xQCo4wuqEG9ZdyBg3h48cN/29cnoA/DQoaRmpXKE0uf+Nv7twfezh2Bd3Ap7RLPrXzuL+9N7Tu1wPtNmjSJxYsXs2LFChYsWAAYK4UNGjSIAQMGmLZ6ly2TJCDKhpOT8dAPjDImq9s2HbZOh1nDwLueUTIIfwA8q5sdqRAOzRJJQCk1EBgYGBhodiiiOLzrQo9XjfWS9/8O0d8Yk9iteNfogtpmJNRrJ3MW2aGCvrlXdKlY4Ps+7j63/OYvSs4SvYNknICdcHY1HvoPzod/bjGqhg7+AVP6wJddjJJCRorZUQo7VNDUzY7OEklA2KHqjaHfe/D8PhjwCWgN858yupkueQ0uHjU7QmFHhg0bxvvvv094eDhHjhwxOxybIlNJC9ugNZzYAJsnw775kJNtNCC3GwUNexptDMLmyVTS5pOppIU1KWWMQvbvCJdPG72Ktk6F7+6CaoHQ9h8Qdi+4VTY7UiHsiny9ErbHq46x3sGze2DwV+BeBRa9CB82gUUvG4viCCFKhSQBYbtc3KDF3fDYMhi5DIL7wZZv4PPWMHMoHF5qVCMJIYpNkoCwBr8IY63kZ/dAt5fh9HajqmhCWyMxZFw1O0KRy0rtjPamOH97SQLCWirXhB6vGMngzi/B1QN+fw4+amKshpZ40uwIHZq7uzsJCQmSCEygtSYhIQF3d/cinecQvYO01iw8sJ1MlxP09u+Nh6tHGUQnTKE1nNwEGycavYpQ0GQgtH8C6rWVAWjlLDMzk7i4ONLSZFpxM7i7u+Pn54erq+tf9hfUO8ghksDuuCTu+vF13HxX4qoqElm/N/c3HUoL3xYoeUjYj8QTsPkrY4qKtCSo29pIBqG3GwPVhHBQNpsElFKhwJtAArBMaz2noOOLmwSupmexYNdpZmxfyeHU5bh47UI5ZeJfqSmzB31HxQrSU9aupF8x1kfe+AVcPAJedaHtKGj9IFT0MTs6IcpdmSQBpdQUYABwTmvd7Ib9fYFPAWfga631uAKu8TywWWu9Rin1m9Z6UEH3LI3BYicSUpgVfZA5+xeQmHYV99TuDGhRm3TvuQwJ6U37Ou1xUtJUYhdycuDQH8b6BsdWg2slCB8O7UZDNZlKWDiOskoCXYErwLfXkoBSyhk4CPQC4oAtwL0YCeHdmy7xSO5/3wBSgI5a604F3bM0Rwzn5Gg2HktgTnQcC/fvwcXvM5RLCl4uNRjSeDD3hQ6hZqWapXIvYQPO7DJKBrtnQ04WhNxmVBX5d5R2A2H3yqw6SCkVACy4IQl0AN7UWvfJff0KgNb65gRw83Wcgbla69vzeG8UMAqgfv36rY8fP17sePOTnJbJrztOMH3nfOKyV+FS6TCgGNPkIx5p3RMXZykZ2I3ks7Dla6NbaepFqBMOHZ6E0DvAWaoFhX0qzyQwBOirtR6Z+3oE0E5r/WQB578KVAK+0FqvLeh+5TF30OFzyXyzKZpFxxaQeLorNb0q0bzJHhrW1DzSchi1KtUq0/uLcpKRkttuMBESDhtrHLR/HMJHgLuX2dEJUarKMwkMBfrclATaaq3HFPsm/GU9gccOHSqfRcwzs3NYtu8cP0WfZH3SF7h4R6OAEO+2/CN8OD3qd8VZVseyvpwcOLgYNoyH4+vAzctYG7ndaGMdBCHsgM1XBxWWWbOInk1KY8qmaOYemkuq+wacXK7gX6EbH/V8l8Y1ZUIzu3FqK6wfD3vngXKCZkOg4xio1ezW5wphw8ozCbhgNAxHAqcwGobv01rvKfZNbmD2VNI5OZrVh84ycctv7DiWQ8bV+jT3z8Kj1hL+2XoEXeq1l3EH9uDScaMRedu3kHkVGvWEjk9Bw+7SiCwsqax6B80CugO+QDzwhtb6G6VUf+ATjB5BU7TWbxfrBn+9V7lXB91KwpV05m47xbQdi0iq9C3KJQVPpzoMCRrCyPC78XaTVdAsL/USRE+BTV/ClXio1Rw6PSONyMJybHawWFGZXRLIi9aatYfPMH7zz8QkL8ap4gmcdEVebvo9g8Ma4O4q7QaWl5UOu36E9Z/DhYNQpT60/ye0GgEVKpkdnRC3ZPkkYIslgbxcuprBFxtWM3//Js7EheHl7kKjkKXcFtyGB1rcQQXnCmaHKEriWiPyuk/h5EZj9HHbUcZWydfs6ITIl+WTwDW2WBLIi9aaTccuMn3jAVZffR2nCudx1pXpXLM/L3Z8GH/pdWJ9JzbCus/gwO/gUhHC74eOT4JPgNmRCfE3kgRMdD45jY/WLmDxiblkusWggB4+L/FS18H4+chsppZ3/gCs/wx2/gg6B5reCZ2fMdoPhLARkgRsQE6OZl5MDF9u/55DByMgx52w4FO0aujMM+2G4ekmdcuWdvm0MfAseipkXIHAKKMROaCz9CgSprN8ErBKm0BhnUpMZdamE8w4Mo6cStGonIq0rNKLsZ0eoWkNmdjM0lIvGVNSbJoEV89D3Qjo8hw07gdOMv2IMIflk8A1Vi4J5CU9M5vJW5bzw75ZJDlvBTQNXfvz325jaeFXxezwRElkpsKOmUa7QeJx8A02qomaD5W1DUS5kyRgAetjj/DhhunsO+FJyqWmtKhfgeaNY3mp8714u3uaHZ4oruwsYwTymo/g3B5jjqKOY4w5iipIm5AoH5ZPAvZWHVSQy2mZ/Lw1jq+2/Uiy10zIrkgzr9680mkkLWoHmB2eKC6tjbUN1nxkdC/18DUmrGszEipKqU+ULcsngWvsuSRws+zsHKZvX8m0Pd9yUW8DFDWcIni78zu0b1BDpqewsuPrjWRw+E9jwro2jxprG3jWMDsyYackCVhc9MkjvLdhCvvPx5J88j5Ca3sRGXaVf7Tvjqebu9nhieI6s9NIBnt/BRc3aPWAMUdRlXpmRybsjCQBO5GSkcW87aeZsmEHZ6v8C5VdmfAqt/F610doXF1WQbOsC4dh3cew8wfjdYthRiOyb5C5cQm7IUnAzmTnZPN19EJm7PuOJPaic1yp69KFl9o/SWSQPDgsKynO6E20bboxX1HTO6DL8zLwTJSY5ZOAIzUMF9Xq2J18sPFrjqWu58qR52hRO4D72lfnjhaNcJPJ66zpynnYOAE2fw0ZydC4L3R9Efzy/AwLcUuWTwLXSEkgf+euJLFk9yWmrovlTMXPcK2QQo/aQ3mt2zBqVJbRyJaUegk2f2WMRE69BA26GclARiGLIpIk4ECys3MYt346vxz5nnR1Fp3pTahnf17t/BBhfnXMDk8UR/oVY12D9Z/D1XNQr72RDAIjJRmIQpEk4IBydA4/7vmTL3dMJSF7D+nx/WjvO4SRnRvQJchXuphaUWYqbJthTGV9OQ7qhBvJILi/JANRIEkCDm5j3E7W7tP8uPkCl9iMT/WDDA8ZwegO3WXRGyvKyoCds2DtR3ApFmo2g64vQJPbZX4ikSdJAgKA9Kxs/r18MotOT0GrNFRaIL3q3s3L3e6gRuWKZocniio7C3bPhjUfQsIhY36iri9A08Gy/KX4C8snAekdVLqS05P5aOMMfj32A5nqEtlXmnF77Vd4tHNDAmvIPEWWk5NtzE+0+gM4txeqNoQuL0CLu2WyOgHYQRK4RkoCpSszJ5PpO+ex7tBV1u+qSXpOKqGN9/B8hwfo2bi+tBtYTU4O7F8Aq/8HZ3dDFX9jGuuW94GLLG3qyCQJiFtKuJLOW8u/Z8Wlj9DZbnhndWF0+EMMa9UcV2epZ7YUreHgElj1HpzeBl5+0OVZY+ZSFzezoxMmkCQgCm1HfAzvrv+CvUlr0DhRITWCfzR9keHtGuHpJvXMlqI1HF5qJIO4LVC5DnR+1pijyFXmnHIkkgREkZ1IOsm7679k+5kjnD14P5XdXRgU4c5TXdtR00seIJaiNRxdaSSDExugcu3cZPCgJAMHIUlAFJvWml1xSXy+eisbM55Hp9UnwucuXu1+J8G1vMwOTxSF1hC7BlaOg+PrwLOWkQxaPwiu0jvMnkkSECWWkpnCNzt+YMa+b0nVCWSn1SLIbSAvd7mHDg1lfQPLObbGKBnErgHPmrnJ4CFJBnZKkoAoNZk5mczeP58vtn9NYmYcV468QMtajRjdrSG9Qmvh7CTJwFJi1xolg2vJoNMzEPGwJAM7Y/kkIOMEbE+OzmH72Rj2xnrx1ZpjxLt9SxXXWowOH8G9EcEyEtlqYtfBqnFwbLUkAztk+SRwjZQEbFNaZgb3L3icA5c3o7PdcLnaiftD72dUpzC83GWwkqXcnAykmsguSBIQ5WJfwj7+t/ELos+vRGtnVPwD3NeyL492akAN6VFkLbHrYOW7udVEtYxBZ9KbyLIkCYhydfzycT7Z/DUp8b34MyYZV484eoXU54UeXQnwlbUNLCV2Lax4F46vNcYZdHnOGGcgg84sRZKAMM3xhKs8sOghLmbvJyu5KW2qDOWlHlE0q+ttdmiiKI6tNpLBifXgVddIBuEPyHQUFiFJQJjqUtolJu+czk8HfiBDXyXrShChFYfyYvfetGtQVbqXWsW1QWcr34WTm8C7njFradhwmajOxkkSEDbhSsYVpu/5nukx35J5oReX4iMIq+/FP7sFEdmkJk7SvdQatIYjy4ySwaloY6K6bi9Bi2EyhbWNkiQgbEpaVhqZWTBv+1k+i57OVdcN+Gb159lOdzKwRV1cZMI6a9AaDv0JK96GMzuMKay7jYXmQ8BJugjbEkkCwmYtOrqEcZs+4mLGabLTauKV3ocn2w5haER93FzkQWIJWsOBRbDiHYjfDb6NoftYCL1TVjqzEZIEhE3Lysli8bHFfBL9JfFpsWQmtcAr+WFGdW3IvW3rU0lmL7WGnBzYP9+oJjq/D2qEQo9XIWSArIFsMptIAkqphsBrgLfWekjuvkrARCADWKm1nlnQNSQJ2LccncPyE8s5eV7xx7ZKbDx+HK9q+3ioxd080ikI74rS+GgJOdmw5xejATnhMNRqAT1eg8Z9JBmYpMRJQCk1BRgAnNNaN7thf1/gU8AZ+FprPa4Q15pzQxIYASRqrecrpX7UWt9T0LmSBBzLe+u/4btDn5CT6YVK6sG9oUP5R5cQfD2lj7olZGfB7p+MieouxULdCOj5GjTsIcmgnBWUBApbYTcN6HvTRZ2BCUA/IBS4VykVqpRqrpRacNNWI5/r+gEnc3/OLmQswkG81OERJveaTLMajVC+v/L96dF0/fIN3vwthrNJaWaHJ27F2QXC7oMno2Hgp5B8FmbcCdNug+PrzY5O5CpUEtBarwYu3rS7LXBYa31Ua50B/ADcrrXerbUecNN2Lp9Lx2EkgkLHIhyHUooOdTrw46BvmdpnKuE1Q6ld8xwzNp6g6/9W8NLP2zh5McXsMMWtOLsa8w89tQ36vW9UEU3tZySEU1vNjs7hFbpNQCkVACy4Vh2klBoC9NVaj8x9PQJop7V+Mp/zqwFvA70wqo7ezW0TGA+kAWvzahNQSo0CRgHUr1+/9fHjx4v0Cwr7kp6dzrmkbN5fvprlSW+SmdiRPn5DeLpnC7mSoAsAABmvSURBVBpV9zQ7PFEYGSmw5WtY+zGkXoTg24wG5FrNbn2uKJZSaRjOIwkMBfrclATaaq3HlEbQeZE2AXFNbFIs72x8nw1nV6Oz3cm81InutQbzXGQ4wbUqmx2eKIz0ZNj4BawfD+mXodlg6P4K+AaZHZndKY02gbzEAfVueO0HnC7B9fKllBqolJqclJRUFpcXFhTgHcDkPhOYPXA23ep3pILvMtalvUyfT5fxjxnRxJyS/1dsnltlY6Tx0zuMuYgOLIYJbWHeE3BJSvzlpSQlARfgIBAJnAK2APdprfeUSaRISUDk7+Clg0Sf3s3Z082Zuu4Yae5raV+rB8/1DCe8vo/Z4YnCuHLeqCLa8jXoHKMdoesLULmW2ZFZXml0EZ0FdAd8gXjgDa31N0qp/sAnGF1Ep2it3y61qP96f1lZTBTazvgD3L94KOS4kn6xA62r3MHzka2ICKhqdmiiMJJOwer3YfsMcHKFdqOMlc485N+vuGxisFhpkJKAKKwjiUeYsH0Sf55Ycj0ZhHvdxbORLWjfsJrZ4YnCuHjUWP94109G1VGHJ6HDE8bPokgkCQiHdTTxKBN2fMHauE1kn3iZC5ehTYMqPBsZTIdG1WQaayuI32tMUrd/AXhUg87PQZtHZcnLIrB8EpDqIFFSVzOv4ow7Mzcd49N9T5J2OZBQj4E8H9WKzoG+kgys4NRWWPYfOLrCWOWs20sQfr+sZVAIlk8C10hJQJRUUnoS/934DotjF12vJmpScSDPRbWia5AkA0s4tgaW/T+I22xMX939VWh2l8xYWgBJAkLc5GjSUSZu/4Ilx5dATgWuxo6mZc0mPB0ZRLfG1SUZ2Dqt4eASWP4fiI+Bms2g579kkrp8WD4JSHWQKCtHEo8w+8DP1NVDmLQylrPpewj1Deb5qDBJBlaQkwN75sLy/8KlY1CvHUS+AQGdzI7Mplg+CVwjJQFRlpLTU4mcHUVqZgbpCZ1oXPE2no9qSXdJBrYvO9PoUrrqf5B8BgKjIPLfULul2ZHZBEkCQhTSgYsHmLjjC5afXAY5FUlP6ExIxf48F9VCSgZWkJkKmyfnzkt0CZoOhp6vQ7VGZkdmKkkCQhTR/ov7Gb99IqviVuBx4Sniz9chvH4VnolqLA3IVpCWBOs/hw0TISsNWo2Abi+DVx2zIzOF5ZOAtAkIsxxLOkbdSv7M2RrHh1s+5XKKM009+/F8r2bStdQKrpyD1R9A9BRwcoa2o6Dzsw43+tjySeAaKQkIs+ToHMYse4rVp1ahsiuTer4bzb368nyvpnSUQWe271KsMfp45w/g5gWdnoL2j0OFSmZHVi4kCQhRSrbFb+Pz7eOJjt+CyvbmatxQWtdsy3O9Gst0FFYQv9foVnpgIXjWNAactXrQ7gecSRIQopRtPrOZCTsm0qriP5i57irnUxJoH1CX53qF0kYmqrN9JzbC0jfhxAbwaWA0HjcdbLcDziyfBKRNQNiytMxs7p43ktikWFLORdK+ZiTP9QqhlUxhbdu0hkN/wNK34NweqNUCot6ARpF2N+DM8kngGikJCFu14sQKPt8+nkOJByGjJqnnouhcpzvP9QqhhV8Vs8MTBcnJht1zYMV/IfEEBHSBXm9B3dZmR1ZqJAkIUQ5ydA5Ljy/l8+3jib18DHXxDi7Ht6dXaE2ejWpMaB0vs0MUBclKh+ipxloGKRcg9Hbo+W/wDTQ7shKTJCBEOcrOyWbhsYW0rt6JOVsS+GrLn6RkZNKnUWeejWpMUE2ZD9+mpScb6x6v/zx3jMED0H2spVc4kyQghIlG/fE4G86sRac2IvVcbwY27sDTUY1p4OsY3RMt68o5o1QQPRWcXIwFbTo9De7eZkdWZJIEhDBRenY6cw7O4cudk7mUfpGcqyGkn+/D4KZtGNMziHpVPcwOURTk4lFY/jbEzIGKPtD1RWgzElzczI6s0CyfBKR3kLAHKZkpfL//e6bsnkqQ691s2NEYjeaeNvV4skcQtbzdzQ5RFOT0Dlj2FhxZDt71oedr0PxuS3QrtXwSuEZKAsIeJGck4+7izoXkLJ5f/AXb43eRfbEXIyLCebx7I3w9rfMN0yEdWWGMMTizA2o2h6g3IdC2u5UWlARsP4UJYWcqV6iMq5Mrtb0r0jPUi4o+u3Fv8D4zD39M1w9/5f0l+0lKyTQ7TJGfRj3gsRUwZApkJMPMu2D6QGP5SwuSkoAQJjt79SyTd01m7qG5aO1Mypn+uKd1ZlSXhjzcuQGebi5mhyjyk5UBW6fCqvcgJQGa3mmsY1C1odmR/YVUBwlhAScvn2Tizom0qNKTZdt8WHrgOD4eFfhnt6bc394fd1dns0MU+Um7nDt19XjIzoCIR6DrS+BZ3ezIAEkCQljSS8vf5s8Tv3Mlvhs+2V15qmcod0fUo4KL1OLarOSzxmyl274FVw9jttIO/zR9tlJJAkJY0O7zu/l026dsOrsJlxwfks/2pJZzJ56LasLtYXVxdrLdhkiHd+GQ0Xi8fwF41jIGm4WPAGdzqvYkCQhhYRvPbOTTbZ8ScyEGz/RunDnaj6AanjzfuzF9mtaStQxs2YlN8Oe/4OQm8A02ehIF9yv3nkSWTwIyTkA4Oq01y08ux79yAAdOevC/ZWuJuxJHqE8EL/YOoYsseWm7tIb9v8PSNyDhMNTvCL3/A355PpPLhOWTwDVSEhDC8J8N/+Wngz/inB7E5dO9aFsnjBf7hNDaX6avtlnZmUZbwcpxcPUchN5h9CSq1qjMby1JQAg7k5GdweyDs5m8azIX0y7ilNKc5DNR9GzUnOd7B9OktsxYarPSrxg9idZ/DtnpEPGoscJZJd8yu6UkASHsVEpmCjP2zmDqnmk0dO9BzK5uJKdnMahlHZ6NakyATFJnu5LjYdU42Drd6D3U+Rlo9zhUKP25pCQJCGHnEtMSjTaBbA/eWvobS44tJf1CD+5u1YSnI4Oo6SXzEtms8weNnkQHfofKdYw5iVreC06lNy5Epo0Qws5Vca+Ct5s33h6uNG94GZcq66kc+D5zj06h2weLeHfRPhJTMswOU+SlemO493t4eBF41YZf/wlfdoXDS8vl9lISEMIOHU08yufbP2fpiaW4UpnkM31wT2vP6G6NeLhTAB4VZCoKm6Q17PnFmK30Uiw07GH0JKrVvESXleogIRzU7vO7+WTbJ4R6t2ffgXCW7jtDNU83no5szLA29WX0sa3KSoct38Dq/0FqolE91PN18K5brMtJEhDCgWmt0WiclBOfbJzJd/u+41JcL2q7teT5XsHc3rIuTjL62DalXoI1H8GmSVC/Azz4W7EuU1ASkDKhEHZOKYXCeMiH1/Xjj9OKdKeppGUF8fyvvfhyVVNe6htMj+AaMuDM1lT0MaqD2oyEzNQyuYWUBIRwMJnZmfx86Gcm7ZxEQloC7ld7cv5Eb9oE+PBy3xAiAqqaHaIoZTZRHaSUagi8BnhrrYfkt68gkgSEKD0pmSlM3zud4CpNOH2mAZ8s301CSjKRQYG82CeE4FqVzQ5RlJISdxFVSk1RSp1TSsXctL+vUuqAUuqwUmpsQdfQWh/VWj96q31CiPLh4erB4y0fp6d/d+5v78/9fY5RpfEHbEr8jr6fL+G5n3YQdynF7DBFGStsm8A0YDzw7bUdSilnYALQC4gDtiilfgOcgXdvOv8RrfW5EkcrhCgzdzW+gzNX41ioF1LVZzOLTnRnwYcduL9dIE/2DKRqpQpmhyjKQKFKAlrr1cDFm3a3BQ7nfpvPAH4Abtda79ZaD7hpK3YCUEqNUkpFK6Wiz58/X9zLCCFuoV7lerzX9T1+GvATbWq3xKX6Aho0WcS09cfo+r8VfLbsEFfTs8wOU5SyknQSrgucvOF1XO6+PCmlqimlJgHhSqlX8tt3M631ZK11hNY6onp121iqTQh71qRaEyb1msRXvb9ifP/nWPJMV1o3hE/Xz6fr+yuYsfE4mdk5ZocpSklJuojm1Zcs31ZmrXUCMPpW+/K80f+tJ1DUGIUQxdS+dvvrPzcL3cU2vsUlO5g3FvfimzWNeaFPMLc1ry3dSi2uJCWBOKDeDa/9gNMlCydvWuv5WutR3t7eZXF5IcQtPNPqGca2HYubxzkqNRhPivd0xsxeyu0T1rH+8AWzwxMlUJIksAUIUko1UEpVAIYBxRvOJoSwaa7OrgxvMpyFgxcyqsUotMceurWJIeFKBvd9vYkHpmxmz+kks8MUxVCocQJKqVlAd8AXiAfe0Fp/o5TqD3yC0SNoitb67TIJUpaXFMKmnEs5h7NyppJLFd5f+Sez9ywn+VwH7mjZgOd6NaZe1dKfE18Un00MFisNMlhMCNszYccEJu2cREUnH5LPRJKd1JoRHRryZI9AfKRbqU2wfBKQkoAQtm1b/DY+3Pohu87vwlP5ceFEb9wzQxndvRGPdGpAxQqlt0CKKDrLJ4FrpCQghO3SWrP0xFI+3fYpnWr25cihdizdd46aXm48G9WYIa39cHGWqavNIElACFFuMnMy0VpTwbkCEzb9zHe7F3D2eE8Cferzct8QIpvIbKXlTaaSFkKUG1cn1+s/+3ilk+Ueg3fQDhKvdmHkzM60rV+PV/qFEF7fx8QoxTWWKAlIm4AQ1hV/NZ6JOycy7/A8XFVFsi8M4FJ8OLc1r82LfYIJ8K1kdoh2T6qDhBCmO3TpEB9t/Yioev05cSKYyasPkZmtub99A8b0DKSap5vZIdotqQ4SQpguyCeIL6K+QGuNClY4+axg9r7fmbmzF3O2xvG49CQyhSWa6pVSA5VSk5OSZESiEFZ3rVE4xDcAj4ppuNefjE/DGXy4YjXdP1jBT1tOkp1jnRoKq5PqICGEadKy0vh+//d8tesrrmamUC19MMeORhBcszJj+4XQPbi69CQqBdImIISwaZfSLvHlri/pG9CX0/E1GbdkGycSMujYsDav9m9Cs7oyeWRJSJuAEMKm+bj7MLatsUJtWA3YfDmaxceWsyc+igGfn+eOMD9e6BOMn4/MSVTaLNEmIIRwLLcHDaJRVT9yfH+iXvNJLD66mp4fruLdhftISs00Ozy7YonqIBknIITj0Vrzx/E/+Hjrx5y6coogl3vYHhOOd0VXnuoZxP3t/angIt9jC0PaBIQQlpWRncGs/bOI8o8i8bInby5cQ3RsMvW9a/Fy3xD6Nasljce3IG0CQgjLquBcgQebPghAXU+oEbCIqm6bybzSkye+T6RVvRq8dlsTWvtXNTlSa5KylBDCUsa2HUu3el1I9lhI7WafEpu+kru+WMfj320l9sJVs8OzHEkCQghLqedVj4+6f8SMfjNoWKUumVV/IKr9IVYdPE+vj1fx1vw9XLqaYXaYliFtAkIIy7rWeNypTidS0lz495KF/BmTiIeqy5iegTzYMQA3F5mGoqA2AUuUBGTaCCFEXpRS9Anog2cFT2p4uXOl0lw8G32Gr/8C3l0STeSHq5i/8zRW+rJb3qQkIISwGxfTLvLFji+YfXA2rk5uVEjuxakTEYT5Vef125oQEeCYjceWLwkIIURhVHWvymvtX2Pu7XNpX7stlz1+ZXjkJc4kpTJk0gYe/24rxxOk8fhGUhIQQtitbfHbaFm9JelZmteWzGLJrqtkptTngQ4BjOkZSBWPCmaHWC5knIAQwiG1qtkKgIqumuP6F5z9DuHn0o6pm7sxZ2scT0UGMcLBRx477m8uhHAYSim+6/cdo1uO5pLeiXfgx/jW/5P/LNxG749XsTjmrMM2HksSEEI4BA9XD/4Z9k/m3zmf/g37cc5pMWPv8MDV2YnR323lni83sisu0ewwy520CQghHFJsUiwB3gFkZefw4h8TWLUHLl5oyB1hdXixbwh1q1Q0O8RSY/k2gRtmETU7FCGEnQjwDgBAq2yOpS8js/oRguuGsehgJItizvJYl4aM7t4ITzdLPCaLTUoCQgiHl5mdyff7v+fLnV9yNSuF2qo7B/Z3pFrFqjzfuzF3R9TD2cm6M5XKOAEhhCiAq7MrDzZ9kN8H/87djYdyQa1l/P2NCajmwStzd3PbZ2tYe+iC2WGWCSkJCCHETS6mXaSqe1W01jyx+E22HaxC/JkgeobU5NX+IQTWqGx2iEUiJQEhhCiCqu7G9BIpWSmczdhFSpUpBId9x5ZTu+nzyRre+DWGi3YyU6kkASGEyEcl10rMHjib19u9TiqnUX6fENJsEd9F76Pb+yv4avVRMrJyzA6zRCQJCCFEAVycXLgn5B4WDF7Ag00f5IrzLn4a1Z7W/j68vXAfvSw+2EzaBIQQoghSMlPwcPUgR+fwwIIniT0eSFxcY9o1qMa/BoTSrK632SH+jbQJCCFEKfFw9QCMxuM0fY6kylMJCf+WA4l7GTh+LS/M3kn85TSToyw8SQJCCFEMvhV9+XHAj7zZ4U1SiSe71ic0a7mQX3cdpscHK/ls2SFSM7LNDvOWyjUJKKUaKqW+UUrNuWHfHUqpr5RSvyqlepdnPEIIURLOTs7c1fgufr/zdx5t9iguFU+zaExPugZV56M/DxL54Up+3XHKptsLCt0moJSaAgwAzmmtm92wvy/wKeAMfK21HleIa83RWg+5aZ8P8IHW+tH8zpM2ASGELcvMycTVyZWUzBTuW/AoiWfbEns8iPD6PvxrQCit6vuYEldptQlMA/redGFnYALQDwgF7lVKhSqlmiulFty01bjF9V/PvZYQQliSq5MrAAmpCVRwySbBYwpNWn3LiSsHGDxxPU//sJ3TiakmR/lXhU4CWuvVwMWbdrcFDmutj2qtM4AfgNu11ru11gNu2s7ldV1leA9YpLXeVtxfRAghbEU9r3rMum0Wb3V8ixQdT0bNjwkLX8TiPSfp8cFKPvrjACkZWWaHCZS8TaAucPKG13G5+/KklKqmlJoEhCulXsndPQaIAoYopUbncc4opVS0Uir6/PnzJQxXCCHKh7OTM4ODBvP7nb/zSLNHqFk1k6XPRtK7aS0+W36IHh+s5OetceTkmNteUKRxAkqpAGDBtTYBpdRQoI/WemTu6xFAW631mNIPVdoEhBDWlaNzcFJOxF+N54GFI8lKiORIbCAt/Krw7wGhRARULbN7l+U4gTig3g2v/YDTJbzm3yilBiqlJiclJZX2pYUQolw4KeNxm5SRRGU3N85V/Iamrb/jTOohhkzawJhZ2zllQntBSUsCLsBBIBI4BWwB7tNa7yn1SJGSgBDCPmTnZDP38FzGbx/PpbRLBHlEErMzCnBiVNeGjO7WiEqluJhNqZQElFKzgA1AsFIqTin1qNY6C3gSWALsA34qiwQgJQEhhD1xdnJmaOOhLLhzAQ+EPkALvyosf6EnvZvW4vPlB+n5Yfm1F8jcQUIIYTKtNUop9iTs4ellL8DFARyO9aelXxX+PbAprf1LNr5A5g4SQggbppSxdGVWThaVKlQgvuIkWkT8yOmUY9z1RdmOL7BEEpDqICGEI2hZvSVzBs1hbNuxJGQeJbPWh3SI2MDimLO8OGdnmdxTqoOEEMIGJaYlMmHHBKp7VKd/vftJy8wu9rKWBVUHlV7zsxBCiFJTxb0Kr7V/rczvI9VBQgjhwCyRBLTW87XWo7y9bW/FHiGEsDJLJAEhhBBlQ5KAEEI4MEskAWkTEEKIsmGJJCBtAkIIUTYskQSEEEKUDUkCQgjhwCyRBKRNQAghyoalpo1QSp0HjufxljdQmAzhC1wo1aBsX2H/NuWhvGIpzfuU9FrFPb+o5xX2+MIc54ifE7Dvz4q/1rp6nu9orS2/AZMLeVy02bHa6t/GnmIpzfuU9FrFPb+o5xXhM3DL4xzxc1Ia/9ZWjcUS1UGFMN/sAGyYLf1tyiuW0rxPSa9V3POLel5hj7el/x9sjS39bcotFktVB5WUUipa5zOTnhDCIJ8Tx2IvJYHCmmx2AEJYgHxOHIhDlQSEEEL8laOVBIQQQtxAkoAQQjgwSQJCCOHAHDYJKKXuUEp9pZT6VSnV2+x4hLBVSqkmSqlJSqk5SqnHzY5HlC67SgJKqSlKqXNKqZib9vdVSh1QSh1WSo0F0FrP01o/BjwE3GNCuEKYpoiflX1a69HA3YB0HbUzdpUEgGlA3xt3KKWcgQlAPyAUuFcpFXrDIa/nvi+EI5lGET4rSqlBwFpgWfmGKcqaXSUBrfVq4OJNu9sCh7XWR7XWGcAPwO3K8B6wSGu9rbxjFcJMRfms5B7/m9a6IzC8fCMVZc3F7ADKQV3g5A2v44B2wBggCvBWSgVqrSeZEZwQNiTPz4pSqjswGHADFpoQlyhDjpAEVB77tNb6M+Cz8g5GCBuW32dlJbCyfEMR5cWuqoPyEQfUu+G1H3DapFiEsGXyWXFAjpAEtgBBSqkGSqkKwDDgN5NjEsIWyWfFAdlVElBKzQI2AMFKqTil1KNa6yzgSWAJsA/4SWu9x8w4hTCbfFbENTKBnBBCODC7KgkIIYQoGkkCQgjhwCQJCCGEA5MkIIQQDkySgBBCODBJAkII4cAkCQghhAOTJCCEEA5MkoAQpUApVUkpNV4p1d7sWIQoCkkCQpSO0RhTLXc2OxAhikKSgBCloy9wENhhdiBCFIUkASFKSCnlDjgDrYBVJocjRJFIEhCi5IIwksB+rXWm2cEIURSOsLKYEGWtOtCY3PV4hbASKQkIUXJ1gJ8BJ6WUj9nBCFEUkgSEKAGllAtGW0AtYBKQbW5EQhSNLCojhBAOTEoCQgjhwCQJCCGEA5MkIIQQDkySgBBCODBJAkII4cAkCQghhAOTJCCEEA5MkoAQQjiw/w9OOS2ySqBBlwAAAABJRU5ErkJggg==\n", - "text/plain": [ - "
" - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], - "source": [ - "# And let's compare the fitted cls\n", - "nz = smail_nz(a=a, b=b, z0=z0)\n", - "loglog(ell, data,label='data')\n", - "loglog(ell, lensing_cls(cosmo, ell, nz2, nz2), label='initialization')\n", - "loglog(ell, lensing_cls(cosmo, ell, nz, nz), '--', label='fit')\n", - "legend()\n", - "xlabel(r'$\\ell$')\n", - "title('Lensing angular $C_\\ell$')" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.8.2" - } - }, - "nbformat": 4, - "nbformat_minor": 2 -}