Fair inverse transform implementation

inverse-transform.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Inverse transform"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "See evaluated notebook at: https://gist.github.com/newkozlukov/9663aca1ea70b03ec52385eb881a30ff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "import numpy as np\n",
    "\n",
    "\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example CDF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAErCAYAAADEyxRmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAwXklEQVR4nO3dfZwU1Z3v8c9vEAQCDgKuioBjUEEeJ4KKiuu4GCPB+JTcJQkmkpuEuDFuTKJrEvcqcWMe7jWum1dismxMSFYSNUZ3NFFjcNUVA/I4AgooyOOiBnAlKhBBzv2jarAZzsx0z9TpUzP9fb9e/Zru6nr49uma/nVVV9Ux5xwiIiJNVcUOICIi+aQCISIiXioQIiLipQIhIiJeKhAiIuKlAiEiIl4qECIi4qUCIdIOZvZxM1tkZm+a2ctm9rCZTTCzGWa2x8zeSG8vmNkPzOzogmnrzGxfOm3j7cGYr0ekkAqESBuZ2ZeB24BvAUcCg4HbgYvSUe52zvUG+gKXAEcBiwuLBLDFOder4Pahsr0AkVaoQIi0gZlVAzcBVzrn7nPOveWc2+Oce9A5d23huOnw54ApwFbgKxEii5RMBUKkbU4HugP3FzuBc+4doB44K1QokSypQIi0TT9gm3Nub4nTbSHZ5dRogJm9XnD72+wiirTPIbEDiHRQ24H+ZnZIiUXiGOC1gsdbnHMDs40mkg1tQYi0zTxgN3BxsROYWRXwIeCpQJlEMqUtCJE2cM7tMLMbgB+a2V7gUWAPcC5wDrCzcVwz6wocD8wgOZLp1rIHFmkDbUGItJFz7lbgy8A/khydtAn4AvAf6ShTzOxN4HXgAZLdUmOdc1vKHlakDUwdBomIiI+2IERExEsFQkREvFQgRETESwVCRES8VCBERMSrU50H0b9/f1dTUxM7hohIh7F48eJtzrkjfM91qgJRU1PDokWLYscQEekwzGxDc89pF5OIiHipQIiIiJcKhIiIeHWq3yB89uzZw+bNm9m9e3fsKBWne/fuDBw4kK5du8aOIiJtEKxAmNlPgQuAPznnRnqeN+BfgA+SXPlymnNuSfrc+elzXYCfOOe+09Ycmzdvpnfv3tTU1JAsUsrBOcf27dvZvHkzxx13XOw4Ih3T7Nlw/fWwcSMMHgw33wxTp5Zt8SF3Mc0Czm/h+UnACeltOvAjADPrAvwwfX448DEzG97WELt376Zfv34qDmVmZvTr109bbiJtNXs2TJ8OGzaAc8nf6dOT4WUSrEA45/6LA3vOauoi4BcuMR/oY2ZHA6cCa5xzLznn3gbuSsdtMxWHONTuIu1w/fWwc+eBw3buTIaXSczfII4huX5+o83pMN/w05qbiZlNJ9kCYfDgwdmnzMArr7zC1VdfzcKFCzn00EOpqanhtttuY8yYMQwbNozdu3fTu3dvrrzySi6//HIAZs2axbXXXssxxxwDwOjRo/nFL34R82WIyKzbYeiC8ixrYzOnJ2zcAPOmHThs9akw7fOZR4hZIHxfL10Lw72cczOBmQDjxo1rvXOLZ2+AnRuLjFiEnoNhzE3NPu2c45JLLuHyyy/nrrvuAqChoYFXX32VIUOGsHTpUgBeeuklLr30Uvbt28enPvUpAKZMmcIPfvCD7LKKSPs8NQBOrinPso6qhpd3+Ie/p0mGpwbAtOwjxCwQm4FBBY8HAluAbs0Mz8bOjQc3bnu8tb7Fpx9//HG6du3KFVdcsX9YbW0t69cfON173/tebr31Vr7yla/sLxAikjMTluL/DhvAVRPhpgdh9553h3Xvmgz35ro48wgxz4N4APikJcYDO5xzLwMLgRPM7Dgz6wZ8NB23Q1qxYgVjx44tatyTTz6ZVatW7X989913U1tbS21tLT/72c9CRRSRYq0q4xF5k0fBDR+Cbl2Sx0dXJ48njypbrpCHuf4KqAP6m9lm4EagK4Bz7sfAQySHuK4hOcz1U+lze83sC8DvSQ5z/alz7rlQOfOkafev2sUkkjPb+pR3eZNHwX2Lk/t3TGt+vEC5ghUI59zHWnneAVc289xDJAWkwxsxYgT33ntvUeMuXbqUk046KXAiEWmzi54A+kQO4XHRE3S2XUwV4W/+5m/4y1/+wr/927/tH7Zw4UI2bDjwCIX169dzzTXXcNVVV5U7oogUq74udgK/+rogs+30l9qIzcy4//77ufrqq/nOd75D9+7d9x/munbtWt73vvftP8z1qquu0g/UInk2+JXYCfwC5aq8AtFzcKtHHpU8v1YMGDCAe+6556Dhu3btanaaadOmMW3atPYkE5Gs9X89dgK/QLkqr0C0cM6CiEiLlgyDsxpipzjYkmFBZqvfIEREinVOTnusDJRLBUJEpFiBvqm3m7YgREQi29ErdgK/QLlUIEREinXRE7ET+AXKpQIhIlKs+rrYCfzq64LMVgVCRKRYQzbHTuAXKJcKhIhIsXrktIfEQLlUIMrk5ptvZsSIEYwePZra2lqeeeaZoMvr1au0H61mzJjBLbfc4n1u8+bN3H333QDU1NS0N5pIx7Xi+NgJ/ALlqrwT5SKYN28ev/3tb1myZAmHHnoo27Zt4+23344dq2iPPfYYzz//PFOmTIkdRSSu988HusdOcbD3z0cX6yuH2bOhpgaqqpK/GXQQ/vLLL9O/f38OPfRQAPr378+AAQMAuPjiixk7diwjRoxg5syZQHLhvmHDhvGZz3yGkSNHMnXqVObMmcOZZ57JCSecwIIFCw4Y7/LLL2f06NF85CMfYWfTPmyBO++8k1NPPZXa2lo+97nP8c477wDJVs3QoUM599xzWb16tTf73Llz+fKXv8y9995LbW0tb731FgBvvfUWkydPZsyYMYwcOXL/FoZIp/Z0bewEfoFyqUAUmj0bpk+HDRvAueTv9OntLhLnnXcemzZt4sQTT+Tzn/88Tz755P7nfvrTn7J48WIWLVrE97//fbZv3w7AmjVr+OIXv8iyZctYtWoVv/zlL5k7dy633HIL3/rWt/ZPv3r1aqZPn86yZcs47LDDuP322w9Y9sqVK7n77rt5+umnaWhooEuXLsyePZvFixdz1113sXTpUu677z4WLlzozT5hwgROOeUU6uvraWhoYOvWrQA88sgjDBgwgGeffZYVK1Zw/vnnt6uNRDqEPTnd6RIolwpEoeuvh6bfwHfuTIa3Q69evVi8eDEzZ87kiCOOYMqUKcyaNQuA73//+4wZM4bx48ezadMmXnzxRQCOO+44Ro0aRVVVFSNGjGDixImYGaNGjTqgu9JBgwZx5plnAnDZZZcxd+7cA5b92GOPsXjxYk455RRqa2t57LHHeOmll3jqqae45JJL6NmzJ4cddhgXXnhhs/lXr17N0KFDDxg2atQo5syZw3XXXcdTTz1FdXV1u9pIpEM4b37sBH6BcuW0HEaycWNpw0vQpUsX6urqqKurY9SoUfz85z+npqaGOXPmMG/ePHr27EldXR27dydHIzTujgKoqqra/7iqqoq9e/fuf87swP5xmz52znH55Zfz7W9/+4Dht91220Hj+mzfvp3q6mq6du16wPATTzyRxYsX89BDD/G1r32N8847jxtuuKGIlhDpwH43Ab7YEDvFwX43AS7Nfrbagig0uJlLdzc3vEirV6/ev2UA0NDQwLHHHsuOHTs4/PDD6dmzJ6tWrWL+/NK/BWzcuJF58+YB8Ktf/YoJEyYc8PzEiRO59957+dOf/gTAa6+9xoYNG/jrv/5r7r//fnbt2sUbb7zBgw8+6J3/unXr9v9eUmjLli307NmTyy67jGuuuYYlS5aUnF2kwxm2PnYCv0C5tAVR6Oabk98cCncz9eyZDG+HN998k6uuuorXX3+dQw45hOOPP56ZM2fSu3dvfvzjHzN69GiGDh3K+PHjS573SSedxM9//nM+97nPccIJJ/B3f/d3Bzw/fPhwvvnNb3Leeeexb98+unbtyg9/+EPGjx/PlClTqK2t5dhjj+Wss87yzn/YsGFs27aNkSNHMnPmTM444wwAli9fzrXXXktVVRVdu3blRz/6UekNIyK5ZknX0J3DuHHj3KJFB172duXKlaX18zx7dvKbw8aNyZbDzTfD1KkZJ83G+vXrueCCC1ixYkXsKM0quf1F8uzT/1H+XUyfnpX8vWNa8+P8Sy3ccXGbZm9mi51z43zPaQuiqalTc1sQRCSyyXOBHF7RdfJcdB6EHKCmpibXWw8inc6jpe8GLotAuVQgRESK1XVv6+PEECiXCoSISLHObIidwC9QLhUIEZFi/SGnu5gC5VKBEBEp1sg1sRP4BcqlAiEiUqxdObySKwTLpQIhIlKstQNjJ/ALlEsFQkSkWBc9ETuBX6BcKhAiIsWqr4udwK++LshsVSDKYP369fTo0YPa2lrgwC48Q8p6Obt27aK2tpZu3bqxbdu2zOYr0mFUvxk7gV+gXCoQZTJkyBAaGhqApI+Gclz9NOvl9OjRg4aGBu/VXUUqwsmrYifwC5RLBaLMmnbhuW7dug69HJGK8rj3mnbxBcqlAlFmTbvwPO644zr0ckQqSoVtQVTe1VyfSG8fTv9uB6YDM4FRJBdqnAd8DPgt8Bfg48AsYGw6j8XANOCXwBlAXWkRCrvwnDVrFv379+eCCy44aLy1a9fypS99ibPOOouTTjrJO06xy2m0b98+fvvb33L44Yd7+4BomqelfCIVZ1sf4JXYKQ62rU+Q2VZegajj3Q/0UQXDZxTc/0D6d2gzz38o/fv10hffXBeeGzZs4Hvf+x7OOYYMGcLVV1/Nj370I/bu3cvDDz/M0qVLAYr+oG66nFmzZjFnzhzGjRtHnz59qKqqYuXKlXzjG99g6NChPPPMMzzyyCMA3HPPPTzyyCMceeSRbNiwgZ1pB0oqElLxNh5FLgvExqOCzLbyCkRkzXXhefvtt9OjRw969OjB8uXLgeQDuaamhl69epX8Ld63nEmTJjF16lRmzZoFwE9+8hO++93vcswxx/CBD3xg/3gf+MAHmDp1KlOmTGHSpEnaghBpdNETQJ/IITwuegL1B9EJFHbh+cc//nH/8H379jF16lRmzJjBHXfcccA0VVUHvk0f/OAH2bJlS4vDfMuprq4+YBrnHGaGmR0wvHE8Mzto2SIVrb4udgK/+rogs9UWRJn16tWLBQsW7H/8wgsvAPCFL3yBr3/96xx99NH07t2bG2+8cf84Y8aM4eabb2bv3r1cfPHFPPTQQwfNt+mw5pZT6LOf/SzXXXcdJ554Ir16+XvJarpskYrW//XYCfwC5VKf1GWwadMmzjjjDPr167f/XIg8eO2117jtttvYvn07EydO5NJLL21x/F27dnH66aezdetWli9fTt++fVtdRh7aXyQzP7sRxlrr42WpmD6pFzv41DfaNHv1SR3ZoEGD2LRpU+wYB+nbty833XRT0eM3nignUrHmvg/GNsROcbC574NPZT/boDuYzex8M1ttZmvM7Kue5681s4b0tsLM3jGzvulz681sefrcooPnLiIVa/ZsqKmBqqrk7+zZ5VnuaTntAz5QrmAFwsy6AD8EJgHDgY+Z2fDCcZxz/885V+ucqwW+BjzpnHutYJRz0udzevqiiJTd7NkwfTps2ADOJX+nTy9PkQh0OGm7BcoVcgviVGCNc+4l59zbwF3ARS2M/zHgVwHziEhncP31kJ6bs9/Oncnw0F7uH34ZbREoV8jfII4BCne8bwZO841oZj2B84EvFAx2wKNm5oB/dc7NbGba6STnQjN48GBvkMbDOaW8OtMBEFJmc2+ALhv9z23c0PzwedOCRQLggm5ADi9WGeg8iJAFwveJ3NwnxoeAp5vsXjrTObfFzP4K+IOZrXLO/ddBM0wKx0xIjmJq+nz37t3Zvn07/fr1U5EoI+cc27dvp3v3nHbRKPn2s5Phi83s4DiqGl7e4R/+npqgsXi0Dv7uibDLaIv6Orgw+9mGLBCbgUEFjwcCW5oZ96M02b3knNuS/v2Tmd1PssvqoALRmoEDB7J582a2bt1a6qTSTt27d2fgwJx20Sj5dnQL/Y1cNRFuehB273l3WPeuyfDQBrwefhlt0VJ7tUPIArEQOMHMjgP+m6QIfLzpSGZWDZwNXFYw7D1AlXPujfT+eUDxx2MW6Nq1q65kKtLRDG7hekeT04uozaiHt9+Bo6uT4jB5VPPTZOXY7eGX0RYttVc7BPuR2jm3l+Q3hd8DK4F7nHPPmdkVZnZFwaiXAI86594qGHYkMNfMngUWAL9zzj0SKquI5MwzI1t+fvIoGD0Qxh0Lj1xdnuIAMG9IeZZTqtbaq42CnijnnHsIeKjJsB83eTyL5GLahcNeAsaEzCYiOTZhKf6fMSM7e3XsBH4TlqKL9YlIZViV093Cz+XwCCYI1l4qECKSP4E6wGm3bb1jJ/AL1F4qECKSPxc9ETuB34cXx07gF6i9VCBEJH/q62In8PvN2NbHiaG+LshsVSBEJH8CHbbZbjrMVUQksrx2zHPEG7ET+AVqLxUIEcmfJcNiJ/BbVBM7gV+g9lKBEJH8OSenXcCc+3zsBH6B2ksFQkTyJ69bEAtrYifw0xaEiFSMHb1iJ/Db0TN2Ar9A7aUCISL5o/MgSqPzIESkYtTXxU7gp/MgREQiG7I5dgK/E16NncAvUHupQIhI/vTYHTuBX4+3YyfwC9ReKhAikj8rjo+dwG/ZoNbHiSFQe6lAiEj+vH9+7AR+k5bHTuAXqL1UIEQkf56ujZ3A78kTYyfwC9ReKhAikj97gnZ22XYVlksFQkTy5zztYipJoPZSgRCR/PndhNgJ/B6ojZ3AL1B7qUCISP4MWx87gd9JW2In8AvUXioQIiLipQIhIvmzqiZ2Ar+VA2In8AvUXioQIpI/k+fGTuB3YUPsBH6B2ksFQkTy59HxsRP4PTwqdgK/QO2lAiEi+dN1b+wEfhWWSwVCRPLnzIbYCfzOfiF2Ar9A7aUCISL58wftYipJoPZSgRCR/Bm5JnYCv9GbYifwC9ReKhAikj+7usdO4LerW+wEfoHaSwVCRPJn7cDYCfxePDJ2Ar9A7aUCISL5c9ETsRP4fXhx7AR+gdpLBUJE8qe+LnYCv9+MjZ3Ar74uyGxVIEQkf6rfjJ3Ar3pn7AR+gdpLBUJE8ufkVbET+J2yPnYCv0DtpQIhIvnz+LjYCfzmDI+dwC9Qe6lAiEj+5HULYtz62An8tAUhIhVjW5/YCfy29o6dwC9Qe6lAiEj+bDwqdgK/Df1iJ/AL1F4qECKSPzoPojQ6D0JEKkZ9XewEfjoPQkQksv6vx07g1/+N2An8ArXXIa2NYGYLgWXA8sa/zrmtQdKIiAAMWwdY7BQHG7EldgK/YeuCzLaYLYiLgF8D3YArgPVmtiFIGhERgLnvi53A78mhsRP4BWqvVrcgnHNbgC3AIwBmdhLwkSBpREQATlsRO4Hf6WtjJ/A7bQVwceazbXULwswGFz52zq0ERmSeRESkkQ5zLU2g9mp1CwK428wGAetIfofYDQwLkkZEBODl/sDm2CkOtqVP7AR+L/cPMttidjGdDmBmxwOjgL7ArUHSiIhAelx/n8ghPHJ9HsTFmc+26MNcnXNrnHP3O+fucM7lsLSLSKdRXxc7gZ/OgxARiezobbET+A14PXYCv0DtpQIhIvkz+JXYCfyO3R47gV+g9lKBEJH8eWZk7AR+84bETuAXqL1UIEQkfyYsjZ3A7+zVsRP4BWovFQgRyZ9Vx8VO4PfcgNgJ/AK1lwqEiORPXjsM2qYOg0RE4lJ/EKVRfxAiUjHq62In8NN5ECIikekw19LoMFcRqRh57TDoiMrqMEgFQkTyZ0lOrwe6qCZ2Ar9A7aUCISL5c86i2An8zn0+dgK/QO2lAiEi+ZPXLYiFNbET+GkLQkQqxo5esRP47egZO4FfoPZSgRCR/NF5EKXReRAiUjHq62In8NN5ECIikQ3JaZ9kJ7waO4FfoPZSgRCR/OmxO3YCvx5vx07gF6i9VCBEJH9WHB87gd+yQbET+AVqLxUIEcmf98+PncBv0vLYCfwCtZcKhIjkz9O1sRP4PXli7AR+gdpLBUJE8mfPIbET+FVYLhUIEcmf87SLqSSB2ksFQkTy53cTYifwe6A2dgK/QO2lAiEi+TNsfewEfidtiZ3AL1B7qUCIiIiXCoSI5M+qmtgJ/FYOiJ3AL1B7qUCISP5Mnhs7gd+FDbET+AVqLxUIEcmfR8fHTuD38KjYCfwCtZcKhIjkT9e9sRP4VVguFQgRyZ8zG2In8Dv7hdgJ/AK1lwqEiOTPH7SLqSSB2ksFQkTyZ+Sa2An8Rm+KncAvUHupQIhI/uzqHjuB365usRP4BWovFQgRyZ+1A2Mn8HvxyNgJ/AK1lwqEiOTPRU/ETuD34cWxE/gFai8VCBHJn/q62An8fjM2dgK/+rogs1WBEJH8qX4zdgK/6p2xE/gFai8VCBHJn5NXxU7gd8r62An8ArWXCoSI5M/j42In8JszPHYCv0DtpQIhIvmT1y2IcetjJ/DTFoSIVIxtfWIn8NvaO3YCv0DtpQIhIvmz8ajYCfw29IudwC9Qe6lAiEj+6DyI0ug8CBGpGPV1sRP46TwIEZHI+r8eO4Ff/zdiJ/AL1F4qECKSP8PWxU7gN2JL7AR+gdpLBUJE8mfu+2In8HtyaOwEfoHaSwVCRPLntBWxE/idvjZ2Ar9A7aUCISL5o8NcS9MRD3M1s/PNbLWZrTGzr3qen2pmy9LbH81sTLHTikgn9nL/2An8tvSJncAvUHsFKxBm1gX4ITAJGA58zMyaXshkHXC2c2408E/AzBKmFZE8mj0bamqgqir5O3t26fPQeRCl6YDnQZwKrHHOveScexu4C7iocATn3B+dc/+TPpwPDCx2WhHJodmzYfp02LABnEv+Tp9eepGorwuRrv10HkRmjgEKe/jenA5rzqeBh9s4rYjkwfXXw84mfSbs3JkML8XR27LLlKUBr8dO4BeovQ4JMteEeYY574hm55AUiAltmHY6MB1g8ODBpacUqWSzboehC7Kb38YNzQ+fN634+Qw4BuiaRaJsHbs9dgK/wa8EmW3IArEZGFTweCBw0FkmZjYa+AkwyTm3vZRpAZxzM0l/uxg3bpy3iIhIM54aACfXZDe/o6rh5R3+4e8pYTnL62DCExmFytC8IVC7qfXxyu2ZkenX5GyF3MW0EDjBzI4zs27AR4EHCkcws8HAfcAnnHMvlDKtiGRgwtJs53fVROje5Jt/967J8FKcvTq7TFnKa66s38dUsALhnNsLfAH4PbASuMc595yZXWFmV6Sj3QD0A243swYzW9TStKGyilSsVcdlO7/Jo+CGD0G3Lsnjo6uTx5NHlTaf5wZkmysrec2V9fuYCrmLCefcQ8BDTYb9uOD+Z4DPFDutiGQsREczk0fBfenhoHdMa9s8tuW1Y5685uoTZLY6k1qkkul8g9LkNVcHPA9CRPKuvi52Ar+8nm+Q11z1dUFmqwIhUskCHR7Zbnk9nDSvuQK9jyoQIpUsrx3zHJHTjnnymksdBolI5pYMi53Ab1FN7AR+ec0V6H1UgRCpZOcsip3A79znYyfwy2uuQO+jCoRIJcvrFsTCmtgJ/PKaS1sQIpK5Hb1iJ/Db0TN2Ar/c5grzPqpAiFQynQdRmrzm0nkQIpK5+rrYCfzyer5BXnPV1wWZrQqESCUbsjl2Ar8TXo2dwC+vuQK9j6H7pL7ezN5Obw97np9kZm+YmTOzB0uZViTXsuh2M4t5tKbH7uznmYUeb8dO4JfbXGHex5B9UncFZgDvBw4H6szsQ01G20By1dan2zCtSD5l0e1mVl13tmbF8dnOLyvLBrU+Tgx5zRXofQy5BTEN2OGce9I59xbwJHBl4QjOueedcz8H9pY6rUhuZdHtZlZdd7bm/fOznV9WJi2PncAvr7kCvY8hL/c9FCjsKHU9cEbW05rZvwOXAvTt27fUjNLZ/NtMGPnHuBm83W6eDRufLL7bzay67mzNU+fA8HXZzS8rT54INfNipzhYXnM9XZt0q5axXPRJ3Z5pnXOfAD4B6nJUgMV94bSauBm83W6eAUc1FN/tZlZdd7amy0AghwViT9CuatquwnKF3MW0Cuhf8LgGeLkM00olOy8Hu0x83W52+3Vp3W5m1XVna/K6y0S5ShNovQ9ZIH4BVJvZWWb2HuBs4PYyTCuV7HcTYifwd7s54ZuldbuZVdedrXmgNtv5ZUW5ShNovQ+2veSc+4uZ/RPwGMkuo8edcw+Y2ez0+almNgpYCnQBMLO9wLHOuf/2TRsqq3Qiw9bHTpBo2u3mE++0fx4hnLQlzHzbS7lKE2i9D90n9U3ATU2GTS24v7y5DL5pRUSkfHQmtXQuq2piJ/BbOSB2Aj/lKk1ecwVa71UgpHOZPDd2Ar8LG2In8FOu0uQ1V6D1XgVCOpdHx8dO4Pdwxj8uZ0W5SpPXXIHWexUI6Vy6Nj0pPyeUqzTKVZpAuVQgpHM5syF2Ar+zX4idwE+5SpPXXIHWexUI6Vz+oF1MJVGu0uQ1V6D1XgVCOpeRa2In8Bu9KXYCP+UqTV5zBVrvVSCkc9nVPXYCv13dYifwU67S5DZXmPVeBUI6l7UDYyfwe/HI2An8lKs0ec0VaL1XgZDOJVDn7e2W187ulas0ec0VaL1XgZDOpb4udgK/vHZ2r1ylyWuu+rogsw3dJ/WstL9pZ2arPc+bmW1Nn99nZl8vdtpOpRx9D7dHR+lfGaD6zTDzba/qna2PE4NylSa3ucKs9yH7pD4UuBz4GEm/0ieY2eebjPYLoDrN8W3gn0qYtnMoV9/DbdWR+lcGOHlV9vPMwinrYyfwU67S5DVXoPU+5BbEDGCPc+4u59zrwBrg6ibjXAD8p0tcD1SZ2cQip+0cytX3cFt1pP6VAR4fl/08szBneOwEfspVmrzmCrTeh7zc9wig8FNhEzCmyTg9gRUFj/cCZxY5LQBmNh84BaBnz55tT7tsRtunbY+W+h5eNqOcSfyyyFfO13hqf3hrfbbzbKvv1yV/31oPI+e3LVfhPEJoa67WtDd3MblCt41PqPZqSTGvc3yYq8zmsU/qd0qZ1jm3/xTCdvVJPXpGmydtl8Gzkl0uBw0/Nl6mQlnkK+drHJ3t7DJzeuwAzVCu0lRYrpC7mFaQbCE0GgT8T5NxdgIjCx4fAswvctrO4eaboemWT8+eyfA8yCJf3l+jiPg554LcSD7gHfC3QB9gH3Blk3FmA2+TbDHcDLxT7LS+29ixY12HdOedzh17rHNmyd8774yd6EBZ5Mv7axSpUMAi18xnqiXPh2Fm/w5clj5c45w7wcyeSQvTaWZmwDagb1oQZrikq1HvtK0tb9y4cW7RokVZvwwRkU7LzBY757y/cgctEOWmAiEiUpqWCoTOpBYRES8VCBER8VKBEBERLxUIERHx6lQ/UpvZVqCZ03Zb1Z/kiKqYfBnykKtRYZa25gr9evLUXvBunvbkCvma8jrvYqct9/sda/1qbbntyXWsc+4I3xOdqkC0h5ktau6X/JgZ8pCrUWGWtuYK/Xry1F7wbp725Ar5mvI672KnLff7HWv9am25oXJpF5OIiHipQIiIiJcKxLtmxg6AP0MecjWa2cz9ts4jhDy1F7ybpz25Qr6mvM672GnL/X7HWr9aW26QXPoNQkREvLQFISIiXioQIiLipQIhIiJeKhAiIjlnZn8VY7kVVyDMrNrMvmNmq8xse3pbmQ7rU8Ych5nZt83s383s402eu71cOYrR1jYzs6PM7Edm9kMz62dmM8xsuZndY2ZHx8oVipn1MrPvmtlWM3vHzPaZ2R4z21JMptCvJ1Z7mdnDWeQysyVm9o9mNiRU1rbkCrDcvk1u/YAFZnZ4+rhsuSquQAD3kHRfWuec6+ec6weckw77dRlz/IykJ73fAB81s9+Y2RFm9h3gf8f+wDOz8wse3gdMJOnZ7/fAcIprs1nA88Am4HFgFzAZeAr4cQYx8/JeNpoNfJDkkMObgBtJ2uop4NQiMoV+PcHmb2YnN3MbC9RmlOtwkh4mHzezBWb2JTMb0J7cGeXK2jZgccFtEXAMsCS9X75czXU111lvwOq2PBcgR0OTx9enb/A3gGUFw48CrgP+UOZ2WlJw/3Xgm8CxwJeA/yiyPZcW3N/Y0uvvyO9lwTKfLVwusDD9WwWsai1T6NcTcv7AO8B/knwRaHrblUWuJuvkWcDtwCvpMqYHeD+jrF/ANcAjwKiCYeti5KrELYgNZvYPZnZk4wAzO9LMriP5plsuh5rZ/vZ3zt0M7AY+AlQXDH/FOfddYHAZszVVBfwZ2O2c+2egpsg2K1y/ftHCc22Vl/ey0VvAjjTTJ4DX0uFHAP2KyBT69YSc/0rgc865c5reaP0iciXncs495Zz7PMk36+8Cp7czfya5suCcuwX4DHCDmd1qZr1JumQue65KLBBTSP5ZnzSz/zGz/wGeIOkX+2/LmONB4G+aDFsOzCP5NgZE/cD7KzP7spl9heSDbn+bASMprs3qzawXgHPuHxsHmtnxwAsZZMzLe9noCqALye6lfwXGp5meItm6aC1T6NcTcv4zaP7z5KqMch20zjjn3nHOPeKc+1TJiVsXbf1yzm12zv0vkq2jPwA9Y+TSmdQ5YmaHA18FLgIavx28AjwAfNc591pz0wbIcmOTQbc757aa2VHA/3XOfbJcWUQqmZn1AIY451aUfdmVWCDMbBjJh/AxJJtuW4AHnHMry5ihG/BRYItzbo4lRzKdQbKpPtM5t6dcWYrRljYrx2vMw3vZQqZxJIX+BeB7xWSy5AidS4BBwF7gReBXzrkdGWfLvL3aM+82rl8TSH78X+Gce7S9+bPKldFyW1wPypWr4nYxpbtr7iI5gmgBsDC9/ysz+2oZo/yM5IieL5rZvwP/C3gGOBeYb2bfN7N/MbPrzOykMuY6iKfNelBcmzX3Gk8BfhIgV6z3sjHPgoJMp6a3rSTv6WOtZTKzvyc5uqs7SRv1IPmAmGdmdRnkC9Ze7Zl3sdOa2YKC+58FfgD0Bm4M8X7HWr9aWw/KmivUL/F5vZF8m+vqGd4NeLGMOZalfw8BXiXZd30d0AC8DFyW3r6aDvtqXtqM9Iik1trM9xrTx0bBkVod/b0sWO7SxkzpP+0R6fD3ACtay0TyG1RjG/UEnkjvD6bgiLA8tld75l3stBx4VFzT9l0e4P2Msn61th6UM9chvqLRye0DBnBw16RHp8+VS1W6C+Y9JCtBNfBp4GSSwyPvbBzRzG4FngO+U65wZras4OFgYJmZ7SH5cG/8faS1NvO9xteAQ0k+RNsrL+9lo8Yt8uEku2+3Ajjn3jIzisx0CMlBCoeSfDvGObfRzPLeXu2Zd7HTVqW/01VxcPvubWvwDHKF0NJ6ULZclVggribZ3H+Rd48MGgwcD3yhjDnuIDk2vgvJORC/Bv6K5ESYXzUZN8YH3pHAB0jOzTib5MicV0i2bt5rZo/Qepsd9BrN7CVgPMkmcntdTT7ey0bVJP/Qi4A96W61XcBxwAkk+4xb8hNgoZnNB/6a5PBNzOwI3j1ktj2uJlx7tWfexU5bTXLimAHOzI5yzr2SHiln7czfnlxZa209+Fa5clXqj9RVJPuHjyFZsTaTfGt/p8UJs88xAMA5t8WSM6WvAS4n2R1x0BvvnHukjNnuAH7mnJubPi5ss6uBaymizTyv8VySXVQLWpquhJy5eC9bybQV2OycW1vEtCOAk0h+eF1VhmyZtVd75t3OaXsCRzrn1rUjfua52rncFteDcuWqyAIByfkFFBwB4Jx7NUKGQ5xze9P7vYBhwDqSb5u5+cBrlEWbmVlfl/Hhunl4L5vkOQIYTXKs+mZgbVszddT2Ss91GQOsdM49n5dcbVHuXGbWxzn3ei5yhfqhJa83kuvCzCc51PIPwByS3SDzgZPLmGMasJ3kB6dJwEvAYyRbDleQ/BbxPpJvRjHb6wiSE3OeTdup6DYDzkzb+TngtHTal9LXeHpneS8L8gxPl72b5J92B8kuph0kP6q2mKkjtxfJCV390/ufSNfrn5D84HpVFrlIiu78tD1mAocXPLcgwPsZZf0iOax1Dslvkn1i5gr6D5PHG8kRQad5ho8Hni1jjuVAf5L9038GhqRv/OL0AybqB176YTcHWEPy+8cKkq2bWUB1MW1GcgjeKJLLIGwDJqTDTwae7izvZcFyG/9pTyPZ/P95OvyzJNcpajFTR24vkl0hjfcXAv3S+z1p5Yi1YnMBc4HzSS7Ydw1JIR2SPrc0wPsZZf1KPxsuILn443agnuR8oh7lzhXkBeb5RsuHZa4pY46GgvtbCt/4pv9QMT7w0g+7oen9jU0+7O4tps048LDElU2eW5JBxly8lwXLfJYDD8ssvLjc861l6sjtRXL45THp/ceB7un9LsBzWeTi4AtcnkNyAtn4LNonL+tXk/WmB8nlM+5Li8Uvy5mrEo9ietjMfkdy8bjGH4IHAZ8kuYJiuWw0s2+THMK2ysy+R7LP+nySI4X2c87NN7P3lDEbJN9WVqf3/wO43MymkHxr+3p6v7U2KzwR82tNnuuWQca8vJeN1gJHm9nTJEebbDGzM4AakoMNZrUyfUdury8Bj5rZb0jWkf9Mj3Q7i+SEySxymZlVu/RsYufc42b2YZJL5vdtZ/725Mra/iOynHO7SC7vfY+ZVQMXA9vKlasif6Q2s0m8e5p64w/BDzjnHipjhsOAK0n2Vf+A5JDSW9I83ybZzIR33/h1zrmyHbppZveRfCt8DLiUZNfHOmAgUAfcSSttZmYXAnOcczubDB8CfNg5938zyBn9vSzI0gf4OslhwV1JDgveB/yJ5ESuf25l+g7dXukH2MeBE0kOod8M1LsijsYqJpcll2p5yTk3v8m0g4H/45z7bHtfQ1tyBVjmNS65omv0XBVZIPIsLx94BR92w0l2nXzHOfdG+iFwUtN/UhHpfCruWkwtMbPpsTM45x4m2Qf5IefcBc65K2J8G3bOve6c+4c0w/XOuTfS4TsKi0Nb2yx0W+fhvSxkZtPbk6kjt1c5Xne53+9Y61dry806lwrEgUKcjdkWB+XI0wdekyxtbbPQbZ2X97KR0b5MHbm9yvG6y/1+x1q/Wltuprkq8Ufqli6V+685zhH1Ay/NegzJ1VgLszS9Hkyz0znn3ix2uhJynQo459xCMxtO8iP/qnK/l63kmeycm9jG+f3COffJEK/HDrxcdqbzL2XeZnYayVFbf7ak74Ovkhza+zzJZSVKGi/D1/D3wP3OuYM66wq5frW2XGvlMvqZZqm03yAsuVTux0iuBbQ5HTyQpMHvcs6V5YJ4peYws08551o7GiSIdIW9kmQFrCUpYn+fPrfEOXdykdN90TlX39p0JeS6keQkw0NIzhs5jaRnrXOB37ukG9eySfN8kWTLfCvJ8frbSc6Mf4Hk8MQLW5j+gaaDSA7l/E+AlqYtMt8C59yp6f3Pkrw39wPnAQ+2Z91vz7zN7DlgjHNur5nNBHYC9wIT0+GXljJeVsxsB0k3smtJro/2a5deIDCk1pZrZrNJ1vmeJP3F9yI5DHYigHNuWmZhsjxmtiPcyMklokvNQXqJ7Uhtthzold6vAf5C8mEPLZyg5JluUTHTlZirS/qP8mfgsHR4DzK4nHgb8ywh+ad+i6R41ZEctrwWOLuV6ZeQHB1WR3IkVOO0Z7c2bZH5lhbcz/Ry2e2ZNwXnfNDkfAYOPF+oqPEyfD+XkhT780guPLmV5DDSy4HeAdejFpdL4MvoF94q8TeIxkvlNlXuK6YelMOSS2w/Cwwys2UFt+W8e4ntsihcPsk34D+m9x9MR5lkyWXIW9r11cWlu5Wcc+tJPvCKma5Ye13SJ/FOkusd/Tld1i7iXO57L0kvcgtILtX8snPuCZLLbbzhnHuylenHkZxJfz2wo3Fa59yTRUxbjCozO9zM+sGBl8tOs8ea9woza+xT+lkzGwdgZicCe9owXlacc26fc+5R59ynSf5fbyfZbfhSgOUVu9zGy+j35t3L6EN2l9HfrxJ/g7iafFwi2pdjKMk3gulA4QeCAX8sYzY48HLfvwS+SbKvtzHLBcBPSc6PaM4rZlbrnGsAcM69aWbFTFest82sZ1ogxjYOTA/FjVEg3iY5g/ifzexe4FYze5Xkn7bVPM65fcA/m9mv07+vku3/aDXhLpfdnnl/BvgXM/tHkkuMzDOzTST/F59pw3hZOSC3S7rIfQB4IP0NJJTWlnsFYS+j/26QdNOkolhOLhHtyfFJ4Bbn3H95xv2lc+7jZcy2/3LfZjaQ5Nv6K02zmNmZzrmnm5nHAdM1ea7Z6UrIeKhz7i+e4f2Bo51zyz2TBePLY2aTSX4T+WmpedJpz3TOfT3DmL7lhLxcdtHzNrPewHtJT7JzzVydtNjx2svMTnTOvRBi3u1drgW+jP7+5VRigRARkdZV4m8QIiJSBBUIERHxUoEQEREvFQgREfFSgRAJzMwGWtJ/hkiHogIhEt5EkmsGiXQoOsxVJKD0onX1JNfMeQO4JMT5BiIhqECIBGZJ15vXOOdWxM4iUgrtYhIJbyiwutWxRHJGBUIkoPTidTvS6+mIdCgqECJhHUfSEZRIh6MCIRLWKqC/ma0wszNihxEphX6kFhERL21BiIiIlwqEiIh4qUCIiIiXCoSIiHipQIiIiJcKhIiIeKlAiIiIlwqEiIh4/X/CFUgP1yhtHwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "torch.manual_seed(42)\n",
    "\n",
    "cdf = torch.tensor([0, 0, 0, .1, .1, .2, .75, 1.0])\n",
    "t_left = torch.arange(len(cdf), dtype=cdf.dtype)\n",
    "sizes = torch.ones_like(t_left)\n",
    "t = t_left + torch.rand_like(sizes) * sizes\n",
    "\n",
    "plt.fill_between(torch.cat((t_left, (t_left + sizes)[-1:])),\n",
    "                 torch.cat((cdf, cdf[-1:])),\n",
    "                 alpha=.5, color=\"orange\", step=\"post\", label=\"CDF\")\n",
    "plt.scatter(t, cdf, color=\"red\", label=\"Sampled $t$'s\")\n",
    "plt.vlines(t, 0, cdf, color=\"red\")\n",
    "plt.vlines(t_left, 0, cdf,\n",
    "           color=\"magenta\", linestyles=\"dashed\", linewidth=.5,\n",
    "           label=\"$[t_{\\\\mathrm{left}}..t_{\\\\mathrm{right}}]$\")\n",
    "plt.vlines((t_left + sizes)[-1:], 0, cdf[-1:],\n",
    "           color=\"magenta\", linestyles=\"dashed\", linewidth=.5)\n",
    "#plt.step(t, cdf, color=\"red\", where=\"mid\")\n",
    "plt.title(\"CDF\")\n",
    "plt.xlabel(\"$t$\")\n",
    "plt.ylabel(\"$u$\")\n",
    "plt.yticks(cdf)\n",
    "plt.xticks(torch.sort(torch.cat((t_left, t, (t_left + sizes)[-1:]))).values, rotation=90)\n",
    "_ = plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first density we encounter is in bucket $[t_3..t_4]$, so the minimal value we'd want to sample would be the left bound of that bucket, that is $3$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Inverse transform implementation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def inverse_transform(t_left, sizes, cdf, min_divisor=1e-6):\n",
    "    \"\"\"Assuming sample dimension is the innermost, as required by :func:`searchsorted`\"\"\"\n",
    "    n_samples = t_left.shape[-1]\n",
    "    \n",
    "    def transform(u):\n",
    "        \"\"\"transform(u)\n",
    "        - :obj:`u` uniform sample of type :obj:`torch.Tensor`\"\"\"\n",
    "        \n",
    "        iright = torch.searchsorted(cdf, u, right=True).clamp(1, n_samples - 1)\n",
    "        ileft = (iright - 1).clamp(0, n_samples - 2)  # we must put $u$ _before_ found indices\n",
    "        \n",
    "        tleft = torch.gather(t_left, -1, ileft) + torch.gather(sizes, -1, ileft)\n",
    "        s = torch.gather(sizes, -1, iright)\n",
    "        qleft = torch.gather(cdf, -1, ileft)\n",
    "        qright = torch.gather(cdf, -1, iright)\n",
    "        \n",
    "        u = (u - qleft) / (qright - qleft).clamp(min_divisor)\n",
    "        t = tleft + u * s\n",
    "        return t\n",
    "    \n",
    "    return transform"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Why `right=True`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With `right=False` if we sample (e.g. deterministically, via `linspace`) $u=0$, the `searchsorted` would return the index $0$. After clamping, this corresponds to the interval `t_left[0]..sizes[0]` which is known not to contain any density (first interval that does have density is `t_left[2] + sizes[2]`)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tensor([0.0000, 0.0000, 0.0000, 0.1000, 0.1000, 0.2000, 0.7500, 1.0000])\n",
      "tensor([0, 3, 3, 6, 6, 7])\n"
     ]
    }
   ],
   "source": [
    "print(cdf)\n",
    "print(torch.searchsorted(cdf, torch.tensor([0, .05, .1, .7, .75, 1.0])))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tensor([0.0000, 0.0000, 0.0000, 0.1000, 0.1000, 0.2000, 0.7500, 1.0000])\n",
      "tensor([3, 3, 5, 6, 7, 8])\n"
     ]
    }
   ],
   "source": [
    "print(cdf)\n",
    "print(torch.searchsorted(cdf, torch.tensor([0, .05, .1, .7, .75, 1.0]), right=True))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Eval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "T = inverse_transform(t_left, sizes, cdf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAloAAAFjCAYAAADlxqUSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAyGklEQVR4nO3de7xcdX3o/c+XhHBJQYxcKpAYihEVudko0dhjvKTiFdqnFWliq6eWp68jtsjBlgpPQQvneI59UnrResBatSAgKvuxiiCtpiqSqJhAQKFc5JJAAQWEck/4Pn/MbNkMs/es2TNrrbl83q/XvLJnrd9av++smcz+7t/6rt+KzESSJEn9t13dAUiSJI0qEy1JkqSSmGhJkiSVxERLkiSpJCZakiRJJTHRkiRJKomJliRJUklMtCRJkkpioqWRExHXRsSK5s8HRMSGiHgwIv6ogr5viYjXl91PWWYbf0R8OiJOLyOmKgx7/GVqPTZT/39VHMf/jIjju2j/vYg4sMSQpEJMtDRwIiIj4vkty06LiHOKbJ+ZB2bm2ubTPwHWZuYumfk3fQ5VfdBM7hbXHYeKafn/1RcR8XsR8Z/Nx6MRsW3K8/sjYi/gd4H/08Vu/xL48JQ+/JypFiZaGnXPA66dzYYRMbfPsWhMtfss+fl6SmZ+JjN/KTN/CfgfwFcmn2fmbjSSrIsz85Eudvtl4DUR8dwSQpYKM9HS0Gn+ZXpiRFwdET+PiAsiYseW9a+PiG8ArwH+rvmX8Qsi4kURsbb5V/K1EfG2Nvv+04i4GngoIuY2l32g2d9DEfEPEbFXRHyteUryXyLi2W3i/EBEfLFl2d9GxJnTvK4/jYgtzX1eHxGvay4/KSJuai7/UUT8Rku8hWNrtv+z5n7ui4h/nHrsWuLZOyK+GBH3RMRPpp56jYjDIuKHzT4uANruY0r7M6a+7ojYtxnvM76DpjsObdp1Oi4zfUYKxx8RCyPiS83j8LOI+Lsp66b9PM3wWWpdNtNxnvZYzPReFohr1scmppxeLrCvl8ZTp+4vbK7vdIr2UOCqlmVvBP6tJY4ZP1OZ+ShwJfDr7Top+jmTepaZPnwM1ANI4Pkty04Dzmn+fAvwPWBvYAHwY+APp7S9BXh98+e1wHuaP28P3Ah8EJgHvBZ4EDigZduNwEJgpynL1gF7AfsAdwM/BA4DdgC+AZzapu/nAg8BuzWfz21u+6ttXvMBwO3A3s3ni4H9mz//dvO1bgcc3dznc7uNbUr7a5qvbwFwOXB667Fr9nUl8OfNY/UrwM3AG5rPbwXe3zymvwU8MXU/bV7fV4F3T3n+ZuB73RyHNm07HZe2n5Fu4gfm0Pil/1fAfBpJx6uKfJ6Y/rP0i2UdjvOMx2K697JgXLM+Njz9M15kX3/c3NdvAo+3O84tx/xm4Ddblt0DvKzbzxTwN8CaXj5nPnz0+nBES8PqbzLzjsy8F/hnGn8Fd7IM+CXgI5n5eGZ+A/gKcEybfd+eTz9N8beZeVdmbgG+DazPzA2Z+RhwEY3E5mky807gWzQSAoAjgJ9m5pVtYttGIzF6cURsn5m3ZOZNzf1c2HytT2bmBcANwMt7iO3vmq/vXuCMNq8f4GXAHpn54eaxuhk4G3hH8zhuD5yZmU9k5heA77fZx1SHAldPeX5Iy/OOx6FVgeMy3Wekm/hfTiOJ+EBmPpSZj2bmd6bsp9Pnqd1naeqymY5zkWPR7r0sGlevx6bIvuY21z+RmV+ikZRNKyJ2pZH0bGxZtRuNZHGqQ+n8mXqwuW2rwp8zqVcmWhpE22h82U+1PY2/rCf9x5SfH6bxi6WTvYHbM/PJKctupTESNNXtbba9a8rPj7R5Pl3/nwFWN39eDfxTu0aZeSNwPI2Ru7sj4vyI2BsgIn43IjY2TwPdD7wE2L2H2Ka+vltpHJdWzwP2nuyz2e8HaYyc7Q1sycxs2U9bEbFHc7uptXKH8MzTQzMehzb77XRcpvuMdBP/QuDWzNzaZl2Rz1O7z9LUZdMe54LHot17WSSufhyb2eyr3fGY6hAaydFPWpbfB+wy+aSLz9QuwP2tnXTzOZN6ZaKlQXQbjb9qp9qPzl/4ndwBLGypC1oEbGlpl/TPBHBwRLwEeAtw7nQNM/NzmfkqGr98E/hfEfE8GiMcxwHPyUZh8DVA9BDTwik/L6JxXFrdDvwkM3eb8tglM98E3AnsExHRsp/pHAjckI2amcki8NfQfkSr7XFobdPjcekm/tuBRdG+cL3I56ndZ6k18ZjuOBc5Fu3ey6Kf83a6fW+73dfC6Ro3HQpc3ZKcQeOz8oIpz4t+pl5Em4Qein3OpH4w0dIgugA4pVncul2z8PatwBd63O96GnU8fxIR20djLqC3Auf3uN9pNX8RfAH4HI36kdvatYvGfF+vjYgdgEdpjERto1EXlDRqVIiId9MYuenFe5vHdgGN0ZML2rT5HvBAs2B4p4iYExEviYiXAVcAW4E/ahZz/yZPP2X3jJcH7Nxsux3wv4E9aJNozXAcWvVyXLqJ/3s0EoaPRMT8iNgxIpY31/Xj8zTtcS54LNq9l73E1e1722lf24Djmvs6ssC+DuWZpw0BLgZePeV5x89U87j9KnBZ6866+JxJPTPR0iD6MPBd4Ds0Thn8b2BVZl7Ty04z83HgbTSuYPop8HHgdzPzut7C7egzwEFMc9qwaQfgI824/gPYE/hgZv4I+H9p/NK6q7mfy3uM53PA12kUHd9Mo4D6aTJzG41fzofSOI3zU+CTwLOax/E3gXfReH+OBr40Q3/fpvEL8Doav/RuAzZn5n1t2rY9Dm3im/Vx6Sb+Kcfh+ZNxN9v35fM003Gm2LF4xnvZS1yzeG+L7Ov3aZy+W02jVuyxGTY7hPaJ1meBN0XETs3nRT5Tb6Mxh167EdtCnzOpH+KZI7SS+ikiFtH4hfDLmflAzbHcQuMqzH+pMw71bhjfy4hYD3wiM/9xFtv+D+DuzDyzi75+v9c/0KReOWGeVKLmaY0TgPPrTrKkqkXEq4HraYwcrQIOBi6Zzb4ys6sRp8w8fDb9SP1moiWVJCLm0zitdSuNqR2kcXMA8HkaVyLeBPxWNqY9kcaGpw4lSZJKYjG8JElSSUy0JEmSSmKiJUmSVBITLUmSpJKYaEmSJJXEREuSJKkkJlqSJEklMdGSJEkqiYmWJElSSUy0JEmSSmKiJUmSVBITLUmSpJKYaEmSJJXEREuSJKkkJlqSJEklMdGSJEkqiYmWJElSSUy0JEmSSmKiJUmSVBITLUmSpJLMrTuA6ey+++65ePHiusOQVJErr7zyp5m5R91x9IPfX9L4me47bGATrcWLF/ODH/yg7jAkVSQibq07hn7x+0saP9N9h3nqUJIkqSQmWpIkSSUx0ZIkSSqJiZYkSVJJTLQkSZJKYqIlSZJUEhMtSZKkklQ2j1ZEvB94D5DAJuDdmfloVf1LqsbEhi2c9uVruf+RJwB49s7bc+pbD+Sow/apLaaI+BTwFuDuzHxJm/UB/DXwJuBh4F2Z+cPmuiOa6+YAn8zMj1QWuKTKTWzYwkcvvZ477n+EvXfbiQ+84YCevr8qGdGKiH2APwKWNr/k5gDvqKJvSdU5/IzLOP6Cjb9IsgDue/gJPvCFq5jYsKXGyPg0cMQM698ILGk+jgX+HiAi5gAfa65/MXBMRLy41Egl1WZiwxbef8FGttz/CAlsuf8R3n/Bxp6+v6o8dTgX2Cki5gI7A3dU2LekEq1cs5bFJ32Vux58vO36J7YlH730+oqjekpmfgu4d4YmRwKfzYZ1wG4R8Vzg5cCNmXlzZj4OnN9sK2kEHX/BRrJlWQInXnjVrPdZSaKVmVuAvwRuA+4Efp6ZX29tFxHHRsQPIuIH99xzTxWhSerBqrOvYPFJX+WGux/q2PaO+x+pIKJZ2we4fcrzzc1l0y1/hll/f63tNtQ+sE/7HKb+Kupz5ZrpO9n6ZGv6VVxVpw6fTeOvwP2AvYH5EbG6tV1mnpWZSzNz6R57jMS9ZaWRdMrEJhaf9FUuv2mmQaKn23u3nUqMqGfRZlnOsPyZC2f7/bW2eNO+sU/7HKb+KuqzyB+Ms1FVMfzrgZ9k5j0AEfEl4JXAORX1L6kPTpnYxDnrbpvVth94wwF9jqavNgMLpzzfl0Z5w7xplksaIavOvqK0fVeVaN0GLIuInYFHgNcB3tpeGhITG7Zw/AUbZ739XrvMq/WqwwK+DBwXEecDh9Mob7gzIu4BlkTEfsAWGhfx/E6NcUoqQafR+SV7zp/1vitJtDJzfUR8AfghsBXYAJxVRd+SenPwqZfwwGPbZr39kj3nc9kJK/oX0CxExHnACmD3iNgMnApsD5CZnwAupjG1w400pnd4d3Pd1og4DriUxtXSn8rMayt/AZJKc8rEpo5tevkOq2wercw8lcaXm6QhsHLN2p5qFvbaZR7rT17Zx4hmLzOP6bA+gfdOs+5iGomYpBHUqRyi12L2yhItScNh1dlXdFXk3mq7gDVvP3TQTxVKUqHRrDVHH9pTHyZakoDeR7AAzjzaBEvS8ChycU+v32kmWtKY63UEC2D1skWcftRBfYpIkspXZDRr9bJFPfdjoiWNqV6mapi0fP8FnPsHr+hTRJJUnSLff/34A9JESxozJliSxl2Rexf2YzQLTLSksdHrXFgwGFM1SFKvTvj8xo5t+lUOYaIljYFe58LadYc5XP2hI/oYkSTVY2LDFjrdunD5/gv61p+JljTCer2S0KkaJI2aEy+8qmObfpZGmGhJI8i5sCSpva0dhrP6OZoFJlrSSOnHVA3OhSVpVK1cs7Zjm35f6GOiJY2AflxJ6FxYkkZdp1KKfo9mgYmWNNT6cSWhUzVIGgerzr6iY5syvgtNtKQhZIIlSd3pVFaxZM/5pfRroiUNEefCkqTuFRnNKut70URLGhLOhSVJs9NpNGu7Evs20ZIGXK9zYQXwV15JKGlMFbl59JqjDy2tfxMtaUD1Y6oGrySUNO6KXJFd5h+iJlrSgDHBkqT+KHLz6DNLHM0CEy1pYPRjLiyvJJSkpxS5eXTZZRUmWlLNvJJQkvqvyM2jVy9bVHocJlpSTZwLS5LKU+Tm0VWUWFSSaEXEAcAFUxb9CvDnmXlmFf1Lg+bwMy7jrgcfn/X2jmBJ0syqvnn0dCpJtDLzeuBQgIiYA2wBLqqib2lQ9KMGy7mwJKmzOm4ePZ06Th2+DrgpM2+toW+pMv24enCSc2FJUnF13Dx6OnUkWu8AzquhX6lU/UyspnKqBkkqrq6bR0+n0kQrIuYBbwP+bJr1xwLHAixa1MWVAGubD6lCX7zydm6/75FfPH8Zz+FlPKdv+z9432fx2hfuBRtpPIbRiuZDkipS182jp1P1iNYbgR9m5l3tVmbmWcBZAEuXLu1wUeYUK/DLXKV7xq1wDiynH68klKTZKVKbVfWFRFUnWsfgaUMNiV6vDOyWVxJKUm861WaVefPo6VSWaEXEzsBK4P+uqk+pqIkNWzjhgo08WUPfJliS1LsitVll3jx6OpUlWpn5MPSxgEXqQZ2J1SSL3CWpfzrVZu26w5xartx2ZniNhX7Mwt6L7QLWvN3pGSSpDKdMbOrYpq45CE20NJL6MTlorxyxkqRqdPq+r6M2a5KJlkZCWXNYFTV3u+Avf/sQR6wkqWJFRrPqqM2aZKKloVR3YrXXLvNYf/LK2vqXJDUUOXtR5x/BJloaCnUnVl4ZKEmDZ2LDlo5tVi/rYgL0EphoaSA9Y3LQijlpqCQNvhMvvKpjm7prZU20NBDqTqwcsZKk4bP1yZlvIlP3aBaYaKkGgzCHlSNWkjTcitxup+7RLDDRUkWqvp1NK6dakKTR0uksyPL9F1QUycxMtNR3gzBiZWIlSaOryGjWoJy1MNFSz+qeHNQ5rCRpvHQazVqy5/yKIunMREtdqzux8nY2kjS+itw8epAubjLRUkd1z2HliJUkaVKn30eDNJoFJlpqo+7EyqkWJEntFLndzqD9/jDRkomVJGkoDPLNo6djojWGnBxUkjRsBv3m0dMx0RoDdc9h5eSgkqReDfrNo6djojViBmEOKxMrSVI/XXfnA/DsmdsMwu122jHRGnKDkFg5OagkqUyXXvsf8KqZ2wzq7yETrSHjHFaSpHEysWELM986enBHs8BEa+DVnViBI1aSpPqc8PmN/BFLZmwzyL+jTLQGTN1TLTg5qCRpUExs2MKTHYazBuXm0dOpLNGKiN2ATwIvARL4r5nZeR79EVf3VAt77TKP9SevrK1/qQoRcQTw18Ac4JOZ+ZGW9R8AVjWfzgVeBOyRmfdGxC3Ag8A2YGtmLq0scGnMnXjhVR3bDPrFV1WOaP01cElm/lZEzAN2rrDvgVLnqJVzWGncRMQc4GPASmAz8P2I+HJm/miyTWZ+FPhos/1bgfdn5tT/pK/JzJ9WGLYkYGuH4axBu91OO5UkWhGxK/BfgHcBZObjQH0TO9Wo6hEsp1qQeDlwY2beDBAR5wNHAj+apv0xwHkVxSZpGivXrO3YZhgGDqoa0foV4B7gHyPiEOBK4I8zs75zZjWY2LCl9CTLxEp6hn2A26c83wwc3q5hROwMHAEcN2VxAl+PiAT+T2aeVVagkp7S6ffloNdmTaoq0ZoLvBR4X2auj4i/Bk4C/p+pjSLiWOBYgEWLurhUc23zMeBu+ZeHOL7DlRPdWvjsnfi/fnXhUwu2AKf1tQtp9lY0H/WKNsumOx/xVuDyltOGyzPzjojYE7gsIq7LzG89o5PZfn9JmpVhGVSoKtHaDGzOzPXN51+gkWg9TfMvxbMAli5d2mnajKesYBC+zGf0wpMv5tFXFX9J03GqBalrm4Epf42wL3DHNG3fQctpw8y8o/nv3RFxEY1Tkc9ItGb9/SWpa8NQmzWpkkQrM/8jIm6PiAMy83rgdUxfHzFyVq5Zy6Pbuv/edXJQqS++DyyJiP1ojPm+A/id1kYR8Szg1cDqKcvmA9tl5oPNn38d+HAlUUtjrFN91jDUZk2q8qrD9wHnNq84vBl4d4V916abuiznsJL6LzO3RsRxwKU0pnf4VGZeGxF/2Fz/iWbT3wC+3lI7uhdwUURA4/vyc5l5SXXRS+OpzmmP+q2yRCszNwJjN//M8Rds7NjGuaykcmXmxcDFLcs+0fL808CnW5bdDBxScniSRth2dQcwyl548sUd2+w4J0yyJEkqaJjqs8BEqzRF67KuO+NNFUQjSdJwmNiwZcb1w1SfBSZapShalzXIdxuXJKkORW67M0xMtErw/oJ1WU7TIEnS03W67c6wMdHqs4NPvWTamRAnWZclSdJ4MNHqo5Vr1vLAY9s6trMuS5Kk7g3LbXemMtHqE+uyJEnqzSkTm2ZcPyy33ZnKRKtPrMuSJKk356y7re4Q+s5Eqw+K1GUB1mVJkjRmTLR6VLQu68yjDy0/GEmSNFBMtHpQtC5r+f4LvH+hJEk9GNYaZxOtHhStyxrG4j1Jkqq0cs3aGdcPa42zidYsWZclSVL/FDlDNIxMtGbBuixJklSEiVaXrMuSJKlaUXcAPTDR6pJ1WZIk9dfEhi0zrn/Dgb9cUST9Z6LVBeuyJEnqvxMvvGrG9S987q4VRdJ/JloFWZclSVI5tj5ZZBhjOJloFWBdliRJmg0TrQKsy5IkqR7L919Qdwg9MdHq4PAzLrMuS5KkkpwysWnG9cM+iGGiNYNVZ1/BXQ8+3rGddVmSJM3OOetuqzuEUs2tqqOIuAV4ENgGbM3MpVX1PRsTG7Zw+U33dmxnXZYkSZpOZYlW02sy86cV9zkrJ3x+Y8c21mVJkqSZeOqwjVMmNlHkSlPrsiRJKs/qZYvqDqFnVSZaCXw9Iq6MiGMr7Ldr562/vWMb67IkSerNqrOvmHH96UcdVFEk5any1OHyzLwjIvYELouI6zLzW1MbNBOwYwEWLeoii13bfPTJ+779/BnXL5g/j6N23Af+v/71KY2kFc2HJLVRpBZ62FWWaGXmHc1/746Ii4CXA99qaXMWcBbA0qVLi08Tu4K+fpn/7WM3si3bd7/rDnO4+kNH9K8zSZI0sio5dRgR8yNil8mfgV8Hrqmi79k45vCFbZfvtcs8kyxJkiowKkXkVY1o7QVcFBGTfX4uMy+pqO+uTZ4TPm/97WzLZE4Exxy+cCTOFUuSNAzWjEgtdCWJVmbeDBxSRV/9cvpRB5lYSZJUkpVr1s64flTmqByVkTlJkjREbrj7obpDqISJliRJUklMtCRJ0kBZsuf8ukPoGxMtSZJUqYkNW2Zcf9kJK6oJpAImWpIkqVInXnhV3SFUxkRLkiRVamuRGwqPCBMtSZKkkphoSZKkgTFKhfBgoiVJkip0ysSmGdePUiE8mGhJkqQKnbPutrpDqJSJliRJUklMtCRJkkpioiVpaETE/IiYU3ccksqxetmiukPoOxMtSQMrIraLiN+JiK9GxN3AdcCdEXFtRHw0IpbUHaOk4ladfcWM608/6qCKIqmOiZakQfZNYH/gz4BfzsyFmbkn8GvAOuAjEbG6zgAlFXf5TffWHULl5tYdgCTN4PWZ+UREPC8zn5xcmJn3Al8EvhgR29cXniTNzBEtSQMrM59o/nhR67qIWNbSRtIQG9WEpOvXZTGqpKpExNsj4iPALhHxopbvnrPqiktS/605+tC6QyhFx1OHEbEd8A5gFfAy4DFgh4i4B7gYOCszbyg1Sknj6nJgR+A9wBrggIi4H7gDeKTGuCR1aeWatTOuP+qwfaoJpGJFarS+CfwLjWLUaybrJCJiAfAaGsWoF2XmOeWFKWkcZeYW4LMRcVNmXg6/+O7Zj8YViJKGxA13P1R3CLUokmi9vl0NhMWoksoWEZENl08ua3733NvappYAJamDjjVak0lWRHw/Iv4hIo6PiNdGxB6tbSSpz74ZEe+LiKfNYhgR85rfQ58Bfq+m2CT1yZI959cdQmm6KYY/ErgQmAf8IXBLRNzaTWcRMSciNkTEV7rZTtLYOgLYBpwXEXdExI8i4mbgBuAY4K8y89N1Biips4kNW2Zcf9kJK6oJpAaF59HKzDtoFKBeAhARLwJ+q8v+/hj4MbBrl9tV7pSJTZy3/na2ZTIngmMOXziSM9ZKgywzHwU+Dny8WaKwO/BIZt5fa2CSunLihVfVHUJtCo9otQ7dZ+aPgQO72H5f4M3AJwtHV5PDz7iMc9bdxrZm2ce2TM5ZdxunTGyqOTJpfGXmE5l5p0mWNHy2Pjm+ZZTdnDq8ICI2R8S3I+LjEbEGeGEX258J/AnwZId2tVq5Zi13Pfh423Xnrb+94mgkTYqIP607BknqVjenDl8BEBHPBw4CFtCY16ajiHgLcHdmXhkRK2ZodyxwLMCiRV3cwXtt89Gj6+58gDdfu/fMjU7rvR9pbKxoPmYhIj4/9SlwKPC/eoxI0oBZvv+CukMoVZEJS5926XRm3gjcOFObNpYDb4uIN9GYfHDXiDgnM592M9jMPIvmbM9Lly4tPs64gll/mU/1xpO+Tb5q+vVzIjj+tBf03pGkIh7IzPdMPomIv68zGEmzs+rsK2Zcf+4fvKKiSOpR5NRhz5dXZ+afZea+mbmYxizz32hNsup28KmX0CmzO+bwhZXEIgmAM1qen1xLFJJ6cvlN93ZuNMKKnDo8AvivNC6v3g+4H9iJRpL2dRqXV28sK8AqrFyzlgce29axnVcdStXJzJ8ARMTumfnT5kSlkjRUiiRap2fmifTp8urMXEtfKqr6Y2LDlkK3BVi9rIuaMUn99CngbXUHIUmzUSTReu3kD80Z4O8sL5zqvf+CjR3b7LXLPEezpPpE3QFIKseZRx9adwil62Z6h5FTpC4LYP3JK0uPRdK0xncCHmnIrVyzdsb1Rx22TzWB1KjIiNYhEfETYBNwzZR/rxvmexwWrcsah2xbGnCOaElDqkhpzqgrMqJ1NY3pGf4O+Bnw68A/Aj+NiGtKjK00Reuylu+/YCyybWlQRcRfAn9WdxySNFuFJiydcp/Dr08ui4gAnl9SXKUqWpc16nN7SEPgtc2LcSSNmCV7zq87hEoUSbQ+1m5hc4LSG/obTvmsy5IkqXwTG7bMuP6yE1ZUE0jNOiZamTnwN4EuyrosaeiMZI2oNA5OvPCqukMYCGNz1aF1WdJQGrkaUWlcbH3SC4ahi5tKDzvrsqThNGo1opLGy1gkWtZlSUNrpGpEJTUs339B3SFUZuRPHVqXJQ2vUaoRlcbJKRObZlw/TmePRjrRsi5LkqTqnbPutrpDGBgjnWhZlyVJkuo0sonW4WdcZl2WJEmq1UgmWqvOvoK7Hny8YzvrsiRJqtbqZYvqDqFSI5doTWzYwuU33duxnXVZkiT136qzr5hx/elHHVRRJINh5BIt67IkSapPkcGOcTJSiZZ1WZIkaZCMTKJlXZYkSYNtZJKOLozEay5al7Vkz/nWZUmSVJKJDVtmXL9mDAc7RiLRKlKXtesOc7jshBWlxyJJ0rg68cKrZlw/joMdQ59oFa3LuvpDR5QeiyRJ42zrk0V+I4+XShKtiNgxIr4XEVdFxLUR8aF+7Ne6LEmSNMjmVtTPY8BrM/M/I2J74DsR8bXMXDfbHVqXJUnS8Fi+/4K6Q6hFJYlWZibwn82n2zcfPY0vfuifr+3YxrosSZKqccrEphnXj+v8lZXVaEXEnIjYCNwNXJaZ63vZ330PP9GxjXVZkiRV45x1t9UdwkCq6tQhmbkNODQidgMuioiXZOY1U9tExLHAsQCLFhW/F9Ky2xaw7LbnPG3ZEQf+MpzWY9CSerOi+ZCkMVVZojUpM++PiLXAEcA1LevOAs4CWLp06YynFnfbaXvuf6QxqrVu0b2sW/RUvdaSPedz/Akv6G/gkiRJXarqqsM9miNZRMROwOuB63rZ52lvO5Dtt4tnLF+y53zrsiRJGiCrlxU/SzVqqhrRei7wmYiYQyO5+3xmfqWXHU5eSfjRS6/njvsfYe/dduIDbzjAKwwlSarYqrOvmHH96UcdVFEkg6eqqw6vBg7r936POmwfEytJkmpWZLqlcTX0M8NLkiQNKhMtSZJUmnFPNMb99UuSpBKtGfPb4JloSZKkWVu5Zu2M68e9ltpES5IkzdoNdz9UdwgDzURLkiSpJCZakiSpFEv2nF93CLUz0ZIkSbMysWHLjOu9U4uJliRJmqUTL7yq7hAGnomWJEmala1PZt0hDDwTLUmSpJKYaEmSpL5bvv+CukMYCCZakiSpa6dMbJpx/bl/8IqKIhlsJlqSJKlr56y7re4QhoKJliRJUklMtCRJkkpioiVJkvpq9bJFdYcwMEy0JElSV1adfcWM608/6qCKIhl8JlqSJKkrl990b90hDA0TLUmSpJKYaEmSpL4xsXg6j4ckSeqbNUcfWncIA6WSRCsiFkbENyPixxFxbUT8cRX9SpKk/vrsFbfMuP6ow/apJpAhMbeifrYC/z0zfxgRuwBXRsRlmfmjivqXJEl9cO9Dj9cdwlCpZEQrM+/MzB82f34Q+DFgyitJkkZa5TVaEbEYOAxYX3XfkiSpPEv2nF93CAOnqlOHAETELwFfBI7PzAfarD8WOBZg0aIuZpVd23xIGiwrmg9JI2Fiw5YZ1192wopqAhkilSVaEbE9jSTr3Mz8Urs2mXkWcBbA0qVLs/DOV+CXuSRJJTvxwqs4jufXHcZQqeqqwwD+AfhxZq6pok9JktRfW58sPgaihqpqtJYD7wReGxEbm483VdS3JElSLSo5dZiZ3wGiir4kSVL1lu+/oO4QBpIzw0uSpI5Omdg04/pz/+AVFUUyXEy0JElSR+esu63uEIaSiZYkSVJJTLQkjbyIOCIiro+IGyPipDbrV0XE1c3HdyPikKLbStJMTLQkjbSImAN8DHgj8GLgmIh4cUuznwCvzsyDgb+gOZ9fwW2lsbd6WReTjI8ZEy1Jo+7lwI2ZeXNmPg6cDxw5tUFmfjcz72s+XQfsW3RbaRy88OSLZ1x/+lEHVRTJ8DHRkjTq9gFun/J8MzPf1P73ga/Ncltp5Lzw5It5dJsTlc5Wpfc6lKQatJvDr+1vjYh4DY1E61Wz2HZ292qVBtTEhi0cf8HGusMYeo5oSRp1m4GFU57vC9zR2igiDgY+CRyZmT/rZlto3Ks1M5dm5tI99tijL4FLdTllYlPhJGvJnvPLDWbIOaIladR9H1gSEfsBW4B3AL8ztUFELAK+BLwzM/+9m22lUbPq7Cu4/KZ7C7e/7IQV5QUzAky0JI20zNwaEccBlwJzgE9l5rUR8YfN9Z8A/hx4DvDxiADY2hydarttLS9EqsDKNWu54e6HCrc/8+hDywtmRJhoSRp5mXkxcHHLsk9M+fk9wHuKbiuNosPPuIy7Hny8cPszjz6Uow7z2pBOTLQkSRpzB596CQ88tq1Q27kR3PKRN5cc0egw0ZIkaUx1e2XhrjvM4bjXLSkvoBHkVYeSJI2hbq4sBNhrl3lc/aEjygtoRDmiJUnSmOn2ysIle8736sJZMtGSJGmMdHtl4fL9F3DuH7yixIhGm4mWJEljotsrC1cvW+R9DHtkoiVJ0hjo5srCAP7K6Rv6wkRLkqQR182NoXfdYY5F733kVYeSJI2oiQ1bWHzSVwsnWV5Z2H+OaEmSNIJOmdjEOetuK9zeKwvLUcmIVkR8KiLujohrquhPkqRxtursK7pKspbvv8AkqyRVnTr8NOBYpCRJJVu5Zm1Xc2StXrbI6RtKVMmpw8z8VkQsrqIvSZLGlTeGHjzWaEmSNOS6vWfhjnOC6854U3kB6RcGKtGKiGOBYwEWLVpUfMO1zYekwbKi+ZBUmi9eeTv//dGrC7d3+oZqDVSilZlnAWcBLF26tNi1qOCXuSRpLK1cs5Y337d34fZ77TKP9SevLDEitXIeLUkjLyJOjojHm4+vtVn/xoh4MCIyIv65m22luhx+xmVd3bNwyZ7zTbJqUNX0DucBVwAHRMTmiPj9KvqVpIjYHjgNWAk8G1gREW9taXYrcBxw+Sy21ZhbdfYVLD7pqyw+6auc+S//zqqzr6ikz26K3p2+oT6VJFqZeUxmPjczt8/MfTPzH6roV5KAdwE/z8x/y8yHgH8D3ju1QWb+KDM/A2ztdluNt4NPveQZUylcftO9pSdbTt8wPAaqRkuSSnAA8NMpz28BXlnBthpxM90/sJtEqCzbBax5u9M31M1ES9KoizbLil5sU3jbiPgn4DcBFixYUHD3GkbdTqVQB2+nMzgshpc06q4Ddp/yfDFwZ7+3zcx3Zub8zJy/3377zSJMDYNTJjYNRJK1/Qy/va3HGiwmWpJG3WeBZ0XEr0XEfODVwMcr2FYjppv7B86UCPXDR3/70LbLrccaPJ46lDTSMvOxiPgL4F9pnAr8ZmZ+OSLOba5fFREHARuAOQARsRV4XmZuabdtLS9EtVq5Zm1XUylMlwj1y2Td1UcvvR6AfXbbiQ+84QDrsQaQiZakkZeZHwY+3LJs1ZSfNzHN92G7bTVeBvX+gUcdtk+jn0fh+JNeUHp/mh0TLUmSpnHwqZfwwGPbCrXdcU5w/OtfAIeVHJSGijVakiS18cKTLy6cZO26wxxv0qy2TLQkSZpiYsMWFp/01WnnyGq11y7zvEmzpuWpQ0mSmk6Z2FT4ykJwvip1ZqIlSRKN6Ru6mdF9+f4LnEpBHZloSZLGXrdXFq5etojTjzqoxIg0Kky0JEljrZsrC6G66Rs0Gky0JElja6YbQ7facU54ZaG65lWHkqSx0+2VhU7foNlyREuSNFa6vbJwr13msf7klSVGpFFmoiVJGhteWaiqmWhJksaCVxaqDiZakqSR103RO3hlofrHREuSNLImNmzh+As2Fm7vlYXqNxMtSdJI+sZ1d3H8oxsLt991hznes1B95/QOkqSRs+rsK7h6888Lt/fG0CpLZYlWRBwREddHxI0RcVJV/UpSRHw6IrL5uL7N+oiIe5rrn4yIDxbdVoNn5Zq1XV1ZuGTP+U7foNJUkmhFxBzgY8AbgRcDx0TEi6voW9J4i4gdgN8DjgGeDSyJiP/W0uyzwLNofCf+T+AvuthWA+SUiU3ccPdDhduvXraIy05YUV5AGntVjWi9HLgxM2/OzMeB84EjK+pb0ng7DXgiM8/PzPuBG4HjW9q8BfhGNpwMbBcRryu4rQbIeetvL9RuxznBLR95s9M3qHRVFcPvA0z99G8GDm9tFBHHAscCLFq0qPje1zYfkgbLiuajXgcCD095fjtwSEubnYFrpjzfCiwvuC0AEbEOeBnAzjvvXDy6FcWb9s0I97ktn5rCYd2in7VtU2rR+4pydjtQfVbdX1199klViVa0WfaMCU0y8yzgLIClS5cWn/BkBUP9JkgqVaHvnza2dbNtZi6b/Lnr76+qjXCfcyJ+kWytW/TMOq3Sb6ezorxdD0yfVfdXV599UtWpw83AwinP9wXuqKhvSePtGhojVpMWAve1tHkYeMmU53OBdQW31QA55vCF066z6F11qCrR+j6NItL9ImIe8A7gyxX1LWm8/QWwfUS8PSJ2A54P/HVLm4uB1zavPjwDeDIz/7Xgthogpx91EKuXPb30JLDoXfWp5NRhZm6NiOOAS4E5wKcy89oq+pY03jLz4Yg4B7iguejGzPxYRKxvrj8cWA0cATxJ49TgaTNtW2X86t7pRx1kkbsGRmUzw2fmxTT+apSkSmXmO4F3tiw7fMrPCTyn6LaSVJQzw0uSJJXEREuSJKkkJlqSJEklMdGSJEkqiYmWJElSSUy0JEmSSmKiJUmSVJLILH5LripFxD3ArQWbvxC4rqRQ9gN+UtK+dwd+WtK+NdjG5b3v5nU+LzP3KDOYqnT5/eVnYbSMw+sch9cI3b/Ott9hA5todSMitmXmnJL2/VBmzi9p3z/IzKVl7FuDbVze+3F5nb0Yl2Pk6xwd4/AaoX+v01OHkiRJJTHRkiRJKsmoJFrfL3HfXypx32eVuG8NtnF578fldfZiXI6Rr3N0jMNrhD69zpGo0ZIkSRpEozKiJUmSNHAGPtGKiE9HRDYf17dZHxFxT3P9kxHxwS62PTkiHm8+vtZm/Rsj4sHm9v/c5bZHRMT1EXFjRJzUZv2qiLi6+fhuRBxSdFsNrgLv+wciYmPzcU1EbIuIBc11t0TEpua6H1QffXER8amIuDsirplmfUTE3zSPw9UR8dIp6/x8N43DsYiIhRHxzYj4cURcGxF/XHdMZYmIORGxISK+UncsZYmI3SLiCxFxXfM9fUXdMfVbRLy/+Vm9JiLOi4gde9phZg7sA9gBSOAdwG7Ak8B/a2nzT8DjQABnANuKbAtsDzwBvBqYDzwCvLVl3y8Gfg/4DvDPRbcF5gA3Ab8CzAOuAl7csu9XAs9u/vxGYH3RbX0M5qPb9w54K/CNKc9vAXav+3UUfK3/BXgpcM00698EfK35/3KZn+/ePy/D+gCeC7y0+fMuwL+P4utsvr4TgM8BX6k7lhJf42eA9zR/ngfsVndMfX59+9CYO3On5vPPA+/qZZ+DPqJ1GvBEZp6fmfcDNwLHt7R5C41fVpmZJwPbRcTrCmz7LuDnmflvmfkQ8G/Ae6fuODN/lJmfAba29Nlp25cDN2bmzZn5OHA+cGTLvr+bmfc1n64D9i26rQZWt+/dMcB5lUTWZ5n5LeDeGZocCXy2+f9yHbBbRDwXP99TjcWxyMw7M/OHzZ8fBH5M45fZSImIfYE3A5+sO5ayRMSuNP7I+geAzHy8+ft11MwFdoqIucDOwB297GzQE60DgYenPL8dWNDSZmdg6umLrcDyAtsewNNnfL0F2LtgXJ223afZ36TNzPzF8vs0/vqfzbYaHIXfu4jYGTgC+OKUxQl8PSKujIhjS4uyGtMdCz/fTxm7YxERi4HDgPU1h1KGM4E/oXH2ZFT9CnAP8I/NU6SfjIhSJvSuS2ZuAf4SuA24k8agytd72eegJ1rRZlmRyyS3Fdh2tvsusm3hfUfEa2gkWn/ah7hUr27eu7cCl2fm1FGh5Zn5Uhqnkt8bEf+l3wFWaLpj4ef7KWN1LCLil2j8YXF8Zj5Qdzz9FBFvAe7OzCvrjqVkc2mUDPx9Zh4GPASMVG1hRDybxsjyfjQGUOZHxOpe9jnoidY1NEasJi0E7mtp8zDwkinP59I4Fddp2+to3Mdo0mIa2WsRnbbd3Oxv0r60GXqMiINpDDMfmZk/62ZbDaRu3rt30HLaMDPvaP57N3ARjVNLw2q6Y+Hn+yljcywiYnsaSda5mVnm3IR1WQ68LSJuoXEK+LURcU69IZViM7A5MydHJL9AI/EaJa8HfpKZ92TmEzTm0nxlT3usu/CsQ1HazjT+wns7TxW0v7elzbm0L4afcVsaxfJPAL/GUwXtb5smjrU8vRh+xm1pJHs308iIJ4tcD2zZ5yIadWOvbFnecVsfg/ko+t4Bz6JR3zR/yrL5wC5Tfv4ucETdr6nD613M9MXwb+bpxfDf6+YYjcNjXI5F8zPwWeDMumOp6PWuYLSL4b8NHND8+TTgo3XH1OfXdzhwbTOHCBrF/+/rZZ9zGWCZ+XDzr4ILmotuzMyPRcT65vrDgdU0al2epJFYnTbTtlP2/VhE/AXwrzQO5jcz88sRcW5z/aqIOAjYQOPqICJiK427c29pt+2UfW+NiOOAS5vbfiozr42IP2yu/wTw58BzgI9HBMDWzFw63bb9OaIqU8H3HeA3gK9n40KKSXsBFzU/C3OBz2XmJdVF352IOI/GL5TdI2IzcCqNq3EnX+fFNK48vJHGqPO7m+v8fDeN0bFYDrwT2BQRG5vLPpiZF9cXknrwPuDciJhH4w+Fd9ccT19l5vqI+ALwQxo13xvocYZ4Z4aXJEkqyaDXaEmSJA0tEy1JkqSSmGhJkiSVxERLkiSpJCZakiRJJTHRkiRJKomJliRJUklMtDQQImJd84azRMQ+EfGDmkOSpK5ExL4RcXTdcWiwmGipdtGYDn0RcGtz0cHApvoikqRZeR2jd+8/9chES4Pg+TRu4jl5mwITLUlDJSJeBawBfisiNkbEfnXHpMFgoqVBcBBPT6yWAlfXFIskdS0zvwN8HzgyMw/NzJ/UHZMGg4mWBsEC4BGAiHgR8GYc0ZI0fA4Arq87CA0WEy0NgkuB10XE54HfBn6WmXfVHJMkFRYRzwF+nplP1B2LBks8VRYjSZJmIyKWAh/OzDfVHYsGiyNakiT17jpg94i4JiJeWXcwGhyOaEmSJJXEES1JkqSSmGhJkiSVxERLkiSpJCZakiRJJTHRkiRJKomJliRJUklMtCRJkkpioiVJklSS/x89FQ5EoaSvwQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x360 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "u = torch.linspace(0, 1, 1000)\n",
    "\n",
    "f, (ax_t_of_u, ax_u_of_t) = plt.subplots(1, 2, figsize=(10, 5))\n",
    "\n",
    "ax_t_of_u.scatter(u, T(u))\n",
    "ax_t_of_u.set_xlabel(\"$u$\")\n",
    "ax_t_of_u.set_ylabel(\"$T(u)$\")\n",
    "ax_t_of_u.hlines(t_left, 0, 1, linewidth=0.5, label=\"t_left, t_right\", color=\"magenta\")\n",
    "ax_t_of_u.set_xticks(cdf)\n",
    "\n",
    "ax_u_of_t.scatter(T(u), u, )\n",
    "ax_u_of_t.set_xlabel(\"$t$\")\n",
    "ax_u_of_t.set_ylabel(\"$T^{-1}(t)$\")\n",
    "ax_u_of_t.set_yticks(cdf)\n",
    "ax_u_of_t.vlines(t_left, 0, 1, linewidth=0.5, label=\"t_left, t_right\", color=\"magenta\")\n",
    "_= f.suptitle(\"Uniformly sampled $u$'s and corresponding $T(u)$'s\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[0.0000, 1.0000, 2.0000, 3.0000, 4.0000, 5.0000, 6.0000, 7.0000],\n",
       "        [0.0000, 1.0000, 2.0000, 3.0000, 4.0000, 5.0000, 6.0000, 7.0000],\n",
       "        [1.0000, 2.0000, 3.0000, 4.0000, 5.0000, 6.0000, 7.0000, 8.0000],\n",
       "        [0.8823, 1.9150, 2.3829, 3.9593, 4.3904, 5.6009, 6.2566, 7.7936],\n",
       "        [0.0000, 0.0000, 0.0000, 0.1000, 0.1000, 0.2000, 0.7500, 1.0000]])"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "torch.stack((torch.arange(len(cdf)), t_left, t_left + sizes, t, cdf))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The $u$-values of $[0..0.1)$ interpolate between $t$-values of $[3,4]$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([3.0000, 3.0100, 4.0000])"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T(torch.tensor([0, .001, 0.1-1e-8]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is because `searchsorted(cdf, .)` query for these values returns the bucket"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([3, 3, 3])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "torch.searchsorted(cdf, torch.tensor([0, .001, 0.1-1e-8]), right=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note how on the right boundary of this bucket we have the jump"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([3, 5])"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "torch.searchsorted(cdf, torch.tensor([.1-1e-8, .1]), right=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The result is that the query for $u=0.1$ skips the interval $[2,3]$, which has the probability measure of $\\mathrm{cdf}_3 - \\mathrm{cdf}_2 = 0$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([5.0000, 5.9990, 6.0000, 6.9818])"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T(torch.tensor([.1, .1999, .2, .74]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Eval stats"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we sample $t$ using inverse transform and overlay the CDFs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA6IAAAILCAYAAADynCEVAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAg0UlEQVR4nO3df5Tld13f8deb3WD4HSUrP/KDiTVQIwrS3Yj1RxFEEhCitrZEK8rRppxD/FGOQrTU1h9tUauHQ0FzUkRU1KAYadAo1uNBpQpkEQTCD13Cwm7Cjw0gQmQbA+/+ce+yN5PZmbu7dz93ZvbxOGdO5t77nXvf90fms8/5fudOdXcAAABglHssewAAAABOL0IUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCrBFVNVNVfW4k91mjtvZX1XfcDLXsRlV1cur6qeWPcc8quoRVfWWqvpkVX3/nF+zrZ63k7k/W+m5BjhdCVGAkzT9B/Onq+pTMx8vXvTtdPeXdvfrTnabk7Hqvn64qn65qu47c9kdVXX2qq95a1V1Va0c43qOfDz0VM09r6r69qraO53ng1X1B1X1Natm/mRV/V1V/UVVPauq7jHz9Yu6X89N8rruvl93v+gYs26r8Bypqm6pqkfPnN4/+/oE4NQTogCL8dTuvu/Mx5Ujb7yqdg68uad2932TPCbJniTPn7nsfUkun5nry5Lca73rmfm49ZRNPIeqek6SFyb5b0kelOT8JL+Q5LKZzZ7a3fdL8rAkL0jyvCS/tOqqFnG/HpbkphP4OjYw/UHJFyZ517JnATidCVGAU2i6p+WHq+ptVXV7Vf1SVT1ouqftk1X1x1X1+au2/5GqemdVfXy6x/HMmcu+YdW2z6uqtyW5vap2zm5TVedV1XVVdaiqPjq7l7aqrqqq905neGdVfcvx3rfuviXJHyR55MzZv5bkGTOnvyvJrx7vdc875/T+/tD08f1EVb1y5vH6iqr6q+nXvjLJmevczgOS/ESSZ3f3dd19e3f/Y3e/prt/eI37/onuvj7Jv0nyXVX1yNXbzHHfvqSqXjfdu3pTVT1tev6fJPn6JC+e7lF9+Bpf+2uZhPJrpts8d3rRo9d6LKZf89Cq+p3p6+F9tc4hv9PX1S3Tx+49VfWE6fnzPB9zvd7Xe62vMc+6s8/7XFfVFyc5kMm/fz46/f9i56pt1rzvACyWEAU49f5lkicmeXiSp2YSbz+a5OxMvg+vDoLvSPKkJP9k+jXPz7FdnuQpSc7q7juPnFlVO5L8XpL3J1lJck6Sa2e+7r1JvjbJA5L8eJJXVNVDjudOVdV5SZ6c5C0zZ78hyf2nkbUjk1B7xfFc7yrzzPmvk1yS5IIkX57ku6vqnklenUkYf0GS387keTiWr8okXn73eIbr7jclOTidcW5VdUaS1yT5o0z2zn1fkl+vqkd09+OT/HmSK6d7VP9mjdv9ziQfyNG9rz8zvehuj8X09u4xvb2/zuS18IQkP1hVT1pjtkckuTLJnune3ycl2T+9eJ7n43he7xu+1jea/Xie6+7el+SHkrxq+rg9sLvv7O6V7t6/wX0HYIGEKMBivHq6Z+vIx7+buex/dveHp3sQ/zzJG7v7Ld39/zIJn69YdV0v7u4D3f2xJP81M4e6ruFF020/ver8i5M8NMkPT/fuHe7u1x+5sLt/u7tv7e7Pdvcrk/zt9Gvmvq9JXp/kTzM5lHXWkb2iT0zy7iS3rHc9049Xr7XBnHO+aLrNxzIJlkcneWySM5K8cLpn81VJblznPj0wyW2zMX8cbs0kgOa+X9P57pvkBd19R3f/SSY/OFjvuZ7HWo9FMjmEeld3/8T09m5O8r+SPH2N6/hMks9LclFVndHd+7v7vcncz8fxvN7nea1vNPvxPtePSvLWY1x2zPsOwGKN/J0igO3sm7v7j49x2YdnPv/0Gqfvu2r7AzOfvz+ToDyWA8c4/7wk7z9WWFXVM5I8J5O9pZnOcPZa265hvfuaTEL0zzLZK7feYbkbXc+8c35o5vN/yOTxemiSW7q7Zy57/zo39dEkZ1fVzhOI0XOSfGzm9Ib3azrfge7+7Kr5zjnO215trccimfzO6UOnP0A4YkcmoXgX3b2vqn4wyX9J8qVV9dokz+nuW+d8Po7n9T7Pa32j2Y/3uX50JntQ72a9+77O9QFwAuwRBdh8zpv5/PxM9rgdSx/j/ANJzl/9+29JUlUPy2SP0pVJHtjdZyV5R5I6oWlXD9T9/kzetOjJSa470es5yTk/mOScqprd9vx1tv/LJIeTfPNxzrgnk3h8/UbbrnJrkvNq5h13p/Mda+/xWo713K/lQJL3dfdZMx/36+4nr3nF3b/R3V+TSQR2kp8+Ra+beV7rG80+93M9fbwfmWPvEV3zvs99bwCYmxAF2HyeXVXnVtUXZPK7da88get4Uyb/QH9BVd2nqs6sqq+eXnafTP6BfShJquqZuesbDi3C9yR5fHfffhLXcTJz/mWSO5N8f03exOlbs86hx939iSQ/luQlVfXNVXXvqjqjqi6tqp9ZvX1V3b+qvimT37t9RXe//bjuWfLGJLcnee70dh6Xye9TXrveF63y4SRfNOe2b0ry99M34rlXVe2oqkdOQ/ouavI3TB9fVZ+XSZx/OpNDVk/F62ae1/pGsx/Pc32v6cea//5Z574DsGBCFGAxjrx76ZGP43rTm1V+I5M3sbl5+vFTx3sF3f2ZTMLmizN5U5uDmbxxULr7nUl+LpN/wH84yZcl+b8nMe9at//e7t57ktdxwnN29x1JvjWTN+v5eCb3fd29s93985kcdvr8TGLrQCZ7/149s9lrquqT08v+Y5KfT/LMOe/S6vmeluTSJLdl8mdintHd7z6Oq/nvSZ4//V3UH9rg9o68Hh6dyd7q25K8NJM3HVrt8zL50zS3ZXKo7xcm+dFT9LrZ8LW+0ezH81xPfzBydZJ3VtXBNTZZ876f4H0DYB1111+pAGCZqmp/ku+d43cMYUvzWgc4vdkjCgAAwFBCFAAAgKEcmgsAAMBQ9ogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADDUzmXd8Nlnn90rKyvLunkAtpk3v/nNt3X3rmXPsZVZmwFYpPXW5qWF6MrKSvbu3busmwdgm6mq9y97hq3O2gzAIq23Njs0FwAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUBuGaFW9rKo+UlXvOMblVVUvqqp9VfW2qnrM4scEAI6wNgOw1c2zR/TlSS5Z5/JLk1w4/bgiyS+e/FgAwDpeHmszAFvYhiHa3X+W5GPrbHJZkl/tiTckOauqHrKoAQGAu7I2A7DV7VzAdZyT5MDM6YPT8z64esOquiKTn8zm/PPPX8BNH4errko+9KGxtwmwVT34wckLXrDsKThx1maAAVYe/G3LHuGU2P+Cp5zy21hEiNYa5/VaG3b3NUmuSZLdu3evuc0p86EPJSsrQ28SYMvav3/ZE3ByrM0AIxxe9gBb1yLeNfdgkvNmTp+b5NYFXC8AcGKszQBsaosI0euTPGP6Dn2PTfKJ7r7boT8AwDDWZgA2tQ0Pza2q30zyuCRnV9XBJP85yRlJ0t1XJ7khyZOT7EvyD0meeaqGBQCszQBsfRuGaHdfvsHlneTZC5sIAFiXtRmArW4Rh+YCAADA3IQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAw1M5lDwAAADDayuE9yx7htGaPKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIbauewBAAAA5rVyeM+yR2AB7BEFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGConcseAAAA2P5WDu9Z9ghsIvaIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEPNFaJVdUlVvaeq9lXVVWtc/oCqek1V/XVV3VRVz1z8qADAEdZmALayDUO0qnYkeUmSS5NclOTyqrpo1WbPTvLO7n5Ukscl+bmquueCZwUAYm0GYOubZ4/oxUn2dffN3X1HkmuTXLZqm05yv6qqJPdN8rEkdy50UgDgCGszAFvaPCF6TpIDM6cPTs+b9eIkX5Lk1iRvT/ID3f3Z1VdUVVdU1d6q2nvo0KETHBkATnvWZgC2tHlCtNY4r1edflKStyZ5aJJHJ3lxVd3/bl/UfU137+7u3bt27TrOUQGAKWszAFvaPCF6MMl5M6fPzeSnq7OemeS6ntiX5H1J/uliRgQAVrE2A7ClzROiNya5sKoumL7JwdOTXL9qmw8keUKSVNWDkjwiyc2LHBQA+BxrMwBb2s6NNujuO6vqyiSvTbIjycu6+6aqetb08quT/GSSl1fV2zM5XOh53X3bKZwbAE5b1mYAtroNQzRJuvuGJDesOu/qmc9vTfKNix0NADgWazMAW9k8h+YCAADAwghRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMtXPZAwAAAJvbyuE9yx6BbcYeUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgqJ3LHgDYHFYO71n2CAu3/8wblz0CAABrsEcUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQc4VoVV1SVe+pqn1VddUxtnlcVb21qm6qqj9d7JgAwCxrMwBb2c6NNqiqHUlekuSJSQ4mubGqru/ud85sc1aSX0hySXd/oKq+8BTNCwCnPWszAFvdPHtEL06yr7tv7u47klyb5LJV23x7kuu6+wNJ0t0fWeyYAMAMazMAW9o8IXpOkgMzpw9Oz5v18CSfX1Wvq6o3V9Uz1rqiqrqiqvZW1d5Dhw6d2MQAgLUZgC1tnhCtNc7rVad3JvlnSZ6S5ElJ/lNVPfxuX9R9TXfv7u7du3btOu5hAYAk1mYAtrgNf0c0k5+ynjdz+twkt66xzW3dfXuS26vqz5I8KsnfLGRKAGCWtRmALW2ePaI3Jrmwqi6oqnsmeXqS61dt87+TfG1V7ayqeyf5yiTvWuyoAMCUtRmALW3DPaLdfWdVXZnktUl2JHlZd99UVc+aXn51d7+rqv4wyduSfDbJS7v7HadycICNrBzes+wRPmf/mTcuewS2EWszAFvdPIfmprtvSHLDqvOuXnX6Z5P87OJGAwCOxdoMwFY2V4gCAABbz2Y6OghmzfM7ogAAALAwQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAy1c9kDAAAAd7dyeM+yR4BTxh5RAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYaueyBwBOzsrhPcseAQAAjos9ogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDzRWiVXVJVb2nqvZV1VXrbLenqj5TVf9qcSMCAKtZmwHYyjYM0arakeQlSS5NclGSy6vqomNs99NJXrvoIQGAo6zNAGx18+wRvTjJvu6+ubvvSHJtksvW2O77kvxOko8scD4A4O6szQBsafOE6DlJDsycPjg973Oq6pwk35Lk6vWuqKquqKq9VbX30KFDxzsrADBhbQZgS5snRGuN83rV6RcmeV53f2a9K+rua7p7d3fv3rVr15wjAgCrWJsB2NJ2zrHNwSTnzZw+N8mtq7bZneTaqkqSs5M8uaru7O5XL2JIAOAurM0AbGnzhOiNSS6sqguS3JLk6Um+fXaD7r7gyOdV9fIkv2ehA4BTxtoMwJa2YYh2951VdWUm77i3I8nLuvumqnrW9PJ1f/cEAFgsazMAW908e0TT3TckuWHVeWsuct393Sc/FgCwHmszAFvZPG9WBAAAAAsjRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACG2rnsAQBOByuH9yzkevafeeNCrgcAYJmEKAAALNCifvgI25lDcwEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMNTOZQ8AwPxWDu8Zc0MP3pNc9ftDbmr/C54y5HYAgM3DHlEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGCouUK0qi6pqvdU1b6qumqNy7+jqt42/fiLqnrU4kcFAI6wNgOwlW0YolW1I8lLklya5KIkl1fVRas2e1+Sf9HdX57kJ5Ncs+hBAYAJazMAW908e0QvTrKvu2/u7juSXJvkstkNuvsvuvvj05NvSHLuYscEAGZYmwHY0uYJ0XOSHJg5fXB63rF8T5I/WOuCqrqiqvZW1d5Dhw7NPyUAMMvaDMCWNk+I1hrn9ZobVn19Jovd89a6vLuv6e7d3b17165d808JAMyyNgOwpe2cY5uDSc6bOX1ukltXb1RVX57kpUku7e6PLmY8AGAN1mY4RVYO71n2CHBamGeP6I1JLqyqC6rqnkmenuT62Q2q6vwk1yX5zu7+m8WPCQDMsDYDsKVtuEe0u++sqiuTvDbJjiQv6+6bqupZ08uvTvJjSR6Y5BeqKknu7O7dp25sADh9WZsB2OrmOTQ33X1DkhtWnXf1zOffm+R7FzsaAHAs1mYAtrJ5Ds0FAACAhRGiAAAADCVEAQAAGEqIAgAAMJQQBQAAYKi53jUXAAA2s5XDe5Y9AnAc7BEFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFA7lz0AAACnt5XDe5Y9AjCYPaIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMNTOZQ8Ao60c3rPsEQAA4LS2LUJ05arf33ijB39bcvjUz7Io+8+8cdkjAACsyw93gRPl0FwAAACG2hZ7RLejRf2E0Z5VAABgs7FHFAAAgKHsEQUA2IDfhQRYLHtEAQAAGEqIAgAAMJRDc7e5RRxK5A2PAACARbJHFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABjKu+YCANvaIt5BHoDFskcUAACAoYQoAAAAQzk0lw0t6pCm/WfeuJDrAQAAtjZ7RAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABD+fMtAMCmtKg/HwbA5mOPKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAof76FLcPb+AMAwPZgjygAAABDCVEAAACGEqIAAAAMJUQBAAAYypsVAQALt3LV72+80YO/LTl86mcBYPOxRxQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoeYK0aq6pKreU1X7quqqNS6vqnrR9PK3VdVjFj8qAHCEtRmArWzDEK2qHUlekuTSJBclubyqLlq12aVJLpx+XJHkFxc8JwAwZW0GYKubZ4/oxUn2dffN3X1HkmuTXLZqm8uS/GpPvCHJWVX1kAXPCgBMWJsB2NJ2zrHNOUkOzJw+mOQr59jmnCQfnN2oqq7I5KeySfKpqnrPcU17Ys5Octv5yTlnJGcMuL1N6xPJmQ9IDi/r9r/4JL9+nhfrvJb9WGwmHoujPBZHfSw54+PJ3464rfrphV3VwxZ2TZuftXmb8H3nKI/FUR6LozwWR223tXmef9vXGuf1CWyT7r4myTVz3ObCVNXe7t498jY3q6rae8hjkcRjMctjcZTH4ijfOzc9a/M24fvOUR6LozwWR3ksjtpu3zvnOTT3YJLzZk6fm+TWE9gGAFgMazMAW9o8IXpjkgur6oKqumeSpye5ftU21yd5xvQd+h6b5BPd/cHVVwQALIS1GYAtbcNDc7v7zqq6Mslrk+xI8rLuvqmqnjW9/OokNyR5cpJ9Sf4hyTNP3cjHbejhRpucx+Ioj8VRHoujPBZHeSw2MWvztuKxOMpjcZTH4iiPxVHb6rGo7rv9uggAAACcMvMcmgsAAAALI0QBAAAYatuGaFWdWVVvqqq/rqqbqurHlz3TMlXVjqp6S1X93rJnWbaq2l9Vb6+qt1bV3mXPs0xVdVZVvaqq3l1V76qqr1r2TMtQVY+Yvh6OfPx9Vf3gsudahqr6D9Pvme+oqt+sqjOXPRPbh7X5rqzNR1mbj7I2T1ibj9qua/O2/R3Rqqok9+nuT1XVGUlen+QHuvsNSx5tKarqOUl2J7l/d3/TsudZpqran2R3d9+27FmWrap+Jcmfd/dLp++8ee/u/rslj7VUVbUjyS1JvrK737/seUaqqnMy+V55UXd/uqp+K8kN3f3y5U7GdmFtvitr81HW5qOszXdnbd6ea/O23SPaE5+anjxj+rE9q3sDVXVukqckeemyZ2HzqKr7J/m6JL+UJN19x+m+0E09Icl7T7eFbsbOJPeqqp1J7h1/d5IFsjYfZW1mLdbmY7I2b8O1eduGaPK5Q17emuQjSf5Pd79xySMtywuTPDfJZ5c8x2bRSf6oqt5cVVcse5gl+qIkh5L88vTQsJdW1X2WPdQm8PQkv7nsIZahu29J8j+SfCDJBzP5u5N/tNyp2G6szZ/zwlibZ1mbJ6zNa7M2b8O1eVuHaHd/prsfneTcJBdX1SOXPNJwVfVNST7S3W9e9iybyFd392OSXJrk2VX1dcseaEl2JnlMkl/s7q9IcnuSq5Y70nJND4F6WpLfXvYsy1BVn5/ksiQXJHlokvtU1b9d7lRsN9Zma/MxWJsnrM2rWJu379q8rUP0iOkhDa9LcslyJ1mKr07ytOnvXlyb5PFV9YrljrRc3X3r9L8fSfK7SS5e7kRLczDJwZm9Ea/KZPE7nV2a5K+6+8PLHmRJviHJ+7r7UHf/Y5LrkvzzJc/ENmVttjbPsjZ/jrX57qzN23Rt3rYhWlW7quqs6ef3yuRJfPdSh1qC7v6R7j63u1cyOazhT7p7W/wU5URU1X2q6n5HPk/yjUnesdyplqO7P5TkQFU9YnrWE5K8c4kjbQaX5zQ99GfqA0keW1X3nr6pzBOSvGvJM7GNWJsnrM13ZW0+ytq8JmvzNl2bdy57gFPoIUl+ZfouW/dI8lvdfdq/PTp5UJLfnfx/nJ1JfqO7/3C5Iy3V9yX59elhLzcneeaS51maqrp3kicm+ffLnmVZuvuNVfWqJH+V5M4kb0lyzXKnYpuxNrMWa/NdWZunrM3be23etn++BQAAgM1p2x6aCwAAwOYkRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABD/X+CuSIQfx2OIgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x576 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "torch.manual_seed(42)\n",
    "u = torch.rand(1000)\n",
    "Tu = T(u)\n",
    "Tu = torch.sort(Tu).values\n",
    "\n",
    "# _1, _2 = np.histogram(Tu, bins=torch.cat((t_left, (t_left + sizes)[-1:])), normed=True)\n",
    "\n",
    "_1, _2 = np.histogram(Tu, bins=np.linspace(Tu.min(), Tu.max(), 20), density=True)\n",
    "\n",
    "f, (axd, axc) = plt.subplots(1, 2, figsize=(16, 8))\n",
    "axd.bar(_2[1:], _1)\n",
    "axc.bar(_2[1:], np.cumsum(_1 / _1.sum()))\n",
    "axd.fill_between([Tu.min(), Tu.max()], 1, color=\"red\", alpha=.5)\n",
    "axc.fill_between([Tu.min(), Tu.max()], 1, color=\"red\", alpha=.5)\n",
    "\n",
    "_ = f.suptitle(\"Empirical PMF and CDF of the sampled $t$'s\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeMAAAHSCAYAAADfUaMwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAax0lEQVR4nO3de5CddZ3n8c/PDhCQYBCipSQx0UkkMQkBQwQvLCzLbWSAZXQGZplBSifDIDprLSpuwS7MYM2F1XFwgAxeFqawkBItNkvFSzkbLyPjmBARAkiMJIQecAgJsFwnBH/7R3d6O52GPiHd+SXdr1dVqvqc8+vn+SYnybue5/R5Tqm1BgBo51WtBwCAsU6MAaAxMQaAxsQYABoTYwBoTIwBoLFxrXZ88MEH12nTprXaPQDscnfcccdjtdZJA+9vFuNp06ZlxYoVrXYPALtcKeXBwe53mhoAGhNjAGhMjAGgsWavGQ/mhRdeSHd3d55//vnWo7CbGD9+fCZPnpy99tqr9SgAI2a3inF3d3cmTJiQadOmpZTSehwaq7Vm48aN6e7uzvTp01uPAzBidqvT1M8//3wOOuggISZJUkrJQQcd5EwJMOrtVjFOIsRsw98HYCzY7WLcWldXV+bPn9/36y/+4i+S9LyeffHFF2fGjBmZM2dOFi5cmG9+85tJet4zPXfu3MydOzezZ8/OJZdckn/7t39Lkqxbty777rvvNtvcvHlzx/OcffbZmTdvXv76r/96m/tvvfXW3HvvvX23jz322D3mfdt70qwAu8Ju9Zrxdi67bJdvb999982dd9653f2XXnppHnnkkaxatSr77LNP/vVf/zXf//73+x5ftmxZDj744Dz99NNZtGhRFi1alBtuuCFJ8pa3vGXQbQ7lV7/6VW6//fY8+OD27xG/9dZbc+qpp2b27Nk7vF0Adi+OjDvw7LPP5gtf+EI+//nPZ5999kmSvP71r8/v/M7vbLd2//33z+LFi3Prrbdm06ZNHW3/+eefz3nnnZe5c+fm8MMPz7Jly5IkJ554Yh599NHMnz8/P/zhD/vW33777VmyZEk+/vGPZ/78+fnlL3+ZJPna176WhQsXZubMmX3rX3zxxXz84x/PkUcemXnz5uXv/u7vttv/M888k/e+97057LDDMmfOnNx8881Jkj/90z/NkUcemTlz5mTRokWptSbpObL92Mc+lmOOOSazZs3K8uXLc+aZZ2bGjBm55JJLkvScETj00ENz7rnnZt68eXnf+96XZ599drt9f+c738nRRx+dI444Iu9///vz9NNPd/RnBjCaiPEAzz333DanlG+++easWbMmU6dOzQEHHNDRNg444IBMnz49v/jFL5Ikv/zlL/u29+EPf3i79VdffXWS5O67785NN92Uc889N88//3yWLFnSd1T9nve8p2/9O9/5zpx22mm58sorc+edd+Ytb3lLkmTLli35yU9+ks997nO5/PLLkyRf+tKX8prXvCbLly/P8uXL84UvfCFr167dZv/f+ta38sY3vjE/+9nPsmrVqpx88slJkgsvvDDLly/PqlWr8txzz+W2227r+5699947P/jBD3L++efn9NNPz9VXX51Vq1bl+uuvz8aNG5Mk999/fxYtWpS77rorBxxwQK655ppt9vvYY4/liiuuyHe/+92sXLkyCxYsyGc/+9mO/owBRpPd+zR1A4Odpr7rrrt2eDtbjyKToU9T/+M//mM+8pGPJEkOPfTQvOlNb8rq1as7jv9WZ555ZpLk7W9/e9atW5ek58jzrrvuyi233JIkefLJJ/OLX/xim7cKzZ07NxdddFE++clP5tRTT+0L/7Jly/JXf/VXefbZZ7Np06a87W1vy2/91m8lSU477bS+733b296WN7zhDUmSN7/5zXnooYcyceLETJkyJe9617uSJOecc06uuuqqXHTRRX37/fGPf5x77723b83mzZtz9NFH79DvGWA0EOMO/MZv/EbWr1+fp556KhMmTBhy/VNPPZV169Zl5syZefLJJ4dc3z/cO2PrKfSurq5s2bKlb9uf//znc9JJJ73k982cOTN33HFHli5dmk996lM58cQT84lPfCIXXHBBVqxYkSlTpuSyyy7b5i1GW/f1qle9qu/rrbe37nvgT0IPvF1rzQknnJCbbrppJ37XAHs+p6k7sN9+++WDH/xgPvrRj/b9JPQjjzySG2+8cbu1Tz/9dC644IKcccYZOfDAAzva/jHHHJOvfOUrSZLVq1dn/fr1eetb3/qy3zNhwoQ89dRTQ277pJNOyrXXXpsXXnihb/vPPPPMNmsefvjh7LfffjnnnHNy0UUXZeXKlX3h3fpDaVuPrHfE+vXr80//9E9Jkptuuinvfve7t3n8qKOOyo9+9KOsWbMmSc9r86tXr97h/QDs6cR4gIGvGV988cVJkiuuuCKTJk3K7NmzM2fOnJxxxhmZNOn/fyTlcccd1/eWp6lTpw76g1Iv5YILLsiLL76YuXPn5nd/93dz/fXXb3O0OZizzjorV155ZQ4//PC+H+AazIc+9KHMnj07RxxxRObMmZM/+qM/6jty3eruu+/OwoULM3/+/Hz605/OJZdckokTJ+YP//APM3fu3Jxxxhk58sgjO/79bDVr1qzccMMNmTdvXjZt2pQ//uM/3ubxSZMm5frrr+97+9ZRRx2Vn//85zu8H4A9XRnqFGkp5ctJTk3yaK11ziCPlyR/k+Q3kzyb5AO11pVD7XjBggV14HtN77vvvsyaNavz6dltrVu3LqeeempWrVq109vy9wIYLUopd9RaFwy8v5Mj4+uTnPwyj5+SZEbvr0VJrn0lAwLAWDVkjGutP0jycm+YPT3J39ceP04ysZTyhuEakD3TtGnThuWoGGAsGI7XjA9J8lC/29299wEAHRiOtzYNdiX/QV+ILqUsSs+p7EydOnUYdg3ASDjuhuN2zY7WrhuRzS7LuTu/keG+JPPLGI4j4+4kU/rdnpzk4cEW1lqvq7UuqLUu6P+TyAAwlg1HjJck+YPS46gkT9ZaHxmG7QLAmDDkaepSyk1Jjk1ycCmlO8l/T7JXktRaFydZmp63Na1Jz1ubzhupYUfaxo0bc/zxxyfp+cSkrq6uvvcS/+QnP8nee++90/t44YUXcumll+brX/969tlnn+y33365/PLLc8opp2TatGl9V/h68cUXc+aZZ+bSSy/NPvvsk3Xr1mXWrFnbXAxkuGYCoK0hY1xrPXuIx2uS7T/9YBgM92sWy85d9rKPH3TQQX3XkL7sssuy//77b3Mt5S1btmTcuJ17mb3FRzECsHtzbeohfOADH8hrX/va/PSnP80RRxyRCRMmbBPpOXPm5Lbbbsu0adNy44035qqrrsrmzZvzjne8I9dcc026urr6trX1oxjXrl3b8UcxTpkypeOPYgRgz+RymB1YvXp1vvvd7+Yzn/nMS6657777cvPNN+dHP/pR7rzzznR1dfVdb3qrkf4oRgD2TI6MO/D+979/myPcwfzDP/xD7rjjjr5rOD/33HN53etet9P73pGPYgRgzyTGHXj1q1/d9/W4cePy61//uu/21k83qrXm3HPPzZ//+Z+/5HZG+qMYAdgzOU29g6ZNm5aVK3s+B2PlypVZu3ZtkuT444/PLbfckkcffTRJsmnTpjz44IPbfO9IfxQjAHsmMd5Bv/3bv51NmzZl/vz5ufbaazNz5swkyezZs3PFFVfkxBNPzLx583LCCSfkkUe2f7v1SH4UIwB7piE/QnGk+AhFOuXvBex6LoeZEbkc5s58hCIAMILEGAAaE2MAaGy3i3Gr17DZPfn7AIwFu1WMx48fn40bN/oPmCQ9Id64cWPGjx/fehSAEbVbXfRj8uTJ6e7uzoYNG1qPwm5i/PjxmTx5cusxAEbUbhXjvfbaK9OnT289BgDsUrvVaWoAGIvEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaG9d6AAB2Q2vX7Zr9PPH4jq2feODIzNGYI2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGisoxiXUk4updxfSllTSrl4kMdfU0r536WUn5VS7imlnDf8owLA6DRkjEspXUmuTnJKktlJzi6lzB6w7MNJ7q21Hpbk2CSfKaXsPcyzAsCo1MmR8cIka2qtD9RaNyf5apLTB6ypSSaUUkqS/ZNsSrJlWCcFgFGqkxgfkuShfre7e+/r72+TzErycJK7k/xJrfXXAzdUSllUSllRSlmxYcOGVzgyAIwuncS4DHJfHXD7pCR3JnljkvlJ/raUcsB231TrdbXWBbXWBZMmTdrBUQFgdOokxt1JpvS7PTk9R8D9nZfkG7XHmiRrkxw6PCMCwOjWSYyXJ5lRSpne+0NZZyVZMmDN+iTHJ0kp5fVJ3prkgeEcFABGq3FDLai1bimlXJjk20m6kny51npPKeX83scXJ/mzJNeXUu5Oz2ntT9ZaHxvBuQFg1BgyxklSa12aZOmA+xb3+/rhJCcO72gAMDa4AhcANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAYx19njEAnTvuhuNaj8AexpExADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNjWs9AAAdWrtu1+3ricd3bP3EA0dmjjHCkTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA01lGMSyknl1LuL6WsKaVc/BJrji2l3FlKuaeU8v3hHRMARq9xQy0opXQluTrJCUm6kywvpSyptd7bb83EJNckObnWur6U8roRmhcARp1OjowXJllTa32g1ro5yVeTnD5gze8l+UatdX2S1FofHd4xAWD06iTGhyR5qN/t7t77+puZ5MBSyvdKKXeUUv5gsA2VUhaVUlaUUlZs2LDhlU0MAKNMJzEug9xXB9wel+TtSd6b5KQkl5ZSZm73TbVeV2tdUGtdMGnSpB0eFgBGoyFfM07PkfCUfrcnJ3l4kDWP1VqfSfJMKeUHSQ5LsnpYpgSAUayTI+PlSWaUUqaXUvZOclaSJQPW/K8k7ymljCul7JfkHUnuG95RAWB0GvLIuNa6pZRyYZJvJ+lK8uVa6z2llPN7H19ca72vlPKtJHcl+XWSL9ZaV43k4AAwWnRymjq11qVJlg64b/GA21cmuXL4RgOAscEVuACgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxrUeAGCPdtllg9y5bmT29cTjO7Z+4oEjMwfDzpExADQmxgDQWEcxLqWcXEq5v5SyppRy8cusO7KU8mIp5X3DNyIAjG5DxriU0pXk6iSnJJmd5OxSyuyXWPeXSb493EMCwGjWyZHxwiRraq0P1Fo3J/lqktMHWfeRJF9P8ugwzgcAo14nMT4kyUP9bnf33tenlHJIkv+YZPHwjQYAY0MnMS6D3FcH3P5ckk/WWl982Q2VsqiUsqKUsmLDhg0djggAo1sn7zPuTjKl3+3JSR4esGZBkq+WUpLk4CS/WUrZUmu9tf+iWut1Sa5LkgULFgwMOgCMSZ3EeHmSGaWU6Un+JclZSX6v/4Ja6/StX5dSrk9y28AQAwCDGzLGtdYtpZQL0/NT0l1JvlxrvaeUcn7v414nBoCd0NHlMGutS5MsHXDfoBGutX5g58cCgLHDFbgAoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGOopxKeXkUsr9pZQ1pZSLB3n8P5VS7ur9dXsp5bDhHxUARqchY1xK6UpydZJTksxOcnYpZfaAZWuT/Lta67wkf5bkuuEeFABGq06OjBcmWVNrfaDWujnJV5Oc3n9BrfX2WuvjvTd/nGTy8I4JAKNXJzE+JMlD/W539973Uj6Y5Js7MxQAjCXjOlhTBrmvDrqwlOPSE+N3v8Tji5IsSpKpU6d2OCIw1hx3w3GtR9gB61oPwCjQyZFxd5Ip/W5PTvLwwEWllHlJvpjk9FrrxsE2VGu9rta6oNa6YNKkSa9kXgAYdTqJ8fIkM0op00speyc5K8mS/gtKKVOTfCPJ79daVw//mAAweg15mrrWuqWUcmGSbyfpSvLlWus9pZTzex9fnOS/JTkoyTWllCTZUmtdMHJjA8Do0clrxqm1Lk2ydMB9i/t9/aEkHxre0QBgbOgoxgB7lLXrdt2+nnh86DX9TTxwZOZgj+ZymADQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQmxgDQmBgDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JgYA0BjYgwAjYkxADQ2rvUAwBhx2WU7sHjdzu3ricd3bP3EA3duf7CTHBkDQGNiDACNiTEANCbGANCYGANAY2IMAI2JMQA0JsYA0JiLfsDuZocujrGTvve9HVt/7LG7Zl/zH3chDsYUR8YA0JgjY+jnuBuOaz1CdvpSkDti/g5eNjIP7sJ9wdjhyBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbGtR4A9ghr1+26fT3x+I6tn3jgyMwB7DKOjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBoTIwBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaEyMAaAxMQaAxsQYABoTYwBorKMYl1JOLqXcX0pZU0q5eJDHSynlqt7H7yqlHDH8owLA6DRkjEspXUmuTnJKktlJzi6lzB6w7JQkM3p/LUpy7TDPCQCjVidHxguTrKm1PlBr3Zzkq0lOH7Dm9CR/X3v8OMnEUsobhnlWABiVOonxIUke6ne7u/e+HV0DAAxiXAdryiD31VewJqWURek5jZ0kT5dS7u9g/7uTg5M81nqIMc5zsJ0nd/W+dtFzsMt/X3vavjp4HvbI39ew7Kfk8p3f1eVDbuOV/Ft402B3dhLj7iRT+t2enOThV7AmtdbrklzXwT53S6WUFbXWBa3nGMs8B+15DnYPnof2hvM56OQ09fIkM0op00speyc5K8mSAWuWJPmD3p+qPirJk7XWR4ZjQAAY7YY8Mq61bimlXJjk20m6kny51npPKeX83scXJ1ma5DeTrEnybJLzRm5kABhdOjlNnVrr0vQEt/99i/t9XZN8eHhH2y3tsafYRxHPQXueg92D56G9YXsOSk9HAYBWXA4TABoT4w4MdTlQRl4pZUopZVkp5b5Syj2llD9pPdNYVUrpKqX8tJRyW+tZxqJSysRSyi2llJ/3/ns4uvVMY00p5WO9/w+tKqXcVEoZv7PbFOMhdHg5UEbeliT/pdY6K8lRST7seWjmT5Lc13qIMexvknyr1npoksPiudilSimHJPlokgW11jnp+cHms3Z2u2I8tE4uB8oIq7U+Umtd2fv1U+n5D8hV3naxUsrkJO9N8sXWs4xFpZQDkhyT5EtJUmvdXGt9oulQY9O4JPuWUsYl2S+DXFdjR4nx0FzqczdTSpmW5PAk/9x4lLHoc0k+keTXjecYq96cZEOS/9n7UsEXSymvbj3UWFJr/Zck/yPJ+iSPpOe6Gt/Z2e2K8dA6utQnu0YpZf8kX0/yn2ut/7f1PGNJKeXUJI/WWu9oPcsYNi7JEUmurbUenuSZJH6OZRcqpRyYnrOj05O8McmrSynn7Ox2xXhoHV3qk5FXStkrPSH+Sq31G63nGYPeleS0Usq69Lxc8+9LKTe2HWnM6U7SXWvdelbolvTEmV3nPyRZW2vdUGt9Ick3krxzZzcqxkPr5HKgjLBSSknP62T31Vo/23qesajW+qla6+Ra67T0/Dv4P7XWnT4ioHO11l8leaiU8tbeu45Pcm/Dkcai9UmOKqXs1/v/0vEZhh+i6+gKXGPZS10OtPFYY9G7kvx+krtLKXf23vdfe68OB2PJR5J8pffg4IG4/PAuVWv951LKLUlWpuddHj/NMFyJyxW4AKAxp6kBoDExBoDGxBgAGhNjAGhMjAGgMTEGgMbEGAAaE2MAaOz/AYZ8h+iaRZJgAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "_1, _2 = np.histogram(Tu, bins=np.linspace(Tu.min(), Tu.max(), 2 + len(t_left)), density=True)\n",
    "ecdf = np.cumsum(_1 / _1.sum())\n",
    "                 \n",
    "f, ax = plt.subplots(figsize=(8, 8))\n",
    "ax.bar(_2[:-1], ecdf, alpha=.5, color=\"red\", label=\"ECDF of the sample\")\n",
    "ax.bar(t_left, width=sizes, height=cdf, alpha=.75, color=\"green\", label=\"True CDF\")\n",
    "_ = ax.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "nerfs",
   "language": "python",
   "name": "nerfs"
  },
  "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}