asr_lab2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# ASR Lab 2 - Computing HMM probabilities\n",
    "\n",
    "To begin with, we'll use the code you wrote in Lab 1 to generate a Word WFST for the word \"*peppers*\", using `generate_word_wfst('peppers')`.  By viewing this as an HMM, you'll be able to sample possible paths through the model and also generate the likelihood of an observation sequence $(x_1, \\dotsc, x_T)$.\n",
    "\n",
    "We have also created [technical documentation](https://openfst-python-documentation.readthedocs.io/en/latest/lab2.html) to accompany this lab.\n",
    "\n",
    "In this lab we only compute likelihoods along single (random) paths.  In Labs 3 and 4 we'll build on this to implement the basics of the Viterbi algorithm, which can later be used for word recognition.\n",
    "\n",
    "First, copy your code from Lab 1 into the space below.  You can use the official solutions if you like.\n",
    "If you want to extract the code-only parts of your previous notebook, on the terminal command line you can type:\n",
    "\n",
    "```bash\n",
    "jupyter nbconvert --to python <notebook-name.ipynb>\n",
    "```\n",
    "\n",
    "where <notebook-name.ipynb> indicates the path of the notebook file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openfst_python as fst\n",
    "\n",
    "def parse_lexicon(lex_file):\n",
    "    \"\"\"\n",
    "    Parse the lexicon file and return it in dictionary form.\n",
    "    \n",
    "    Args:\n",
    "        lex_file (str): filename of lexicon file with structure '<word> <phone1> <phone2>...'\n",
    "                        eg. peppers p eh p er z\n",
    "\n",
    "    Returns:\n",
    "        lex (dict): dictionary mapping words to list of phones\n",
    "    \"\"\"\n",
    "    \n",
    "    lex = {}  # create a dictionary for the lexicon entries (this could be a problem with larger lexica)\n",
    "    with open(lex_file, 'r') as f:\n",
    "        for line in f:\n",
    "            line = line.split()  # split at each space\n",
    "            lex[line[0]] = line[1:]  # first field the word, the rest is the phones\n",
    "    return lex\n",
    "\n",
    "def generate_symbol_tables(lexicon, n=3):\n",
    "    '''\n",
    "    Return word, phone and state symbol tables based on the supplied lexicon\n",
    "        \n",
    "    Args:\n",
    "        lexicon (dict): lexicon to use, created from the parse_lexicon() function\n",
    "        n (int): number of states for each phone HMM\n",
    "        \n",
    "    Returns:\n",
    "        word_table (fst.SymbolTable): table of words\n",
    "        phone_table (fst.SymbolTable): table of phones\n",
    "        state_table (fst.SymbolTable): table of HMM phone-state IDs\n",
    "    '''\n",
    "    \n",
    "    state_table = fst.SymbolTable()\n",
    "    phone_table = fst.SymbolTable()\n",
    "    word_table = fst.SymbolTable()\n",
    "    \n",
    "    # add empty <eps> symbol to all tables\n",
    "    state_table.add_symbol('<eps>')\n",
    "    phone_table.add_symbol('<eps>')\n",
    "    word_table.add_symbol('<eps>')\n",
    "    \n",
    "    for word, phones  in lexicon.items():\n",
    "        \n",
    "        word_table.add_symbol(word)\n",
    "        \n",
    "        for p in phones: # for each phone\n",
    "            \n",
    "            phone_table.add_symbol(p)\n",
    "            for i in range(1,n+1): # for each state 1 to n\n",
    "                state_table.add_symbol('{}_{}'.format(p, i))\n",
    "            \n",
    "    return word_table, phone_table, state_table\n",
    "\n",
    "\n",
    "# call these two functions\n",
    "lex = parse_lexicon('lexicon.txt')\n",
    "word_table, phone_table, state_table = generate_symbol_tables(lex)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABqQAAABNCAIAAABPKZIfAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydZ0ATSRvHN/RepZcUescCiJ6FooCIiAqiIiqK3bP3O/HsHT17Q7HgqYCIWAGxixRFQcACShHpECCQEJJ9P+TeXAwBZpcESNjfp5SZ3Sf/zDPl2Sk4GIYhDAwMDAwMDAwMDAwMDAwMDAwMDOFHrLcNwMDAwMDAwMDAwMDAwMDAwMDAwOAPWLAPAwMDAwMDAwMDAwMDAwMDAwNDRMCCfRgYGBgYGBgYGBgYGBgYGBgYGCICFuzDwMDAwMDAwMDAwMDAwMDAwMAQEbBgHwYGBgYGBgYGBgZGv6Xi/iYPM3UFJd1BM068pfS2NSIN89vNhU46iko6DvP+KWD2tjV9E7jxc/y+kKHuu/J72xJRBHN2wYF5d58DTbAPq4AEA1b1CIJ+Vem0/Xz296yB4/4u7SId5sLdA1BnzKO7CaDO/crH0YN5vQDAfFxAYE79C5jz9gC0xzfSnQ6nlpalHbJNXTZlR0bHSZn1H+MOhLrMjSIjvUnl+Vnui8MTPjXA3TFV6CmJv90Yeq+g+NmGATGL9z7pJGU/lroiauOSQ/+8KW/tIl0/lgg1mLMLEMy7+yAwcn6emqSvJQtZheV1no5Rl3Nr/7zRIVfrkd6h4lyw26JDd/LJTBTmCSnU5L//isutbSbnXZtFEiNuSO8gHaYqIopvhZ/PqKbUvjvsqaIamtxJSlEQNmOTuYwY5BRe0kU6zIW7B5DOoB4NYzp3BFh57l8+jh4gr+/fEiEF83FBgTn1r2BNds9SdXQUNOZ0TfsvmA2f7x5eOsbSdHjwtn8yK+nIL91amXn1ryAnU6txy48lFjR131ahpK2tjfXi085BjvsKeCXBpIZhGE5ZMqBjr8ck4geYs/MbzLv7HmiCfTCMVUAChXfVg6mKhv5W6bTFT5ftOtgHw5gLdw9wnWEY60ygB0Tn/ubj3aFjr8ck6h6Yj/MVzKnbgzXZiKHmRy0eSVJT1LCeeuJ9C0wve3Z0zmD30E1+JFWTDS8ZHWesOjZab9UrJgzXXJtqGvhPNQzDlG/JJ1b62JoPC94WlVnR2sEtWgrv7ZlqFRx+a4Wjpryaue/Rd81wc9qesURlHbeDN6/GZf8/H63sZeQf/kNM7CeuPfekqEWgMvQEKKRurf54d+dEn93ZrNLaz6S+ORX336wbiSF78nkne71a7xev/1elDiSCf1XpTa5QS4QaQG3/BXN2hADJy+XdMFZ0exfOYB+vgtsRv1RAmId0AhJV/4W76kHvGzDMQ1hRUBVGIqwIVzpcngXDcOJcZae9iRGBFqqqxHHHcjruzmIujIBu6AzDWGcCGNQ6i7CPdw7C9oWX14u6ROhA2nBjPg4I6oZbxITltDazECQIhQ37kZK7y8P/TG5dU/nLzYOkhx7+lrrOWEYM0g++9bOt84w/Tk+cdP7nLx/9DB8upusV/rax81vcmaMsBUF2y+7n15MLrk4zkHDcVxAXZLPiSU196etDK49k/5odJqfuG6uJcztb1/1f26ugkLrlsjcEQRAka7PuBYXzi34hdeuNQ4cKqEwYLrsV7DDhXGFHMnF7PYsOJIK5VIIgSIglQg2otiwwZ0cIkLxIvRvGiq5g4Qz28Si4hR3lQ1QB9WsPQaIqC+6qB71vFMIws0NhhVpVGIGwIlzpcHsWDMOJc5WJ/iefldb/uOinYhXGbdx/YC4MTnd0hrHOBDCodRZhH+8ChO0LD68XeYnQgVBYzMdBEVjDLVzC/mItTloCIAiFDfsRUnr4N44ti8Sm/EODGXdmyHU5Jb8yes2quHLuT5sKk46v8LY2cZwRdiW9vLXDW1SfGg35XmMFqHO2mov5Xc875Yq3mrT26L1PHH8Srez15bBpQ80H+q45myLkQVW0UreSS9JvrnFWkph87RcB+o/UTWlhI0b+kdbJPFreXs9bIri9Ssf/FnaJUNO1tjCMOTtqupQXkXfDWNEVML8u421fcFt5ZYKRVUD93UPAVYVhnlUPet9oheGv7YQVDVVhBMKKaqXDw7PgxLnK/3awPmyxJKxJ7TAz5sLAdEtnrDMBTHd0FlUf7xpE7Qsvrxd9idCBQFjMx5EgmIZbuITltjYqpusgFDbsR0jRQWe77VyC/degdEBL1sXd1/KoHX3NbPicEL7E3cL0t9k7o7OqC9vfouWyNzQhivU/fNltpzjnIQy3Vb+P3hfqrG+86H4d9cfryD+nOxhbef1+9JGILJdGJ/W/PJqvHXCD1v5zkZea8e2Sn9PM2+0ajl/owOtZcEnU1v6PEHKJUAOkLebsaAGTF4aBvRsruoLm12AfL3F5g6QC6u8eAq5qp1UPOt+AOYQNOX9fdFSFEQkLw6JY6XTewcoOs8Kvft1hZsyFgemGzlhnAgHdKs8wDIuij3cNomqwY68XZYnQASos5uMIEUzDLVzCtqvrACIj2LAfIcysMHNN1+2JX+obq/MTDl/OZNAqIrylB+7luQckDMMwNeefg9c/UmEYZlArnt1MKu0oIaMuO3bf3FFzrzzlvgXcctkbws/+J7+moTh+vpnu7IS6m3/Mv5ZZUlf3dssQ/Oq4MzNdFx6MF7GDUJBLzYJOzr85Z8SshE7OkxFVqetSfh/ueSy3FYaZrVV3I2J5lrU2cuxMFYNFKZRON4r5v0RX69v9ES8jhVgi1ABpizk7WoDkhRF5N1Z0BU37YN8v4rZfGQDDMNIKqL97CKiqYFUPEt9ohGH6f8IaQJCc6KgKgwsLi2ql076D1Vp5fpyU6brXFAaD8mCBnur0mw28fRRzYQSg1RnrTCCjG+UZFlUf7xoE1SCA14ukROgAEhbzceQIpuEWLmG5nOv8/XNdRUawYT8KWnIiQhz1lOQHWE8+lNoIv1tvBEGQDMHE78j39onbco66qXNMhjTd9BaGa6ICTGfcrAW+Bat4u4auHKiqoGHlu+NFDUy/d+t0zHIHDUUdh9Ab3zrf3ld4QSQ1THmwlCgvISmnaT5mcWQ2a5V+/5I6LkiVo6xB7qeqeSS6OfXfr7WWPP/3I4QqCbNEqAHQFnN29HQpLy/vhrGi27u0D/b9Ii5vuCsgzEM6BUhVflQ9vG7Uz4UV9UqHy7PywqxYxcdj+/bhrFejjv7kkQ9zYWSg0BnrTKAATXkWdR/vCsBWG+621wuvROjoWljMx1GBsuEWMWE5rT21pNPICAxjw37hoeWyN3u5NIZAwaTuEkwigYLJKzgwbfkKDobh//qq1CvjZaOD6HGBEpBg6bEb9QUwVQUEJmwvgmnSM2A69wzCqzNWDQoITFgBgQnbK2Bq9AQw7eu+kSbx04qfrjDAhBYomNRdgkkkUDB5BQemLZ8R43gN00p/VEFV5T/bBHzTHrtRXwBTVUBgwvYimCY9A6ZzzyC8OmPVoIDAhBUQmLC9AqbGr5QeHorjRmXeg+5etinS12RDGvRqpWFgND/MFAUwqbsEk0igYPIKDkxbYeG/SX6NFzz+/WzyTdYnJeFO7dIrz73f3dmE7W8kwvD8sYIQtl+pCmPFFQxMk54B07lnwHTmAqsGBQQmrIDAhO0cTI2+R0lJyZ49e8aOHaurq6uhoaGtra2trW1jY7N8+fLExEQmU5R2Oex9Pn36tGLFCktLSw0NDR0dHV1d3TFjxpw9e5ZG43FKT3+moKBg7dq1NjY2rDKpo6Pj4uJy4sSJ5ubm3jZN6ElPT1+5cqWDg4Ompqa2traGhoaBgcG4ceOOHz9eWVnZ29YJPXV1dTt37hw8eDDbxydNmnT37t3etkuU+XUZb7dhMpmvX79OTU399u1bfX29srKyhoaGjY3NmDFjlJSU+Hij/gaZTE5JSXn37l15eTmFQmF5yMiRIx0dHcXExLrOj9Exr169SkxMLCoqEhMTMzAw8PDwcHJywuFwvW1XL/D8+fPk5OTi4mJxcXEDAwMvLy8HB4feNkrU+PLly4sXLz5+/FhXV8dgMFjdiFGjRtnZ2fW2aSIFDMMvX7589uxZUVFRS0uLpqampaXl+PHjNTU1e9u0PgeTyXzx4sXz58+LioqoVKqWlpaVldX48eMHDBjQ26YJN9nZ2S9evMjPz6+vr5eSktLW1jYyMvLy8tLS0upt04QbMpkcHx+fnp5eX1+vpqZmZWXl5+fXP4trXV1dfHx8RkYGmUxWV1e3sbGZOHGimppab9slCpDJ5M2bN9+/fz8oKGj8+PF2dnZSUlIQBMEwXFxcnJycfOPGjdra2v37948aNaq3jRV6mEzmli1bbty4sWLFiokTJ+rq6kIQxGAwUlNTL1y48OrVq4iIiKFDh/a2mb0PDMN79uw5d+7csmXLJk+ebGBgAEEQk8lMT0+PjIxMSko6c+bM6NGje9tMoaSwsHDZsmWVlZXz5s1zc3MzNjZmfU6lUtPT0+/cuXP9+vW5c+du3LhRUlKyd00VUp49ezZnzpyAgICQkBATExMIguh0emJi4v79+5WVlS9fvqyoqNjbNoogfAv2tba2Hj9+/OjRo6ampiNGjDAxMVFVVW1sbCwvL09PT09KSho5cuS2bduMjIz4crv+Q0lJyZYtW5KSktzd3Z2cnLS1teXk5Kqqqn78+JGUlJSfn79s2bKlS5fKysr2tqXCR3Z29rx582RlZSdOnGhsbMxkMgsLC2/fvk2hUE6fPj1w4MDeNrDnePfuXWhoqLKysq+vL4lEYjKZX79+jYuLYzAYp06dsrGx6W0DRYH4+PidO3dSqVR3d3d7e3s1NTUcDlddXf39+/ekpKTS0tJVq1bNnz+fNZzA6A7x8fHr16/X19f39PQ0MjKSkZGprKx8+/ZtfHz8iBEjdu/ezRpIYEAQFBsbu3HjRiKR6OHhQSKRpKSkKisrMzMz4+Pj3dzcdu3ahUWmkALD8M2bN3fu3CkrK+vm5mZtba2qqtra2lpZWfnhw4f79+8bGhru3LkTG7iigMFgHDx48OjRoz4+PsOGDdPS0qqsrMzIyIiNjQ0KCvrzzz9lZGR628Yeoq2tbe/evadPn54wYYKzs7OmpmZFRUV6enpcXNycOXM2bdqENSXdIS8vz8/Pb/r06evWreukUL148WLBggVz585dtWpVT5onYjCZzKCgIAiCzp49Ky8v3z5Bamrq9OnTz5w54+7u3uPW9S3mz59fWVkZGRmprKzc/tusrKwpU6YcOnRowoQJPW+bUJOUlBQaGrpz585p06Z1NOGjrq5uw4YNnz9/jo6OVldX55kGoyOePn0aEhJy69YtW1tbrq9gGN61a9fdu3cfP37cfxrxnoMv8wM/f/5sa2sbGhr648cPngna2trOnz9PJBJPnz7Nlzv2Ey5fvozH40+dOtXRDPaKioply5aZmZnl5ub2sG3CzqNHjwgEwsOHD9t/lZiYSCQS79271/NW9QoJCQkkEunx48ftv7p79y6BQEhOTu55q0SJ2tpaX19fV1fX9PT0jtIUFxfPnz/f0tIyPz+/J20TMRgMxooVK4YOHZqdnd3+WzqdfvToUSKRmJaW1vO29TXa2toWL17822+/8Ww+Wltbw8PDSSTS27dve9424aW2ttbb29vDw4NnCWRx//59e3v7FStWtLW19aRtwg6dTvf39/f396+t5T6FtrGxcfHixS4uLo2Njb1iWw9Do9F8fX1nzJhRX1/P9RWZTJ43b56HhweFQuGZF6NLCgsLSSQSz/5he+rq6lxcXPbv3y9oq0SY7du3+/v7d74mOjc3F4/HFxcX95hVfZAjR454eXl13nAUFhYSCIRPnz71mFUiwKtXr4hEImD3+/Dhw87OzlQqVdBWiRL19fVEIvH9+/edpFm5cuWyZct6zKT+Ax+Cffn5+QQCIT4+vsuUlZWVw4cP37ZtW/dv2h84duzYkCFDSkpKukyZmJhIIBAyMjJ6wCrR4PPnz3g8Pi8vr5MEBAKhkwQiQ05ODoFA+Pr1a0cJPn78iMfjCwoKetIqUaKmpsbOzm7fvn0gm/s8fPiQQCC8efOmBwwTSdavXz9p0qTOO2Fv3rzp3P37CStWrJg6dWpra2snaV6+fInH4798+dJjVgk19fX1AwcOBBn2U6nU4ODgSZMmMRiMHjBMNFizZs2sWbM6qUg3b94cGBjYkyb1FkuXLl2wYEEnUqxevXr27Nk9aZLI0Nra6uTkdOPGDfAsdXV1lpaWPJ+YYnTJp0+fSCQSmUzuMuXFixd9fX17wKS+SUlJCR6Pr6qq6jJlbGysi4tLD5gkGtTX1xsbGyN6BrxixYqlS5cKziTRIywsbNOmTZ2noVKp5ubmOTk5PWNS/6G7wb6mpiZzc/MHDx4ApqdQKAMHDoyNje3mfUWex48fW1pa1tTUAKZ/+fKlsbExtnUoIBMnTrx69WrnaeLj44cPH94z9vQWTCZz+PDhd+7c6TzZ9evX3dzcesYkEYPBYLi5uR05cgQ8y5s3b0gkUkezpDE6ISkpyd7eHmSD6rt379rZ2XUe5xJt7t275+DgAPJo+tatW4MHD8bmoIHg4+Ozd+9e8PSBgYF//vmn4OwRJdLT0y0sLFpaWjpJw2AwRo8efevWrR6zqld48eKFra1t50cW0Ol0Z2fn+/e7fbhH/+P48eNBQUFIc6WlpdnY2GD1JArmzJlz7tw5kJQMBsPOzi4zM1PQJvVNli9ffujQIcDEw4cPT0lJEaQ5osPGjRs3b96MKAuNRjM1NcXCUoDU1dUZGBiAxDQuXbqEovrF6JzuBvvWrFmzYcMGRFkKCgoMDQ0bGhq6eWsRhkajmZiYIG3PDhw4MHfuXAGZJEpkZGTY2tqCzKdwc3MDD2QLIw8ePHB3dwdJOWzYsKdPnwraHtHjxIkTU6ZMQZrr9OnTKHL1c5hMpq2tLfiz2eDg4DNnzgjUpD4Lg8GwsLDIysoCTD916tTIyEiBmiQCREVFeXh4IDqds6mpCRswAOLj4xMXF9dlsoyMDDs7O9E+I3XMmDEgPZOXL186ODj0gD2iBJ1OJxKJ379/R5E3MDDwypUrfDdJtGlubtbR0QFfERkREbFo0SKBmtQ3aWtr09PTAx87R0dHz5gxQ6AmiQZkMllfXx9kYikXV69enTZtmiBMEj0uXboE6LYtLS26urr9ZDuOHqNbwb7y8nJDQ0MUf8natWu3b9/enVuLNqdOnQoJCUGai0ajGRsbY9s0dElwcHBUVBRIyoSEhMmTJwvanl4kICDg9u3bICmjo6Oxhy1IaWlpIRAIRUVFSDMyGAx7e/vU1FRBWCWqPHr0yNPTEzz9ly9fTE1NBWdPX+bOnTuIFkN9/PjRyspKcPaIAAwGA93yk6ioKD8/P0GYJEqUl5eTSCTAJc+jR49+/fq1oE3qLYqKiszMzACjmUOHDn337p2gTRIlbt++jbrX9+zZM8CnpxhsEhISAgICwNPX1dXp6+v3w90Pnj17hqiH09zcrKuri+0r1yVnz55dsmQJiow0Gk1PTw9FlLAfEhgY2OUaMjZBQUEiPz2/hxHrzuEerMmWCgoKSDOuWrXqwoULMJ8OAhY9Tp8+jeJgLykpqaVLl164cEEQJokMMAw/fvzY29sbJPHYsWNfv37d3NwsaKt6BQaD8eLFC8CjzcaNG5ecnEyn0wVtlSgRGxs7atQoQ0NDpBnFxMTWrl17+vRpQVglqsTHxwcGBoKnNzY21tLSysjIEJxJfRakWllaWsrLy2dnZwvOJGEnJSWFSCRaWVkhzTh16tSsrKyKigpBWCUysEL5YmJAXdaAgICYmBhBm9RbPHz40Nvbu6PDIrkQbSkEQXx8fEBAALq8v/3229evX2tqavhrkmhz9+5dHx8f8PQqKiqGhob5+fmCM6lvglQoWVlZW1vbd+/eCc4k0eDu3buTJk1CkVFKSsrDw+Phw4d8N0n0ePbsmaura+dpqqqqCgsLi4qKBg8enJycXPd/sGBR9+lWsC8mJmb69OkoMmpra+Px+LS0tO7cXVQpKChgMpkoxgwQBE2ZMiU6OprvJokSWVlZJBJJSUkJJLGkpOTgwYMzMzMFbVWvkJaWZm1tLScnB5JYVlbW2to6KytL0FaJEjdu3AgKCkKX19fX9+HDhzQajb8miTAPHjzw9PRElGXixIn37t0TkD19FhiGHz16NHbsWES5+qdW4MTExEybNg1FRjExMT8/v9jYWL6bJEog8m5XV9fnz58L1J5eBJNCoDx9+tTFxQVdXhwO5+zsnJ6ezl+TRJsPHz44ODiApAwKCjIyMjI3N//x48eUKVPc3d3HjBkTHh4uaAv7COBCLVy40MjIyMzMLD8/f/bs2Syhdu7cKWgLhZS0tLRhw4YhylJUVPT+/fvi4uIhQ4akpqYKyDCRoaqqSllZucvBZlRUlKWl5cCBA7dt2xYVFUUgEPT19e3t7XvGSNFGAnVOCoXy48cPS0tLdNnd3NxevHjh5OSE2gBRBST+3RF6enpSUlI/f/7U0dHhr1Uiw7t374YOHdplsqlTp75580ZcXJzJZC5cuFBZWVleXn7EiBFbtmzpASN7BkApJk6c+P79ezExMSaTOXv2bBUVFTk5OXd39/Xr1/eAkcILDMNpaWm//fYbuuzy8vKWlpbZ2dlDhgzhr2EiSXFxsby8vJaWFkji0tLS79+/q6qqkkiks2fPCtq2vsbXr1+1tLTU1NRAEpeUlBQXF6uoqBgbG1+7dk3Qtgkvz54927x5M7q8bm5uUVFRixYt4q9JokRWVtbhw4c7T5OXl+fp6SkhIYHD4Wg02uDBg5WVlWVkZPbu3WtjY9MzdvYAWVlZXTbcWVlZvr6+EhISYmJiLClUVFRkZGQOHTpkZmbWM3YKI83NzS0tLRoaGkgzhoaG3r17V01NDYfDHThw4Pnz5zo6OgQCYfz48YKwU5QoKioCXP1gbm5+48YN9vqSvLw8aWlp1NMwhY7v378TCASQlNbW1pGRkVQqlfX206dPUlJSSB/v9RNaWlrExcVlZGQQ5crIyAgKCpKVlaXT6XQ6/ezZswoKCng8PiUlBeml+gPfv38nEoldJgsICNiyZUtdXR37E0lJyenTpwNOY8foBPTBvk+fPllYWCD9D+Li4iIjI42NjSkUSn5+vru7u56e3oABA1CbIXp8/PjR1tYWaa7ff/+9trbWxMREWVn50qVLEydONDQ0lJWVFYSFQk1JSYmBgUGXyZycnG7fvs05r0paWnrixImCNK2nKS0tJZFIXSZzcHC4f/9+a2sr+xNpaWl0c1j6FSUlJdra2iga/r/++isnJ8fc3FxSUvLatWusPgTmy51TWloK0plgUVRU5OrqqqSkxGAw6HS6vLy8uLi4pqZmWloaYAhMqEGk1devXz09PRUVFTm10tbWTk9PV1ZWFqidQgSdTq+pqdHT00OU6/3791u2bNHX19fW1n7+/HlaWpqBgYG2tjbWtW1PVVVVlx1Fc3Pzpqam2tpa1tvS0lIIguTl5U1NTQVuX08BwzCFQlFUVOw8ma2tLZlMJpPJrLclJSUQBCkpKYG0+P2Z0tJSPB6PIqOjo+PVq1d//vwJQVBOTk5ycrKEhMTcuXOxYF/ntLa2wjAM2L2ZOnXqwYMH6+vr2Z9ISUkh2o9CqKmpqQEMQ0+ZMmXTpk3sYB8EQdLS0rNnzxaUZcJMUVERYAiVEzc3NwkJCXZYikajNTU1OTo6YpE+nhQXF4ME9HV0dAwNDXNyctifyMvLz5kzR5Cm9RfQL+MtKioCHzCwUVdXf/To0YEDB86fP//27VtXV1cikainp4ftBcYGXdVTW1t79erVrVu3fvjwYe/evcOGDVNRUVmzZo0ADBRuSktLQYJ9kydPlpaW5vxESkpqxowZArOrFwCMewYEBHDNvpaUlOw/T1NRA/gsqz0MBiM2NnbHjh1Pnz69cOHC8OHDVVVV+0+PFh0/fvwAD7U4OTnJyMjU1NTU19dTKJTm5ubGxkY9Pb3+EOmDEGo1fPhwKSkpLq2IRCIW6ePkx48f+vr6SHMRCIQHDx6cOHFi27Zt9fX1Xl5eFhYWMjIy379/F4CNQgyZTFZWVu4yBorD4fz8/DiTiYmJBQYGcjXlQk1lZSXIgF9MTGzcuHGcUoiLi8+aNUtSUlKQ1gk99fX1qqqqKDK6urpyFTN5efl169bxyS6RBSRyzcbExIRzj3gcDufh4QGeXahhMBgSEqCzc7S1tXV1dTk/cXR0RDFftT9QVlbGpRUIKioqXC2+srIy6qn9Ik9LSwvghlHz58/nTKmtrS1Kz+oAIZPJhYWF2dnZnz9/Li0t5Zxqgxr0wb6GhgYU3f0hQ4aw+h+tra1kMrm2tratrW3hwoVYF4QNOmE9PT1ZrSCFQqmrq6utrZWRkVm4cKEADBRuAOXF4/Fc8whcXFxUVFQEZlcvAChF+w6Wt7c3imN5+hsNDQ3oCsyoUaNY/0tzczPbl5cuXcpvA0WKyspKTU1NwMQSEhJce9+oqKiI0gr9zkGklZSUlJ2dHecnKioqf/zxhwDsEmIaGxtRtNrKysra2toQBLW1tTU0NNTW1jY0NLi4uKB42ifagHeKZsyYwVnrKikpiVjNCS7FzJkzOVMqKipi68S7RFJSsq2tDUVGIyMjzlgMa/M+bB5ll4iLiyMSfNq0aeLi4qzXKioqv//+u2Ds6nOIiYkxGAzw9MHBwezos4qKysqVKwVjl9CD2uUnTZrE6fK6urqAOypidEJAQADbwWVlZRcvXty79vQMdDo9MTFx48aNI0aM0NbWdnR0nDp16qpVq+bNm+fr64vH401NTQMCAk6ePFlYWIjuFuiDfaw9vJDmkpWV5QqHS0lJ9Z/6GgQxMTEUR88MHz6cM2AqJibm6upqbGzMV9NEAUlJScBppEFBQWxJVVVVV6xYIUi7egFwKQIDAzk7WMuXLxekXSICuhoSgqChQ4dyPckxMjJCvfdfP0FCQgJRj41rvqqCggLqnVKFDgkJCURT6QMCAjiXWSkpKY0YMUIAdgkxqJ191AZJfr8AACAASURBVKhRnG+VlJR27drFJ6NEB3DvHjlyJOcfoa6uLmLbe4NL4erqyplST0/PwsJCYHaJCEpKSg0NDejycvqysrJyWFgYn4wSZcTFxRHFsAIDA9nH68nIyPSffhEOh0M0MORstVlTIAVjl9AjJydHoVBQZBw3bhy7KCoqKv755598taufoqWlxQ5cSEpKivyGUQUFBcuXLyeRSCdOnDAwMDh16lRJScmnT5/S09MTExOfPXuWmZn58+fP1NTUefPmFRUV+fj4/PbbbxcvXkS6HBZ9sG/AgAHV1dUoMrq5ubEXF0hLS7NOP0BthugxYMCAqqoqpLmIRCJnS6CgoLB9+3a+2iUiSElJAc6JDQgIkJeXZ72WlJTkGpWJAOBScHawZGVlQY71wNDQ0EDhyBAEycvLcy64UFFR2b9/P//sEk0UFBSamprA03t7e0tJSbFey8nJrV+/vv9slIZUKy8vL/YcAQUFBWxaX3tQd4fGjh3LOUvaxsZm0KBB/LNLRABvqsTFxd3c3Fiv5eTkRO+5FLgU0tLS7KC8vLy86D2tFATdCfb5+vqyfVlPTw/rJoEgLS3NubVclwwaNIg1nYo1eOw/rTYEQTAMg8f7jI2NWQucJSQkZs+eDb4EuL8hLy+PLtjn5OTEroolJSUnT57MV7tECkVFRfYGsl0yf/581tDb3t5ehE90KCkpmTlzpqenp4mJyYcPH27durV48WIrKyue61zV1NTGjh27Z8+ejx8/hoeHp6SkWFhYnD9/HvwZM3r/RxeTgiDIzc3t6tWrrD9eWloa21eOCw0NDXTDhkGDBj1+/BiCIBwO5+joaG1tzW/TRAE5OTnAga6VlRVrs1UpKakFCxaIiaGPjPdNwKUYNGgQ6+ez1pP2qw4WalAH+yAIGjVqFHvfLh0dnf4z6Qw18vLyiAJY+vr6cnJyrK2+JSQkgoODBWZanwNpsM/ExITd/xAXFxexrUv5grq6ek1NDYqMnFPylZWV9+zZw1e7RATwpgqCoODg4OTkZDKZLC4uHhQUJFDDeh5EUsyaNevly5cNDQ1iYmJTp04VqGGigYqKSk1NDQzDXJ2cwsLC3Nzcr1+/Njc3U6lUeXl5XV1dU1NTW1tb9vwpNzc31gIIJSWlbdu29YL1QoiEhISUlFRjYyN7673i4uKPHz9++fKFtUWsgoKClpaWmZmZra0tKwTg5+d37tw5KSmp/rZzv6amZkVFBWvnBwiCSktLP378+PnzZwqF0tTUpKioqKGhYWpqamdnxxJzxowZBw4ckJeXnz9/fq8a3qfR09MrKiri+pBGo3348OHz588/fvygUChSUlJycnIkEsnKyoo970xCQmLIkCFPnjyRlpZesmQJthdZJxgaGrIOiWLR3sc1NTXNzMzs7Ozk5eWnTJmydu1aRUVFEduCgw0Mw4cPHz58+PDmzZsvXLjAGYhvaWlhFbyysjIKhSIjIyMnJ2dkZGRlZcXeFMLBwSEyMvLbt29r1qw5f/78+fPnQebsow/2aWpqsk6e4qSqqurNmzf5+fkVFRUNDQ3S0tIqKipGRkb29vY2NjaseMHw4cNZiSUlJYOCgrBNQ7ngKSyL5ubm2tra1tZWOTk5dXV1rspl3LhxL168aG1tVVZW3r17d48YK3wQCIT2i97r6+tZR/ipqqpybs/s7+9//PhxGRmZefPm9aiVPQJPKerq6lgnTKmpqbE3P8LhcBMnToyIiJCUlOxvHSzUaGhodOTIEASxtuhiMplKSkpqampcoWR3d/eYmJimpiZlZeUDBw4I3lihR19fv32PDYIgCoVSW1tLp9MVFBTU1NQ4m1Vvb+9z586xhsHsiav9gY60Yp1k2tbW1l4rT0/Py5cvS0pKBgcHA2603K+QkJBgMpl0Or19jx+G4draWjKZjMPh1NTUuNYxcJ7hY2RkxO4dYXAiKyvLYDCoVCrnWYc0Gq22tralpUVaWlpNTY0dcxk7dixrYeCoUaPQHbbQl1FVVW1oaGAwGOyNNaBfpVBXV2er5O3tzWQycTgc60DtXjJZmJCUlDQwMCgoKGAN6V+8eHH58uWHDx+qqanZ2dmRSCRlZeUBAwbU1dW9e/cuKirq/fv35ubmgYGBAQEBurq6srKyZDJZVlbW19e3t3+K0EAikb59+0aj0S5dunTv3j0FBQVbW1sTExNlZWV1dfX6+vqPHz/GxMS8ffvWyMjI39/f29v7woULNjY2IOfLiRJGRkYFBQXl5eWRkZF3795lbadramqqrKyspqZGJpPz8/Pj4+MzMzMNDAz8/f3d3NzCw8MNDAzMzc172/a+i6KiorS0dHV19YABAxobG6Ojo69du5adnW1nZ2dubq6jo6Orq8s6YyA+Pn7nzp0/f/4cM2bMzJkzR48e7e/v/+bNGwkJiSVLlvT27+jTGBoaFhUVZWRksHxcTk7Ozs6O08dzc3Nv3bqVmZlJIpH8/f1NTEzy8/N9fHx623D+U11dHRwcLCsrm5mZyZ63WF9ff/369evXr+fn59vZ2VlYWGhpaenr61OpVDKZHBsbGxYWVltb6+HhMXPmTNbeBUQiMSYm5s6dO+PGjfvrr7+6nK+APthnaGhYWVlJoVDk5eWrq6svXbp05cqVrKysjqYZq6urT5w4MSQkZNiwYaw1QbKystjhNe0ZNGjQ2bNn2W/Lysqio6OTk5PfvXvHZDK1tbVxOFxra2tlZaWSkpKzs/PYsWN9fX3l5eVHjBihoKBQW1trbm4+ZMiQXvwJfRkTE5P4+HgIgl6/fh0fH//06dOvX7+qqakpKCjgcLimpqbq6mojI6NRo0b5+PgEBgaeOnXK0tIS5NRwocPExOTZs2cwDL969So+Pv7Zs2cFBQXq6uqsx6eskb+xsfGoUaMmTJgwY8aMS5cuDRkyhP1cEaNz5OXlFRUVOY+cz87Ojo2NffLkSV5enoKCgoqKCg6Ha25urqysNDAwGDFihLe3N2t2wLBhw1ihFnV1dS8vr179HcKBlZVVdnY26/X3799v3rz55MmTrKwsCQkJDQ0NHA5Ho9EqKioGDBgwbNgwDw+P8ePH+/n53bx5E4bh1atX967xPYytre2HDx9YrwsLC6Ojo1NSUt6/fy8lJTVgwAC2VhoaGsOHD/f09Bw3btzkyZNv374NQRC2z3dH2Nravn371snJCYIgJpP55MmThISE58+fFxUVaWhoyMnJwTDc0NBAJpNNTU1dXFwmTpzIWrE7aNCg5ORkZWXlffv29faP6LuwIgJ6enqxsbFJSUmpqanNzc2ampqSkpJMJrOiokJcXHzIkCHu7u6TJ092dHR8/vy5qC5cNTAwKC4uVlVVjYmJSU5OfvPmTUtLi5aWloSEBIPBqKiokJSUdHBwYElha2ubkZEhesuZBYejo2NaWlpubu62bdsUFBTmzJmzZ8+ejqLGTCYzPT396tWr27dvnz59+siRI2NiYjZs2MAZisXoHGVl5enTpysqKoaEhISFhXW0cA+G4aysrCtXrhw4cEBMTCwkJKSH7ex1lJSUFi1aJCUlFRIS8vz5cy0trY5SZmdnR0VFzZs3j3UGd08aKYzY29u/fPkyIyPj4sWL48ePDwsLGzp0aEcu3NDQkJCQsHv37hUrVoSGhtJotHHjxnXyX2BAEJSVldXS0rJ06dKQkJA///yzozleMAy/f//+6tWrRUVFBgYG9fX1IjbeLCoqGjduXGhoKLtzUltbu3fv3qioKD8/v927dzs4OHS0iLCuri4+Pv6PP/5oamoKCwtjRUJ9fHwGDhw4ZcqUgoKCv/76q7N7w93Az8/vzp07a9asYUUHcDicvb3977//fvr06fv37z958iQxMfH69es7duwYP348e+rEb7/99ttvv+FwuKCgoO7cXVRpbGzU19eHYfjNmze+vr5EInHjxo0pKSkUCoUrZWVlZWxsbEhIiJ6e3u+///7t2zcFBQVlZWVWBAeDJx8+fHBycrK0tBw5cuSxY8dyc3Pb2to4EzAYjLy8vJMnT7q6upqZmcnLy0dFRfWWtQIlLS1t2LBhZmZmLi4uJ0+ezMvLYzAYnAna2tpyc3OPHz8+atQoc3NzWVnZW7du9Za1wsjixYsvX77MZDKjo6MdHBwGDx68f//+t2/ftra2cqUsLCy8cOGCj48PHo/ft28f63BPBQWF6OjoXrFcGLG0tIyOjh47dqyZmVlYWNiLFy9aWlq40pSVlV2/fj0oKEhPT2/VqlWSkpJDhw7tFWt7F2Nj45iYGDc3NwsLi23btr169YpKpXKl+fHjx7Vr16ZNm6anp7d27VpxcfFRo0b1hrHCwfbt2/fv39/c3Hz48GESieTl5XXu3LkvX74wmUzOZHQ6/cOHD+Hh4c7OznZ2dlevXt2/f7+4uLi5uXlvWS4UzJo1y8fHR09Pb9GiRfHx8bW1tVwJGhsbExMT16xZg8fjnZyc1NXVuZozkSEwMNDHx0dfX3/p0qUJCQl1dXVcCRoaGh4+fLhy5UpDQ0NHR0dtbW2uQojRCQcPHiSRSC4uLqmpqeC56uvr//rrLw0NDWlp6cbGRsGZJ0qUl5dPmjRp6NChT548Ac/V1NQ0evRoAwODkydPCs62PkVNTc2MGTMGDRr06NEj8FzNzc2enp76+vqHDh3CaoBOWLlypY6OzubNm1lL+AF5+/atp6enjIzM7du3BWebsFNRUeHv7+/k5JSSkgKeq7S0dOnSpYaGhsePHxeZovv9+3djY2POMV1kZCQej9+2bVt9fT34dVJTU0ePHu3l5VVcXMz6pKmpycvLa9WqVZ3k6lawLzQ0lLU0QEVFZePGjV+/fu0kcUtLyz///MN67o3D4XA4XHZ2dnfuLsLY29v7+vra2dnFxMSAdFibmpoOHjyIx+O1tbUtLS17wEIh5cWLFzY2NkFBQbm5uSDpP3/+bGtra25ujqgvIhQ8efLE0tJy9uzZnz59Akmfm5trZWVlbW394sULQdsmMly/fn3y5MnDhw8fP378mzdvQLKUlpauWLGCSCRaW1vr6emJ6niV75SWltrZ2Q0aNOjevXsg6evq6rZv3y4tLb1o0SKR6UwAUlRUZG1t7eDgADhyqKmp2bp1q7S09PLly/ubVuA8efLEycnJxMRk8eLF379/B8ny7t27SZMmsfZbwQYMHUGlUv/44w8ikXj8+PH24fv2tLW1nT17VlNTc+HChWQyuQcs7DGam5vXrVtnZGR05syZ9tH59tDp9GPHjmloaCxbtgyLQIFw9+5dPB5/6dIldNnfvXuHx+ODg4NBCmo/58mTJ0Qi8cSJEyg6OVQqtayszNvbe9KkSSJfsFNTU42NjQ8ePMg1LwEEKpVaVVU1ZcoUT0/P9g9IMGg02oIFC5ydnTsPX3TC8uXL8Xh8TEwMfw0TDZ49e0YikY4ePYpuIFNeXj5hwgRfX18RaMerqqosLS3Z5aS5uTkoKMjFxaWoqAjF1ZhM5sWLFwkEwv3791mfUKlUV1fXPXv2dJQFfbBv69atOBxOTExsyZIl7R8tdsKtW7d0dHQgCLK1tWUHJjHYvHz5Ul9ff+fOnXQ6HVHGmpoaa2tra2vrsrIyAdkmvDCZzB07dpibmyONVVGp1Ddv3lhZWW3ZskU0BrpMJnPLli1WVlaA4Sc2VCr1xYsX5ubmO3bsEA0pBM3Ro0f19PTu3LmDNGNeXh4ej3dxcaHRaIIwTMR48OABgUA4ffo00mJ54cIFX1/fcePG9Z9+cEJCAh6Pj4iIQKrVuXPnxo8fP2HCBEQPIfsJdDp96dKl5ubmHz58QJr33r170tLSR44cEYRhws63b9+GDBmyZMmShoYGRBkbGhq2bt1qZmb27t07AdnWw3z58sXe3n7lypVNTU2IMtbX12/evNnS0hJ7vt45Fy5cMDMzA3z82RFkMnn16tWjRo3C6slOiImJIZFIrH2fUMNkMrdt2+bg4FBZWckvw/oa9+/fJxAIiCaZ8uTgwYO2trY/fvzgi1WiAYVC8fT0XLhwYfulNoj49u2btbX1sWPH+GWYaBAXF0ckEt++fdudizCZzJ07dw4aNKiiooJfhvU8DAbDw8Pj77//Zr2tr68fOXLkmjVrUITvOcnPzzczM7t48SLrLZlMtra2fvjwIc/EaIJ9TCaTdUiKhoZGYmIiiivU1NSMGTMGgiADA4MvX76guAIgzIZPt/fOcXLbmSe4e/CV27dvk0ikjIwMdNm/fPly/vx5U1PTz58/89ewXym/t3GsqZq8os7A6cczkXU8e4O2tra5c+f6+Pigfj7Q2Ng4ZcqUGTNmdLNV6ApG4Y0FjtoKitpD5l77Kog5XXQ6PSgoaPLkyagfh5LJ5AkTJsyZM6eb9RTborKnR4LtvY6UdJFO6Hw5LCzM2dn558+f6LJ//vx5yZIlbm5ufHpwDaizkLk2DMORkZFmZmaAc3Xbw2QyDxw4YGNjU1payg9zAHVG5OmgPtIl586ds7CwQN06MJnMPXv22Nvboy7VoDcSKmdvaWnx8fEJCQkBmWzFk/T0dBcXl5UrVwr4IYqQefeHDx9IJBJrY010PH36lEAgJCUl8dGqdgi81YZhODMzk0gkdmf6Z1JSEoFA6Mk9XoTLi6OiomxsbPhVs+3atWvkyJE9Nb9PyPz67t27pqamgNOfu+TkyZODBw9G+jAALT3h7GxSUlJIJFI3o89sLl++bG1tjWilKlKEyOVbW1s9PT03btzIl6vV1NQMHjz47NmzfLkaL4TMxx88eGBiYlJYWMiXq505c2bgwIGCfHwiWL/et2+fv78/63Vzc/PIkSN3797Nlyv//PnT2tr6n3/+Yb3NyckhEAg8A6Nogn1//PEHBEF4PL47ESUGg7Fs2TIIgkgkUnl5OerrdM7PU5P0tWQhq7Auqh5GXc6t/fNGh1xFXJYqzgW7LTp0J5/Mhz76o0ePjI2Nu+8eCQkJxsbGJSXdHhp2ADX577/icmubyXnXZpHEiBvSO07aN4RdsGDB9OnTkc6U5ILBYMydO3fmzJkCHI8V3wo/n1FNqX132FNFNTSZ35dnMpmzZ8/ufpyOFTEMDQ3lhxQZm8xlxCCn8K5Kq3D58q5du1xdXZFOwWhPWFiYu7s7P+b3AeksdK5948YNKyur7o/Qrly5Ym1tXV1d3W2LwMozMk8H9ZHOuXLliq2tbfcnQVy4cMHW1hbRdH6kCJGzt7W1+fr6rlq1qpuVIZVKnTBhwrp167ppT2e3APTuPqAqDMNfvnwhEond30Dj48ePRCLx1atX3TeJNwJutWEYzs3NJRAI3f8J79+/JxAI6emdVOv8BMiL+0ZhS01NNTIy4tPznn9ZvXr1zJkz+XjBjhAuv2YNR/k7yWPHjh0+Pj49sdxE8M7OpqCgAI/H83c27t9//+3i4sKn5/Q8ECKXX7hw4cKFC7t/HTYVFRVmZmaIdqYDR7h65nl5eQQCgV9BahZ79uwZN26coPY1EqRfFxUVGRoaVlVVsd5Onz6dv928kpISIyOjtLQ01tsDBw7Mnj27fTLEwb7bt2/jcDgtLa3uV9ZMJnPevHkQBLm5uQlwa6qUJQM6rnqYDZ/vHl46xtJ0ePC2fzIrUUSDWiszr/4V5GRqNW75scQC9GP7r1+/4vH4vDz+PBG5evWqg4OD4NcAVh0dBY05zeNJUd8R9siRIx4eHnyZkcca1/ErJM/7BqwXn3YOctxXwO/L792718fHhy8tPZ1O9/T0DA8P7/6l4Lb46bJggQwh8eW4uLiBAwfy61HzkiVLFi1axIcLgesMw0Lh2u/evSMSifzaC+Lvv/92d3fnQ0sEojNST0f23/EgLS3NyMiIX6t4Dhw44OXlJdjBlZA4+7p16/j1BKi5uXnYsGGRkZHdv1RX8PbuvqNqY2OjlZVVQkIC6itwkpmZSSAQBLWETcCtdn19vZmZGaKN+Tvh9evXJBKp55ZEdezFfaewNTQ0mJiYvH79GvUVeNLW1ubq6speY9Uj9HW/plKpNjY2/CrMbJhM5sSJE/nTHe0cATs7Gzqd7ujoKIid4IKDg//66y++X/Y/hMHlo6Ojhw4dyvd1WllZWSQSiR/PjDuhr/fMaTSanZ0d4N7ZiJgyZcr+/fv5flkYFqxfBwQEnD9/nvU6IiKCP6OMX3n58qWpqSlryEmn021sbNo/0kMW7KuurtbQ0BAXF+dX9Lq1tdXZ2RmCIIQb1rQU3tsz1So4/NYKR015NXPfo++aO0z7erXeL1VPzbWppoH/VFO+JZ9Y6WNrPix4W1RmBafPU/OjFo8kqSlqWE898Sa3/Y2a0/aMJSrruB28eTWO/dCFVvYy8g//ISb2E9eee1KEdO4+g8FwcnKKjY1FmK8zli5dijx+jERYGIbhqmOj9Va9YrJVhWG4A2E5VX3fwvNGAhA2OzvbyMiIj5VvfX29qakpwmfjiFRtrf54d+dEn93ZrLqaJWwpwv+FB5mZmSYmJnycj1NdXW1kZIR8jyqukgDDcOJcZae9iRGBFqqqxHHHcjquB3/x5T5a5H7+/InH4/n44Lq1tdXJyQn5Aq5u6AzDfd+1qVSqpaUlf5+jBgYG7tu3D3k+LgUAdebydJhD6m7+d9xQKBQzM7OXL18i/2kdMnny5MOHDyPMhLrh7qOFMDk52dbWtrkZcW3cEaWlpXg8HuHsfqStNtzeu/tUdwiG4dDQ0D///BN5vg6JiIgYM2YMwkzgwvJstatR/TXczJw5c+fOnUhzdcKJEyfGjx+PPB9nScgsfHZ0zmD30E1+JFWTDS+RNdnoCxsM8yhv3S5sK1eu5G9hY1NcXIzH45H3P1FJDcP89GtYIFJv27Zt6dKlyPN1TXV1NR6PR/7kD4XUHTk7VxtEL0Pwx3Fz+PDh4OBghL8FiIaGBjRLg29OxUFsJIbsye8wZZ93eTKZTCQSBbSB2J49e+bPn48sD7i2MNz3e+a7du3i75RJNrW1tXg8HvHyf1B5O/BrZP8ON/n5+SYmJqxlhZWVlXg8nr/zx9ls3rx59erVrNc3b9709fXlSoAs2Ldo0SIIgtasWcMf62AYhuGvX7/KyMgoKyuzZzkCcGeOshQE2S27n19PLrg6zUDCcV+HvWPuYB+Ln+HDxXS9wt+23w4rd5eH/5ncuqbyl5sHSUMQxH2jgrggmxVPaupLXx9aeYR7hjU5dd9YTZzbWYQRlRMnTkydOhVZnq6gUCjGxsYfP35EkgmJsDAM/zg9cdL5X1fQdSDsL6oOPfyNx40EIeyIESP4fsphYmKig4MDkqkcSFRtuewNQRAEydqse0FBdQVeMJlMR0fHBw8eIMnUNfHx8SNGjECWh7skwDCcOFeZ6H/yWWn9j4t+KlZhHS9bQOTLvVTkZsyYwd6HlV/k5OQYGxsjiyl0R2dYCFx727ZtS5YsQZanK1gDBsQbIHArAKYzb0/neUEY2X/Xjs2bN69atQrZj+qKiooKQ0NDhIdBCazh7o1C2NraamFhgfSYoy65ePFi+45apyBvHbi9u291h169emVtbc33aRdjx45lb20DBrCwHfpydxvux48fDxo0iL/L7phM5qhRo+Li4pBl+6Uk4KQlZMQg/eBbP7uwjIcXoy9shTDM7LC8oSxs3759I5FIfIzXc7Fr166VK1ciy4NOapiffi0IqSsrKw0MDAS3/8PJkyd5rl/rDBRSd+Ts3G1Q6jpj4D/uV8hksoGBgeC2t7py5crkyZMRZWm9cehQAZUJw2W3gh0mnCvs5Cf1eZffsmULv7bqaw+NRkO6hTQCbeG+3jOvrq42MDAQ3L6QZ86cQbo3Aqi8Hfg1sn+nHQsXLjxz5gzr9YoVKwS3LpBCoRCJRFYklMFgtN+bG0Gwr6SkRFpaWktLi+8nnW/atAmCoM2bNyPIU31qNOR7jRWBzdlqLuZ3vaOeIe8xA9xUmHR8hbe1ieOMsCvp5f9lLj38239hXEhsyvG/uW+Ud8oVbzVp7dF7nzh0oJW9vhw2baj5QN81Z1OQxcKpVCqBQEB3AHPnXLlyBXEMEVxYuDJ6zaq4du0RT2G5Vf2HxuNG/Bb24cOHLi4uSHKA4u3tfevWLQQZEKgKt5JL0m+ucVaSmHyN4+ciuUJ7bt++7eXlhSADMK6urohiiDxKApw4V/nfJYoftlgS1nR87hgSX+6VIpebm2tubi6ILVEWLVqEaPpzt3Tu865dX1+vr68viC7F4cOHFy9ejChLOwXuAurM29O7+99xw+p+oT6bqBP27du3YsUKZHkE03D3SiGMiIgICgpCkgMIJpM5aNAg9g4sQCBrHXh4d9/pDsEwPGbMGEEsAsrJybGwsEBWOQML25Evd7PhHjFixOPHjxFkACMzM9PW1hbR2nPukhAVM0MOYGMBXl6MvrC1wvDXduWte4VtyZIlfH84x0lDQ4OBgQGiyX0opeavXwtA6j/++EOgC0hpNBqJREI0kkInNU9n59EGMe6A/XHc7Nu3jz1JRxC0tbWZm5uj2DOqKS1sxMg/0jpf39m3Xb6xsRGpPyIlIiICxaxMIG37fM9869atApolzaK1tdXIyOjbt29IM4LI22EjDvrvcNPS0qKnp8daXcvqh/M9esbJkSNH2BMgdu/evWXLFs5vxSBgzpw5Q6PRli9frqCgAJ4LhFWrVsnLy586daq1tRU0j7y8PATDMARBECQtLS2voiKJ7KbyRLfF4QkfMq9MU3m90cV6xJxdMe9rGBDEYDDstrNrKsbNEFXuG5kveJQZNV3jw3a3gYsf1NPKUi9tmeFoNWZb+oA5VzPexu2fN9pQBoklN27cGD16tKGhIbIfAEBgYGBGRkZxcTGCPKDCUt9HnqfP3eWr1e4CvIRt5VJ1qhSPG/Fb2PDw8I0bNyLJAcrGjRsPHTqEIAOS4iqppD9kyv6/AgeIi3O4Z/cKfN+Rgtu/pkpxfovD4f7/IxHQd4rckSNHVq9eLS4ujvQndMm670p3LwAAIABJREFUdetYYxLA9N3QWQhcOyIiIjAwUE1NDUkmIEJDQ+Pj4+vr68GzdCJ15+WZt6fz20fOnj0bHByspKQEngWQRYsW3bx5s7GxEUEewTTcvVIIDx8+vGHDBmTWA4DD4davXy+w9oW3d/ed7lBOTk5VVZWXlxeSTEBYWVmZmpreu3cPQR5gYTvy5e4U+IyMjNbWVhcXFwQGgzFo0CAdHZ2kpCTwLNwlYZqSDELPZYO+sElCkNF/5W1uxINuFraWlpa4uLg5c+ag/CUAKCoqBgQEXLt2DTwLKqn57df8lprBYFy+fHnhwoVIMiFDSkoqJCQkMjISkVUoSjVPZ+fRaovJoPARGIbPnz+/dOlSxDmBERcXX7RoUUREBKJczO+XZy4rXHNju4M84jv2nXIYExPj5eWlrq6O+DcAM2PGjOTkZDKZDJ4FTNu+3jNnMpmRkZGsBaACQlJSct68eRcvXkSUC7DodtSIoy75iYmJI0aMUFRUhCAoKioqICCA79EzTkJCQuLi4qhUKgRBM2bMuHHjxi9fw8AQiUQJCQnAqcVtRbGLh+rIKxLGhqXUAqSfOXMmBEEIVly2XPaG8LP/ya9pKI6fb6Y7O6GDeGkbOXamisGiFEqnOyYw6rJj980dNfdqPczMCjPXdN2e+KW+sTo/4fDll5FcN6q7+cf8a5kldXVvtwzBr447M9N14cH47pxf4+7u/uzZM5CUKM41DwsLQzZxFEhYas4/B69/pMIwzKBWPLuZ1OEa9P8Le+Xpr6pmMtrfiL/ClpWVkUgksI0w0Rxqbm5ujuDxAmhxZUEn59+cM2JWAudpSsiu8Avfv383MzMDe4yP+AByBoNBIpHAV/Nx+1cmo7Xy/Dgp03WvKQwG5cECPdXpNxt43xiZL/d8kaNSqezHOIC05e1zNFoOuDE4eEUBo9dZCFwbhmF7e3vghRLl5zxZ7auM5zmgQ3uXLFly9uxZcGO4pI54cBasPMO8Pb1bPsIDS0vLr1+/gqZuK03cO3ukyQBp260grUxoaCiyAyUE03D3fCF8+/bt0KFDgZK2vfqd80GeXtf+zqpJEEzGBFUVzLt7tTu0Zs0awJlWP+5vGGOiJq+o6zDz/Eewo8hu377t7++PwBoEzS5vX+5Ow71s2TL2CqDOKX+02Z2kLKdu4rnpYTmY9tevX0e0HoqrJJy/f85beuDeLnYyB/BiJIWtEYbp/5U3AwiS62Zhi4mJAZmc277jjWiA8/bt2+HDh4NbhVxq/vs136V+/Pixt7d3l8naS41o1PPt2zcrKytwq1CVapins7dvtWkVEWBX+4WMjAyQvXHayYJsFFNRUUEgEBCYVZfy+3DPY7mtMMxsrbobEdtRx7Cvu7yHh0eXvWheRQ7ZyGjBggWXLl0CtQlIWyHomT979szDw6PLZB15NOBoqKioyNzcHIFZoEUX5t2II8jOzfLlyyMiIlivnZ2ds7KyOk/fgTIIhjDTp09nLzq0sLDg3JIINNj36dMnCIJAF0Uy83cOVhsTnlNb8XytvW7IPUqXOeLi4iAIQrB+quWyN+QaunKgqoKGle+OFx0t57o59d8ettaS5/9+VBMVYDrjZicNdEtORIijnpL8AOvJh1Ib292Ifu/W6ZjlDhqKOg6hN75191CVpqYmPT09wPUUQOea/0pWVhayXdW6FrYt56gb52MR001vEavK60b8FfbixYu///47SEoEh5pzsH79+hMnToBaA1ZcKQ+WEuUlJOU0zccsjsxmec3/hQUs8Lw4derU2rVrgZKiOoD8999/v3DhArA5XCUhL8yKVY48tm8fzno16iivGo3bl/tckUtJSUG22Rbj87GppgMABv8sjh8/vmHDBuCro9BZOFz7x48f4I09JXH96hhkK1iTkpL8/PyQ5OBU4OJqgPLMy9NhDqnR+kg7CgsL7ezsQH9Ha+4xTy0ly6n77nwoawI6zO3evXvIdopA2XD3uUK4c+fOPXv2gKSkpJz++9mPFgYMw3Bb8iLLeQ9AlA0KCkL47LNLVfnh3QJWFYZhKysroF32W+8F647YkVrR3PDlwhT9EUeANuanUqna2tqs3bKBABC2s1Yb7AodYWxsDHRsbu0/Ew3cw7NqGgqjZxKVfCKBnsc3NTXp6uoiWcnLWRJOLTGCIEiGYOJ3pJPd0rvpxbyk4295W7x48eXLl7tMxt3xRjjAYTKZ+vr6SDZSQCS1QPwa5rfUmzdvBgnitx/jIB31mJmZIdkCH1mp7tTZudqgd+tBfKQde/bs2bFjR5fJuGRBMYoZMmQI+DEdcUGqnJOE3E91sAy2b7s8jUbT1tbucjdYHkUO4cgoNjZ21qxZgFYBaCscPfOwsLBDhw51mYy3RyMZDVlaWgLvrN0MUnQ79usy0JLPi8GDB7MetDc0NIDEeXgqg2gIExkZyT4BadGiRdeuXWN/BRrsu3DhAgRBoBsuvFypqzzrfhsMw3D1KXdJt5OVXeWoq6sTExOztbUFtAduuezNXnwuULp9o4qKCiMjoy1bthQU8HjEk5ycPGnSJASX6/hcc54wGAxtbW0aDezBNyxMwh47dmzYsGGXLl3iuQx+3rx50dHRCC/Z4aHm7bl///60adNAL9x9VTu9wvHjx1lS8JxWFhQUlJCQAHQXVAeQx8bGhoSEgKXtVbr9L7x//97CwmLPnj08W5rt27cfPHgQ+GLMzyfWHL+11Qo42JeTk4NogkCv0W2dqVQqHo/fuHFjfj6Po6+uX78Oet4ZI3ebvYad/8qDdz6B7zzS3Nyso6MDnLyXYTKZBgYG69at4znV8dKlS8uWLQO7EvX5ShP5gRteIZiZ+u8OOAgyCE/7snHjRnd395iYmJYWHmXHy8vr9Wsgx+WYXc5MWWK7MAnIpnPnziE4CU14VIVh2M7ObtGiRRkZGe2/qqmpIZFIQFepPjVad9HTNhiG4coTbkMPgJ574erqmpmZCWqrgBvuTqQoKyuztLQEuUPjxfHifldZd/h+YCg0PBzoOQAMDxs2DOEBbj0LPwrb6tWrx44dGxsbS6VS2387ePBg0JOvOTveyAc4kyZN4u/B8XyGH1InJCTY29sfP36c54Zo7u7uoH7XfoyDZNQzd+7c2NhYsLS9w8mTJ52dnSMjI3l21xEUFd6ygI5iVqxYgWxWfg/Aj3JoaGi4du1anjVbZmamu7s70FW4tEU4MqqqqjI1NQWztwfptrxbt251cXG5ceMGzxONPD09Qc8r4y66yEZD8+fPv3nzJpjJgqWTXiKTydTQ0GB1/1JSUkDjPFzKIBzCFBYWDh48mPX66NGjnNv2ge7Zl5eXB0GQjY0NSOKSN2/KCMbG4hAEQerm5srZ2fldZVFRUdHX1//8+TODwQC4A0wr/VEFVZX/bAOxpxvw50YVFRW7d++2s7MzNTUNDw8vLy9nf5WXl2dtbY3gWjIy0khuLSYmZmRkVFhYCJZcmISl0+mpqamLFi3S1tb28PC4ffs2jUZjf5ubm4tMWAiCIAiH07OxUe06HQRZWVmxnAKA7v/YLq7AlkJHR6dbUoiLi0P0mtx7u9an+p9fSQIzDokUvQh/ynZpaenWrVstLCxsbGxOnTpVW1vL/gpRkYMLz55pm73AHMHufmZmZl+/fkVmbi/AH50rKyv379/v4OBAJBL37NlTWlrK/gpBnSlmtvbJuyvzzb/t9HBa/qC26wwQBEGysrJKSkrV1dVo7O4NqqqqDh06NHToUDwev2PHDs5NWhFoVX3lr2PVg82/zicoag6effUr0N+noKAgKSkJvMWhMLUvDQ0NSUlJc+bM0dDQmDx58qNHjzg7J3l5eZaWliDXERP7fy8Lfh372nnqaAmQXD3bvgDCnxt9//799OnTLi4u2traGzdu/Pz5M/ur/Px8KysroKuo+85zvBU6+8rnitdn37uGLyAC3t3Kyio/v8sOKQuBN9xFRUWcUrAW0LDIy8sDlKKutpZRU1MHQRAE4S0sZIF7e0ik6Hn4U9jIZPKjR49mz549YMAAf3//pKQkthfDMFxcXIzH44EuxNHxRjHAMTc35/xz+xj8kZrJZH758mX16tUEAmHYsGFXrlxpampif/v582dTU1OgC7Uf4yAZ9fRtqSEIguh0elpa2uLFi3V0dNzc3OLi4ji7658+fTIzMwO6UAeyAI5i+p5QfBtuh4eHOzs7t+8Offr0ydzcHOgqXNoiHBkNGDCgoaGBTqcjNV6Q8EFeCoWSkpIyZ84cTU1NPz+/Bw8etLX9dzXURRfpaKjvFN1OeollZWU6Ojqs7h/6godwCIPH41kLIyAIMjU1/fLly39XAvxJZWVlEAQBHiJRV1cHsTfAlJaWriovZ3adi0AgUKlUzsFzhzRF+ppsSINerTQMjGZ9Unp4KI4blXkPQKxFdiNUSEhI0On0pqamL1++bNq0ydTU1NLSMjw8vLq6+sePHwYGBt21s1PweHxJSQlQUmETFofDUSgUCoXy6NGjWbNmaWho+Pn53blzp62tDY2wZbFJatvWOuNA0urr67Ocomu6ryqAXJ1IUVZWpq+vD2QqBEHU67OsvDfH3dk2ffPLZqAchoaGoAUMhL5d5KSkpKhUalNTU05OzqpVqwwNDR0cHCIjIxsbGxEUOfj72ZMtcxZbITrIQ0JCQkZGprkZ7F/pkr6ts4SERFtbW2Nj4/fv37du3WplZcWqMysrK5G4tpiMsp61+/yjKbfGPw396zlAOwRBEKI6s0sEpDMHYmJibW1tDQ0NxcXFO3bssLGxMTExOXToUEVFBbhW9JRHT1RHTl9/4tWn5FVKcfMWR1SA3d3Q0PDHjx9ASYWtfYEgqKGhoampKTY21t/fX1VVdebMmS9evIBhmEKhID7z5E3sG+fAkWA+j6BSFUJVmUxmY2NjRUXFwYMHnZ2d9fX1N2/e/O3bt9LSUmDX1p5x+c6U3Lk2phsklq0ZCrztdY92hwAU45Ri2LBhKKTQc3DQfX5sy92SpsaiJynZLRoaGkA/j0+1XJ8vbND/vTg6OnrSpElsL25qalJQUPgvFg8MigGOnp7ez58/UVjOjSDU5p/UkpKSrN7R69evWcGskSNH3rx5s7W1taWlRaCb07Pp01L/H9z/u+uPHz+ePXs2Z3e9qqpKW1sb/aWBRzHoherbLi8uLs7uDm3bts3a2prdHSovL9fV1UV5XYQjI01NTZQPjPu2vBAEUSiUpqamuLi4qVOnqqqqTpo0KSkpCYbhpqYmZWVlxJdDPhrqa0WXZy+xtraW3Up0o+AhG8KIiYkpKChQKBQIgtTV1evq6thfgZ7r5+fnFxcX9+nTJ5CHMz+OjNC/5luSukYfgqBH85X8GTHk82O6yuXj45OQkDBu3Dh5eeSH/fRhaDTa06dP2x/NIy8vT6PR9PT0ZsyYsXPnTtDLpa7Rn6eQlLMVLEoMQRC0cOHC7OxsPT094BzCQW1t7ZMnT7imguJwOHl5eXFxcTqdnpeXh+SM46qYtXsk1hxsf9RRR2hqao4ePRrcYIESGxvbiRS5ubmgT7AhiN5Q+v7Rkd/nHtY93RgdCHQWk6qq6pgxXbq40EMmk9PT0zkrUBYKCgo0Gk1eXj4iIsLPz6/L63w/v/a24+7lNhJQ/g5r9+pzpYeHghlgZmZmZmYmI4Po/DHhg8lk3r9/v31YU05OrrW1VV1dfcmSJX/++Seia1afcvduPPtmLdAMoClTptTU1ACPmnuZhISElpYWrg9ZWqmqqq5YsWLTpk1dXqTyhKtWpFfxm7UGEMRMWaI9ofly4wUPgLt7eXmx/hRUtvddSktLX79+zfWhmJiYrKysnJwcjUYrKCgYMGAA+AXT1v121efZkRFAwYWmpiZjY+ORI0cisFhIuH//PueUHxYyMjIMBkNPT8/KyurOnTs4XJeD1da8s5ui9WcMfbcl9GjzuscPF1sATZk8e/bsvXv3JCXRHiXLVx4+fNjQ0MD1IUsKXV3dgQMHxsTEAESjah7/OX3e0Vdtg+Z5wWeShjz+ut8J5Inl0aNHL1++TCAQUNkuHBQUFLx9+5brQ7YX43C43NxcoLqLo+ONYoBz+fLl48ePI+mOCh9lZWXZ2dnty7OiomJbW5uYmFhqairQNPP2Yxwko56EhIRjx44J4vR5flFdXf38+XPOKVEQuu46D1kQjGKePn26Zs0aIhF0TrSwcOfOHdaBpJzIysrS6XR9fX13d/ezZ892fRVeRQ7RyMjZ2VlFRYV1EqvIUFZW9vLlS64PcTicnJycjIwMjUb7+vWrlhZA4eOQF8Vo6P79++vXrwedKydIysvLnz9/zvUhq32RkpJSV1d//fr1gAED1v6vvTOPa+Ja//8AssomWKmsCQTCEnYEEQUEhOIGigIVBNSCK9Vba0Vbb+tC1at1aZFqsYoiKraIUXABWUQEXKi4gKwCFiggSwJJgKzfP+b3y+UCyoTMhCSc91/DnMznPJzXLOc85znP2bHD2to6IiJibMUPv+sQDmHs7Ozu3Lkzc+bMmpqatWvXFhUVwecRdY8gCJoyZQoEQcNeTx9ippOT7g8VFWxIfwrU8uYN3TncEcFVcMhreHi4pAyxEEKhUAoLC4edhB2abm5u6urqpqammBrAZrP9/f1nz0boT5AYbt++XVBQMOykmpralClTQkJCrl27Joh/8/9tah6E2NMHExMTI9gF2JCbmzvyJNwUwcHBaWlpyD19EH8D8pxLZ+UQjUt5PJ68vLyYNAWmvH379tmzZ8NOws+yv79/VVUVskQHTWmnj8Z9cXTb///bVabqfN/dKGTT3mvWrJGyDsRIWCzW3bvDZ9umTp0qIyMze/ZsGRkZW1tbQTUHBnh6eki/LGw2e9myZQjXaU4sPB4vKytr2EklJSVZWVlnZ2cul+vg4IBER8vAYGpddTUHMpCDZGfM0Nbq1UJmAJfLXbFiBdZfMdFz9uzZYc4+GRkZVVVVRUXF8PDw1NRUgTx9EFR+/alDyGGkYUQsFktFRUUqX6r37t0bdkZRUXHKlClmZmaOjo6wC2Zsldaz6xOmnnlhb+F/s1Qv0HYfOfpKEBIHHpPJNDQ0XLJkyXhMR5vs7OxhZ/hNYWNjo6OjgyzuTNtr/723+yFe028L5s3bcxGRpw+CICaT6eDgsGLFCkHNliBOnTo1zNkHxz4oKyuHhIRkZWWNY5ZiHAMcFotlbm4eHh4uaF0SRGlp6evXr4edVFdX5/F4AQEBOTk5ovmeslis6dOnR0VFiaCu8XHnzp2RDgI1NTU5ObmQkBBBu+tDEGwUw2KxdHR0pO8Tk5mZOeyMkpKSnJycvb29np6eMB0VgUZGXC43LCxMqCBN8SMlJWWYs09GRkZdXV1eXj48PDwlJQWRp+9/GM9oiMVi4fF4cbh1z507N+wMv5fo5+fX1NQE9xLhhUpC1oVwCMNiseC5TDabPXRSE6mzDw7OHBnSMiqybtHrDeYdT3ztvLLjcFJVyL4QJMMGeAHvokWLxHlOZhx0dHTwj5WUlOTl5Q0NDWNjYz///HN1dfW9e/eOnOX+MJzefzrovYw2BtdcBemwobu7OyIiQvpiBIZ2LDQ0NLhcrre399q1a/39/adMmZKbm0ulUrW0kNx6gxVpidWz/hVqrsgd7Hh065XxCu8x3YQDAwMKCgo+Pj7C/AtoUVFRwT8e2RS3b99mMBgqKiqI9di91TeuvPGL+I8Ckl9TqVQNDQ0xaQpMefnyJX8gOnXqVFlZWWtr65iYmKCgIFVVVV9f376+PgQyRjuf8nbChwJG9vX19S1evFhMwlKwY3BwkD/KVVJSUlBQ0NPTi46OhqeCYmNjka5l5vT10JWnqXM6nvy+NXvuv8lIlxF1d3f7+PhIirOPf6yoqKioqAi3VVhY2IwZM2JiYkYG/Y3KlAXrIlV9d+8Js/qGUJbwh0JUqhMyA6hU6oIFCwgEwrjMF19u3LjBP4bd6wEBAevXr3dzc5ORkbl06RKPx0PkloJ5eb3MPvQg4p/39PR88sknUvlS5TeavLy8ioqKmpraxo0bo6KidHV1yWTyzZs3EanwuNyG7CuFMXHOii0djKnq6ggXAfX09ODxeDFp2I80RVpa2shJ4tHhDfY0v3318PKPB+6Qkm5HIU5e0tPTY2dnJyZNgRFpaWn8Yw0NDTk5ueXLl0dGRrq5uTGZzPT0dGQy/9PxHscAh0qlmpmZSXdT9/f38+9nVVVVLpfr6OgYGxsbEBCgoKAwY8YMZJ7rkWMcwUY9VCp15syZ4tzUQ5Oxjuyuk8lk/nD9owxrFoFHMVQqFY50E/4/Eiv4N+HI7lBCQkJXVxcCjQ/dcgKMjKhUqr+/v5SteBg6O8X34/M7Rampqcg6RUObdzyjISqVamJiIg637kd6iTU1NV988QW/aOTiztEYceMJPoShUqmwJb29vUPjQpA6+4yNjSEIqq2tdXNzQ/Bzy7hrR9cGe+p/P91nZ8bF5ZpIqqirq5s+fbqUefpgGAyGurq6np7ehg0bQkNDZ8yYwS8yNjYeucrgg/wZprEyDYKg+cayD9sS5iK7qLa2VvpGYjBycnLKysru7u4bNmzw8/NTUPjvK9jY2Liurs7Z2XksDU5FwiKP2NwuaDu8abzZ7r+qEUx119XVmZiYjN90tBmzKWxsbMYUYdyLJQWdbuZpmbit2JWYuAhZ+oXa2lqxagpModPpampqBAJh06ZNQUFB06b9NxEy3M729vYYVd3b2wvPFmCkL1YMDAxoaGjA7qpVq1YNTXgBtzMSEd7znxesSOaoKypYRx67+IMd4parr6+XoOUtTCZTQ0NDW1s7JiYmLCxsaIJO5G0FKbkfuXlybWQI4aSsxcofr35ri9Ax1djYON5IBHEHdjT7+fmtX7/e09NTTu6/DqVPP/20tbUVefB4RUaZXchexL4+aX6pysjIaGpqTp06dc2aNZGRkUM7JyYmJvX19YhU9GKSTr2MjrQ+2ilP8Npy/vcFCCc/a2trV65cOR67sQFuiqioqKioqGFNcf78eUQS5GiDNQVEj9Dt1wpWkQRIi1ZbWys+eUiwA36KFy5cGBMT4+HhwXc5KSoqysjIIJoHHd7xFniAU19fj2zoJNkMDAyoqqra29tv3Lhx6dKlQzMyzZgxo62tbexAp5FjHAFHPRLx7YbDS+fNm7dhw4bPPvtsaHcdh8M1NDSMnS/rf5rl5DTBRzES0VDjAO4OaWlpxcTEhIeHD+0O4fH4p0+fji0x4pYTdGTE5XKHpmyTJhQVFRUUFHx9fWNiYry9vYd2inR1dVtaWsZOED8uP8ZQxOrW/VAvEe4iwsd4PL64uHhsrREtI+gQpr+/X1ZWVlFREYKg1tbW/3nfItoMmMeDp1s3bdqE8PeCAu/I5uXlhZH+BEKhUH788ceGhoZRS1+8eDFv3jzsaqfRaLq6ulwuF7sqJor79++npqbSaLRRS+Pi4hITE7Gr/cKFC7GxsdjpC8THm2Lbtm3nzp3DrvbTp09/88032OmLD01NTceOHWtpaRm1NDExcefOndjVnpeXFxAQgJ2++MBkMvfu3VtTUzNqaW5u7rJly7CrvaWlxcLCAjt9dOFyud9//31VVdWopbdv3w4NDcWu9oaGBhsbG+z0J5Dr169nZGQMDAyMWrp69WoymYxd7T/++OOhQ4ew059ADh48+Pz581GLBgcHP/30UzabjV3tJBKpqakJO32BOHTo0F9//TVqEZ1O19PTw7TnZmpq2tbWhp2+OPDnn3+SyeTBwcFRS318fJ4+fSoCMzw8PMrLy0VQ0QTy/PnzM2fOdHV1jVq6cuXK7OxsEZgRGBiYn58vgorGTW5u7qVLl/r6+kYtXbdu3Z9//ikCM1avXn3z5k0RVCRivv/++zdv3oxa9PbtW3t7exHY8ObNGxcXFxFUJGLIZHJ6enp/f/+opZ9//vnt27dFYEZQUBC8JciE8/Feoq6uLtxW5eXlHh4eIrDnyZMnPj4+8HF8fPzRo0f5RUidfT09PXJycqampuhbx+PxeLyEhAQIgvbt24eRvtgCZ2Lu7e3FSD8rKys4OBgjcXHm7t27QUFB2OmHh4dnZGRgp48iN27cCAsLw04/KCjo7t272OlLCm/evHFycsJOf/fu3SdPnsROX1JgMBgzZ85ksVgY6V+4cAG7aS0R09vbq6enh533JCkpadu2bRiJizO///77l19+iZ3+/PnzHz9+jJ2+2OLj41NSUoKReGtrK3adWNSZO3fuh7yiwtPQ0EAikTASlxT27t175MgRrGvp7+/X1dXlcDhYVyTO/Pzzz99++y3WtbDZ7JkzZzIYDKwrwo7k5GTRhBEYGRl1d3eLoCKxAo/Hi+C/xnriXzwRzX8N+0w+FNoiVvj5+RUXF/NE+F46fPjwgQMH4OOgoKCcnBx+EdK8b5qamh4eHrW1tU+ePEF4iUBcvnwZgqDAwEAsxMUZWVnZxYsXX79+HSP9K1euBAUFYSQuzsCDJQqFgoV4f39/Xl6epOw/6+vrW1BQIEhqSAHo7e0tLS2dDKuBxsTc3HxwcLC6uhoj/fT0dCRb/Uo9ysrKrq6uI9P8o4U0vTPV1NTs7Ozy8/Mx0pemthKIgICAjIwM4ZMuj0pra2tjY6OTE8KsiVLF8uXLr169ipH41atXJeh2xbQpJu2TO5SFCxcOTbqEEXfu3Bm6fHhy4u/vL4KmLiwstLW1VVZWxroi7PDz88vMzORyuZjW8uzZM319/aFZaCYJPj4+t27dwroWMpns7++PdS3ihr+/P5lMxrqWoqIiKyuroSkCxBZPT0+47y0nJ+fu7o7dmIUP/8bjcrklJSWurq78IgE+P/D2RidPnkTbPOjJkyfFxcUODg7I9rKUNtauXXvq1CkslNvb2wsLC8Vk4zkRo6CgEBQUlJycjIV4SkrKwoULJeJ1A0GQsrLy4sWLL126hIX4+fPnly9fDucIAKxZs+b06dN07auyAAAVzUlEQVRYKN+9e9fIyMjAAHH2dakGu3dmXV1dfX29NDmv161bh1FbvXnzprW1dTLkohqJtra2i4tLRkYGFuJnzpyJiIiYnN6B0NDQ69evI9vpSDA4HE5SUlJkZCTqyhgRFhZ25coVpJsRCQKbzT537lxERATqypKFk5NTT08PnEQIO5KTk6V7H14kEAgEDQ0NjIJF+Fy4cEHSm/rTTz8lEol5eXmY1iIFDTU+Vq9ejdGokE9zc3NNTc28efMwrUUMweFwM2bMQJScTggk6Nb19/fne5ZFcOPV1NRQKBQHBwcIgh49ejTcJYo8PpDJZBoaGsrJyb18+RLdyENfX18Igq5evYqurAQxf/78rKws1GW/+uqr+Ph41GUlhaamJhwOh3q4b39/P4FA+FBCMfEETvf+oVQL44ZGo+FwOPFJgTTh9PX1GRgYfCip37jhcrmzZ88W80w0ooTL5dra2paWlqKuHBERcebMGdRlJxAOh2NpaYnFesDQ0NDz58+jLispPH361MbGBvUl0p2dnQYGBp2dnejKShBfffUVfykKiqSkpGCa2QMLNm3aNDTtDlqcPXsW08weEkRiYmJ0dDR2+pWVlWZmZpimoZQUrly5smLFCuz03717Z2RkJNFreGGysrJ8fX2x0+/s7DQ0NMQueZSY4+jo+OTJE+z0t23b9p///Ac7fXHmjz/+CAwMxE6/tbVVsp5xGxsbOKE2m802Nzd//fo1dnWtW7fu119/hY+jo6MvXLgwtFQAZx+Px4Mdk3PmzEHx0wUv4HV0dJzMKS3KysrMzc3RvYNfv35tYmIiESvbseObb77ZsWMHupo//PDDli1b0NUUAV9++eWePXvQ1dy5cyfqzSvpnDp1KiQkBF3Nc+fOLV26FF1NSScnJ8fJyQndzH0PHz60trZmMpkoaooDWVlZrq6u6H5h8/Ly7O3tJ/kgNiws7MSJE+hqRkdHHzx4EF1NyaKzs9PIyOhDe5qNDyqVisfjJWuKjsfjtbW1GRoaNjc3o6jZ1dWFevNKLv39/SYmJh/K6C88AQEBKSkpGIlLFmw228rKCjs/S1RU1M8//4yRuCjhcrmzZs3Ky8vDSD82NnYyR4FkZGTwNzFAnYaGBiMjow/tviL1cDgcGxsbOFEdFqxbt+748eMYiWPB8ePH+cmdL1y4gN3WghUVFSYmJvBWIfCE8bCbUDBnH5fLnT9/PgRBu3btQsW+2tpaDQ2NKVOmYOpolwh27Nixfv16tNQYDIatra1odsYRZ2g0moWFBYpfzeLiYlNTUwqFgpagyKBSqWZmZkVFRWgJ5ufnm5ubT9qv2ofgcDjz5s27ePEiWoK1tbVGRkYgfHIka9asQTHnd09Pj5mZGXbdlIll1apVKO5/1dnZSSAQRLORpTjT1taGx+NRXOtw48YNBweHD20eOnm4fPmym5sbiq78kJAQCR3fnjt3bv78+Wh51blcbmBg4E8//YSKmnTwxx9/zJ8/H4toAzKZ7ObmhumWypLF/fv3HR0dsZhOKygosLGxwW7bLhFTWlpqZWWFRQTTs2fPiEQi6qt8JAtfX9/U1FTUZblc7uLFi5OTk1FXliAKCgrs7e2x6MMUFRWRSCTJmoyn0WgGBgbt7e08Ho/D4bi6ut66dQv1WjgcjoeHR3p6Ovznnj17Rkb2CObs4/F4zc3Nn3zyiYyMTFJSkpD2tbe3EwgECIIOHTokpJQUwGQy586dm5CQILwUh8NZuXJlXFyc8FJSQHl5OQ6Hq66uFl6qoaHB2NhYcjdJfPr0qbGxcX19vfBSNTU1OBwOu70CJZq///4bj8c/evRIeKmuri4SiUQmk4WXkj5oNJqtre3ly5eFlxocHFywYAEWK+bEhN7eXisrqz///FN4qYGBAU9PT+mIoRCe+/fvE4nEf/75R3gp+FNVW1srvJQUsH79+nXr1qHiKImPj//ss88kd+1IZGTk5s2bUZHas2dPQEAAcD8NIzg4GPWV4/CqUngNF4DP+vXrt2/fjq5mR0eHiYnJs2fP0JWdWHbs2IH6AnMKhUIkEgsLC9GVlTjg+DvUP7XHjx9fvHgxeLtu3rx569at6GrCE8ySOPo+cuTIunXr4OPKykocDvf333+jW8XevXs///xz+Bj+7nR1dQ37jcDOPh6PV1JSoqKiIisre/LkyXEbV19fb2ZmBkFQVFQUeDZgOjs7SSSSkBmjmExmRETEqlWrJLdrizp37twhEAhCrtSoq6sjEokZGRloWTUhkMlkMzMzIT9yVVVVpqamWGSZlBrKysqMjIyE9Pd1dHTMmjXrl19+Qcsq6aO5udnU1DQtLU0YEQaDERAQEBsbi5ZV4kljY6OJiYmQbzA6nb5w4ULUR2sSTXJysrW1tZBrLf/66y88Hv/w4UO0rJJ0WCzWkiVLtmzZImRP5siRI05OThKdoGpwcNDX1/frr78Wsqu8f/9+V1dXOp2OlmFSA5VKtbKyEvI7MhQKhWJra3vlyhW0BKUGBoPh5OT022+/oSVIp9NdXV1Pnz6NlqCYwGQy3d3dUQzCHRwc9Pb2nrTp5IZBJpMtLS07OjrQEszMzDQzM5vMyXb5DAwMuLi4JCYmoiXIYDDQioUSPUwmk0Qi8bt2qampdnZ2KK4OvHr1KolE4vdwAgMDRw3FG4+zj8fj5eXlqaurQxAUERFBpVIFvfz69etaWlrw5VITd40KHR0dzs7OX3/99fhCVdva2ry9vWNiYiZ5KqWRZGdn43C4ca9rzsnJweFwWMTfip7MzEwcDpeTkzO+y2/fvo3D4e7du4euVdIHPHofd0h/WVkZkUhEsU8srbx7945EIu3fv398ToHGxkYXF5e4uLjJMOfU0NBgYWFx+PDh8f2z9fX1Tk5OqKf+lAJSU1NNTExKSkrGd/m1a9dwOBwqscDSxODgYGhoaEBAQHd39zgu7+/vj4mJ8fT0lMS0G8Po7+9ftmzZypUrx9Hf5vF4dDo9KirK19cXpN34EO/evSMQCNeuXRNeqrOzc/bs2WCt9Ifo6OiwtLQ8e/as8FJUKtXb2/u7774TXkoM6enpsbOzQyUtLJ1OX7JkCT99GIDH4505c8bOzq61tVV4qczMTDweD6Ly+bx//1740CWYvr4+X19fiV6q+PjxYwKBwI+2O3r06OzZs1HxC6elpREIhHfv3sF/JiYm+vr6jjoUGqezj8fjlZeXm5qaQhCkq6t7/vx5hM6p169fBwQEQBAkJye3f//+yTC+EhQGgxETE+Pg4CBQwCqXy7148aKRkZGUbSWJIpWVlTY2NjExMSMDXD9CT0/P5s2braysUN+EegJ5+fKllZXV5s2be3p6kF/V1dUVExNjY2NTWVmJnW3SRHNzs7u7e1BQEP9djAQGg/Hvf//b2Nj4wYMH2NkmTVCp1NDQ0Hnz5gn0kLLZ7MTERENDw0m1EXxPT09QUJCnp2dFRQXyq1gs1s8//2xoaIjKQmCppLS01NTUdOfOnQJtivXPP/+Eh4e7uLiADRNGhcvlHj161NjYWNAbLz8/38rKatzzpmIIh8OJj48nEAg3b94U6MKcnBxzc/Pdu3eDmfWP09TUZGVlFR8fL8zA5OXLlxYWFihGtUgl7e3tjo6OcXFxwtyTtbW1dnZ20r2dUXd3t5ubW2xsrDBJ0JqamlxcXHbv3g1G3MNITk42NTUVZgE4l8s9ceIEkUgEnr5hwCuTduzYIcwzXldX5+DgsH//fhQNmxCOHTvm6+vL740kJCRYWFi8evVq3IJwf4BEIvFHlw8fPjQ2NobzA45k/M4+Ho9HpVI3b94sJycHQZChoeGePXueP38+qk+xo6Pj0qVL/v7+srKyEASRSCQUNwqQSnJzc62trZcsWZKbm/vxiBU6nZ6SkmJraxsQEIBKOjYpZnBw8ODBg/r6+nFxcWO2VWNj43fffWdgYLBv3z54jxtpYmBgYN++fQYGBt9++21jY+PHf/z27dtdu3bp6+sfPHgQJI8XCA6Hc+bMGUNDw02bNo35Zm9razty5IiRkdG2bdvGF8ExmYGXqIeEhBQVFX28U0ulUpOSkszNzVetWtXS0iIyC8WH9PR0AoEQFhY2ZjAahUI5ffq0mZlZREQEKpnppBgajbZz505DQ8MDBw6MeV9VVVX961//MjAwOHnyJIjE/ziVlZULFixwdna+du3ax1PLs9nsO3fu+Pr6Ojk5SeVmOy9fvvT09JwzZ056evrHv8UsFiszM9PLy2v27NlgLx2EUCiUFStW+Pj4jGPjZiaTeeTIETwef//+fSxskzLodPqaNWvmzJkzjnl0Npt96tQpHA43GdIZDwwMbN682dHRcRw+KS6X+/vvv+NwOLCi/EM8fPjQxMTkwIED4xjlvX37dtGiRUuWLBEogmTywGAwvvjii9mzZ5eXlwt6LYfD+fXXX42MjK5fv46FbaJn48aNq1at4nf2srOz8Xj80aNHxzEfWV1d7eXlFRwczB8nvnz58uNp9GV4PB4kHK9evYqPj09PT2ez2RAETZs2zcLCQl9fX1NTc2BgoKurq6ampq6uDq4InvqOiIiQl5cXsl6ph8vlksnkpKSkioqK+fPnu7q6mpmZaWtrq6qqUiiUtra2ysrKhw8fPnv2zNfXNzY21sHBYaJNlgzg4evFixenTp3q5eVlZ2eHx+M1NTVlZGQoFEpDQ8OLFy/y8vIoFMrq1as3btyora090SZjRXd396+//pqSkqKhoeHl5WVraws3BQRBPT09DQ0N5eXleXl5NBotMjJyw4YNcBFAUBgMxvnz58+fPz8wMODj42Nvb29qaqqpqSkvL0+hUJqamioqKvLz8xsbG0NDQ7ds2aKvrz/RJkskbDb72rVrZ8+ebWpq8vb2njVrFpFI1NLSUlFR6enpaWlpqaysLCwsfPHixdKlS2NjYy0tLSfa5AmDzWZfuXLl7NmzLS0tPj4+Tk5ORCJx2rRpysrKFAqlubm5srLywYMHr1+/DgwMjI2NJRKJE22yZNDe3p6QkHD58uWZM2d6eXmRSCQcDqepqcnhcCgUSn19/fPnz3Nzc7lc7po1a6Kjo1VVVSfaZMmguLg4MTExPz/f1dXVzc3NwsJCR0dHQ0Ojr6+vs7Ozurq6uLj4wYMHDg4OGzZs8Pf3n2h7MaSwsDAxMfHhw4dubm5ubm7m5uY6Ojrq6up9fX3v37+vqqoqLi4uLCx0dnbeuHHjggULJtpeCSMtLW3Xrl3+/v5ff/01Ho8f8/dMJvPq1asHDx6cM2fO4cOHp0+fLgIjpYNbt25t37597ty5O3fuRPKJ4XA46enp8fHxlpaWP/30k66urgiMFAdycnK2bt3q4OAQFxdHIpHG/D2Xy71169b+/fsNDAxOnDhhZGQkAiMlFAqF8u233+bk5HzzzTerV69WVFQc85Lm5uajR4/euHFj7969ERERMjIyIrBTQsnKyvrqq6/mzJmzc+dOc3PzMX/P4XAyMjIOHDhAJBKPHTump6cnAiNFAIfDWbt2LY1GS01NVVJSgiDo/fv3cXFxJSUlcXFxoaGhCgoKY4o0NDQcOXLk7t27hw4dCg4Ohk8+evQoPDz8woUL7u7uH7oQBWcfTHt7+7Vr17Kzs4uLi7u7u4cWycvLW1paenl5LV++3M3NDTwVgtLd3V1QUFBWVlZTU0OhUGg0mqqqqq6uLpFIdHV1dXV1he8bgKBUV1c/ePDg1atXb9++pdFoEASpqqricDgSieTh4TGpHAFv3ryBm6KxsZHfFHg83tra2sPDA8kLGoCEpqam/Pz8V69e1dTU0Gg0FoulpqZmaGhoaWnp7u5ua2sLxz4DhKSjoyM/P//58+e1tbVUKpXBYKiqqurr6xOJRDc3NxcXFzDbxKe9vT0vL6+8vBxuq/7+fjU1Nbit5s6d6+zsPGXKlIm2UfLg8XivXr168OBBZWVlU1MTjUaTlZVVVVUlEAjW1taenp4mJiYTbaNEwmAwHj169Pjx46qqqvb2drg7pKWlZWZm5uzs7O7urqGhMdE2igg6nV5UVPT48ePq6uqOjg64KbS1teGm8PDwUFNTm2gbJZX+/v6kpKSEhAR9ff1ly5Z5eXkRicRhb8LOzs5Hjx7du3fv5s2bnp6eCL0wgGEwmczk5OQTJ05oaWktX77cy8vLyspq2Ae6p6enuLg4Ozv7xo0bzs7OcXFxjo6OE2XwRMFisS5fvnzs2DFlZeWgoCAvLy9ra+thDgIqlVpaWpqdnX39+nUbG5u4uDhXV9eJMliyqK6uPnz48L179xYtWuTn5zd37lwdHZ2hP+BwODU1Nfn5+WQyuba2dtOmTRs2bABzdUhgMpkXL148fvy4pqYm/IyTSKSRz3hJSUl2dnZGRoaTk9OuXbucnJwmymCM4PF4O3bsKCoqSktL4/vfX79+fejQoYKCgqVLl/r5+bm5uQ2bLmKz2dXV1bm5uTdu3Ghubt6yZUt0dLSysjJc+ttvvx0+fDg9Pd3Ozu4jVaPm7BvK+/fvW1tbaTSakpKSpqamoaEhGFkBAAAAAAAAAAAA8YfH45WUlNy+fbugoKC2tlZdXV1DQ0NJSam3t7ejo0NDQ8PFxcXHx2fp0qVg0YPwPHv2LDMzs6CgoLq6WllZWUtLS0lJiUajtbe3q6qqzpo1y9vbOzAwUIqX2iDkxYsXt27dKigoqKioUFJS0tbWVlJSotPp7e3tysrKjo6OcEMN81UBkNDX13fr1q2cnJzS0tLu7m4dHR01NTUmk9nb20uhUExMTDw8PPz9/efOnQsm5sdBWVkZ/IxXVVUNfcY7OjpUVFT4z7h0B0dfvXo1Li7uhx9+iIyM5Ee/USiUmzdv3r9///Hjx1QqdcaMGerq6gMDA1Qqta+vz9TU1MPDY+HCha6urvxLWlpatm7dSqVSL126NObDjomzDwAAAAAAAAAAAACkgN7eXjg5l4qKira2NghiwI6+vr6uri4ul6usrKytrY1kgdvkhEajdXV1cTgc2OuHZAkqACFsNrurq4tOp0MQpK2tPXlixkVDX19fd3c3h8OZhM94Y2Pj1q1bOzs79+3b5+3tPayUxWJ1d3fT6XQZGRltbW11dfVhP6BSqb/88ktSUlJcXNz69euR+J2Bsw8AAAAAAAAAAAAAAAAAAADAkPv37x84cIDBYKxduzY4OFhLS2vMS548eXLp0qWMjIzw8PDt27cjD4EEzj4AAAAAAAAAAAAAAAAAAAAAc8rKyi5evJiRkaGrqwvvFmBqajp9+nRVVdXBwcG+vr6WlpaqqqqysrIHDx7g8fjw8PDg4GBBE0cAZx8AAAAAAAAAAAAAAAAAAAAgOiorK4uKiqqqqmpqauC98iAImjZtmo6Ojrm5ua2trYeHx8glvQgBzj4AAAAAAAAAAAAAAAAAAABASgC7yQAAAAAAAAAAAAAAAAAAAABSAnD2AQAAAAAAAAAAAAAAAAAAAFLC/wEoWeqTA4oG8gAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def generate_phone_wfst(f, start_state, phone, n):\n",
    "    \"\"\"\n",
    "    Generate a WFST representing an n-state left-to-right phone HMM.\n",
    "    \n",
    "    Args:\n",
    "        f (fst.Fst()): an FST object, assumed to exist already\n",
    "        start_state (int): the index of the first state, assumed to exist already\n",
    "        phone (str): the phone label \n",
    "        n (int): number of states of the HMM excluding start and end\n",
    "        \n",
    "    Returns:\n",
    "        the final state of the FST\n",
    "    \"\"\"\n",
    "    \n",
    "    current_state = start_state\n",
    "    \n",
    "    for i in range(1, n+1):\n",
    "    \n",
    "        in_label = state_table.find('{}_{}'.format(phone, i))\n",
    "        \n",
    "        # self-loop back to current state\n",
    "        f.add_arc(current_state, fst.Arc(in_label, 0, None, current_state))\n",
    "        \n",
    "        # transition to next state\n",
    "        \n",
    "        # we want to output the phone label on the final state\n",
    "        # note: if outputting words instead this code should be modified\n",
    "        if i == n:\n",
    "            out_label = phone_table.find(phone)\n",
    "        else:\n",
    "            out_label = 0   # output empty <eps> label\n",
    "            \n",
    "        next_state = f.add_state()\n",
    "        f.add_arc(current_state, fst.Arc(in_label, out_label, None, next_state))    \n",
    "       \n",
    "        current_state = next_state\n",
    "    return current_state\n",
    "\n",
    "f = fst.Fst()\n",
    "start = f.add_state()\n",
    "f.set_start(start)\n",
    "\n",
    "last_state = generate_phone_wfst(f, start, 'p', 3)\n",
    "\n",
    "f.set_input_symbols(state_table)\n",
    "f.set_output_symbols(phone_table)\n",
    "\n",
    "def generate_word_wfst(f, start_state, word, n):\n",
    "    \"\"\" Generate a WFST for any word in the lexicon, composed of n-state phone WFSTs.\n",
    "        This will currently output phone labels.  \n",
    "    \n",
    "    Args:\n",
    "        f (fst.Fst()): an FST object, assumed to exist already\n",
    "        start_state (int): the index of the first state, assumed to exist already\n",
    "        word (str): the word to generate\n",
    "        n (int): states per phone HMM\n",
    "        \n",
    "    Returns:\n",
    "        the constructed WFST\n",
    "    \n",
    "    \"\"\"\n",
    "\n",
    "    current_state = start_state\n",
    "    \n",
    "    # iterate over all the phones in the word\n",
    "    for phone in lex[word]:   # will raise an exception if word is not in the lexicon\n",
    "        # your code here\n",
    "        \n",
    "        current_state = generate_phone_wfst(f, current_state, phone, n)\n",
    "    \n",
    "        # note: new current_state is now set to the final state of the previous phone WFST\n",
    "        \n",
    "    f.set_final(current_state)\n",
    "    \n",
    "    return f\n",
    "\n",
    "f = fst.Fst()\n",
    "start = f.add_state()\n",
    "f.set_start(start)\n",
    "\n",
    "generate_word_wfst(f, start, 'peppers', 3)\n",
    "f.set_input_symbols(state_table)\n",
    "f.set_output_symbols(phone_table)\n",
    "\n",
    "# We need special code to display the higher-resolution WFSTs inside Jupyter notebook\n",
    "from subprocess import check_call\n",
    "from IPython.display import Image\n",
    "f.draw('tmp.dot', portrait=True)\n",
    "check_call(['dot','-Tpng','-Gdpi=200','tmp.dot','-o','tmp.png'])\n",
    "Image(filename='tmp.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that the WFST has been constructed, we can traverse over the states and arcs.  This example shows how you can do this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0 19 0 0 0\n",
      "0 19 0 0 1\n",
      "1 20 0 0 1\n",
      "1 20 0 0 2\n",
      "2 21 0 0 2\n",
      "2 21 7 0 3\n",
      "3 22 0 0 3\n",
      "3 22 0 0 4\n",
      "4 23 0 0 4\n",
      "4 23 0 0 5\n",
      "5 24 0 0 5\n",
      "5 24 8 0 6\n",
      "6 19 0 0 6\n",
      "6 19 0 0 7\n",
      "7 20 0 0 7\n",
      "7 20 0 0 8\n",
      "8 21 0 0 8\n",
      "8 21 7 0 9\n",
      "9 25 0 0 9\n",
      "9 25 0 0 10\n",
      "10 26 0 0 10\n",
      "10 26 0 0 11\n",
      "11 27 0 0 11\n",
      "11 27 9 0 12\n",
      "12 28 0 0 12\n",
      "12 28 0 0 13\n",
      "13 29 0 0 13\n",
      "13 29 0 0 14\n",
      "14 30 0 0 14\n",
      "14 30 10 0 15\n"
     ]
    }
   ],
   "source": [
    "for state in f.states():\n",
    "    \n",
    "    # iterate over all arcs leaving this state    \n",
    "    for arc in f.arcs(state):\n",
    "         print(state, arc.ilabel, arc.olabel, arc.weight, arc.nextstate)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alternatively, we could begin at the start state, and traverse in a depth-first manner.  **Warning**: the code below specifically handles self-loops, but won't work if your WFST has larger cycles in it!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0 19 0 0 0\n",
      "0 19 0 0 1\n",
      "1 20 0 0 1\n",
      "1 20 0 0 2\n",
      "2 21 0 0 2\n",
      "2 21 7 0 3\n",
      "3 22 0 0 3\n",
      "3 22 0 0 4\n",
      "4 23 0 0 4\n",
      "4 23 0 0 5\n",
      "5 24 0 0 5\n",
      "5 24 8 0 6\n",
      "6 19 0 0 6\n",
      "6 19 0 0 7\n",
      "7 20 0 0 7\n",
      "7 20 0 0 8\n",
      "8 21 0 0 8\n",
      "8 21 7 0 9\n",
      "9 25 0 0 9\n",
      "9 25 0 0 10\n",
      "10 26 0 0 10\n",
      "10 26 0 0 11\n",
      "11 27 0 0 11\n",
      "11 27 9 0 12\n",
      "12 28 0 0 12\n",
      "12 28 0 0 13\n",
      "13 29 0 0 13\n",
      "13 29 0 0 14\n",
      "14 30 0 0 14\n",
      "14 30 10 0 15\n"
     ]
    }
   ],
   "source": [
    "def traverse_arcs(state):\n",
    "    \"\"\"Traverse every arc leaving a particular state\n",
    "    \"\"\"\n",
    "    for arc in f.arcs(state):\n",
    "        print(state, arc.ilabel, arc.olabel, arc.weight, arc.nextstate)\n",
    "        \n",
    "        if arc.nextstate != state:   # don't follow the self-loops or we'll get stuck forever!\n",
    "            traverse_arcs(arc.nextstate)\n",
    "\n",
    "s = f.start()\n",
    "traverse_arcs(s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a more readable table, you could find the indexes of the input and output labels in your symbol tables and print the string instead."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises\n",
    "\n",
    "1. Write code to randomly generate (sample) a path through your word HMM for \"*peppers*\".  You should output the sequence of input and output labels along the path.  To sample from a list of arcs, you can use code like\n",
    "\n",
    "```python\n",
    "import random\n",
    "\n",
    "arc_list = list(f.arcs(state))\n",
    "sampled_arc = random.choice(arc_list)\n",
    "```\n",
    "\n",
    "  Notice that if you repeat your random sampling by running the code multiple times, you'll get paths of different lengths due to the self-loops"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import random\n",
    "\n",
    "def sample_random_path(f):\n",
    "    '''Given an FST, randomly sample a path through it.\n",
    "    \n",
    "        Args:\n",
    "            f (fst.Fst()): an FST\n",
    "        \n",
    "        Returns:\n",
    "            input_label_seq (list(str)): the list of input labels from the arcs that were sampled\n",
    "            output_label_seq (list(str)): the list of output labels from the arcs that were sampled\n",
    "        '''\n",
    "    curr_state = f.start() # start from beginning\n",
    "    weight_type = f.weight_type() # type of weights used in the fst\n",
    "    input_label_seq = []\n",
    "    output_label_seq = []\n",
    "\n",
    "    while f.final(curr_state) == fst.Weight(weight_type, 'inf'): # the .final method returns the probability of a state being final\n",
    "                                                             # it's infinite when the state is NOT final\n",
    "            \n",
    "            # your code here\n",
    "            \n",
    "    return input_label_seq, output_label_seq\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**You may check your answer with the code below, where you can see the expected results and change your code accordingly.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "input_label_seq, output_label_seq = sample_random_path(f)\n",
    "print('\\n'.join(['{} {}'.format(input_label_seq[i], output_label_seq[i]) for i in range(len(input_label_seq))]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get a similar result as below after running the above cell,\n",
    "```\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_3 <eps>\n",
    "p_3 p\n",
    "eh_1 <eps>\n",
    "eh_1 <eps>\n",
    "eh_2 <eps>\n",
    "eh_2 <eps>\n",
    "eh_3 eh\n",
    "p_1 <eps>\n",
    "p_1 <eps>\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_3 p\n",
    "er_1 <eps>\n",
    "er_2 <eps>\n",
    "er_2 <eps>\n",
    "er_3 er\n",
    "z_1 <eps>\n",
    "z_1 <eps>\n",
    "z_2 <eps>\n",
    "z_3 <eps>\n",
    "z_3 <eps>\n",
    "z_3 z\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2. Now it's time to add probabilities to your WFST.  As mentioned at the end of Lab 1, probabilities in WFSTs are traditionally expressed in negative log format, that is, the weight $w$ on an arc transitioning between states $i$ and $j$ is given by $w=-\\log a_{ij}$, where $a_{ij}$ is the HMM transition probability.  Remember that you can add weights using the third argument to `fst.Arc()`.\n",
    "\n",
    "  You should now modify your code above to add weights to your word and phone recognition WFSTs from Lab 1, corresponding to transition probabilities.  Assume that the probability of a self-loop is $0.1$, and that when transitioning *between* separate multiple sets of phones (or words), the probabilities are uniform over all transitions.\n",
    "\n",
    "  Remember to set your fst to use log probabilities and use log weights:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```python\n",
    "import math\n",
    "f = fst.Fst('log')\n",
    "\n",
    "s1 = f.add_state()\n",
    "s2 = f.add_state()\n",
    "weight = fst.Weight('log', -math.log(0.1))\n",
    "f.add_arc(s1, fst.Arc(0, 0, weight, s2))\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "\n",
    "def generate_phone_wfst(f, start_state, phone, n):\n",
    "    \"\"\"\n",
    "    Generate a WFST representing an n-state left-to-right phone HMM\n",
    "    \n",
    "    Args:\n",
    "        f (fst.Fst()): an FST object, assumed to exist already\n",
    "        start_state (int): the index of the first state, assmed to exist already\n",
    "        phone (str): the phone label \n",
    "        n (int): number of states for each phone HMM\n",
    "        \n",
    "    Returns:\n",
    "        the final state of the FST\n",
    "    \"\"\"\n",
    "    \n",
    "    current_state = start_state\n",
    "    \n",
    "    for i in range(1, n+1):\n",
    "        \n",
    "        in_label = state_table.find('{}_{}'.format(phone, i))\n",
    "        \n",
    "        # Your code here\n",
    "        \n",
    "    return current_state"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_word_wfst(word):\n",
    "    \"\"\" Generate a WFST for any word in the lexicon, composed of 3-state phone WFSTs.\n",
    "        This will currently output word labels.  \n",
    "        Exercise: could you modify this function and the one above to output a single phone label instead?\n",
    "    \n",
    "    Args:\n",
    "        word (str): the word to generate\n",
    "        \n",
    "    Returns:\n",
    "        the constructed WFST\n",
    "    \n",
    "    \"\"\"\n",
    "    f = fst.Fst('log')\n",
    "    \n",
    "    # Your code here\n",
    "    \n",
    "    return f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Here is the code for debugging\n",
    "f = generate_word_wfst('peppers')\n",
    "f.set_input_symbols(state_table)\n",
    "f.set_output_symbols(phone_table)\n",
    "\n",
    "\n",
    "# We need special code to display the higher-resolution WFSTs inside Jupyter notebook\n",
    "from subprocess import check_call\n",
    "from IPython.display import Image\n",
    "f.draw('tmp.dot', portrait=True)\n",
    "check_call(['dot','-Tpng','-Gdpi=200','tmp.dot','-o','tmp.png'])\n",
    "Image(filename='tmp.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get an equivalent result as below\n",
    "\n",
    "![generate_word_wfst](figs/generate_word_wfst_weight.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "3. Modify your answer to exercise 1 to sample a path through the word HMM *and* also compute the negative log probability of the path.  This gives you $-\\log p(Q)$ in the lecture notation.  (Recall that $\\log ab = \\log a + \\log b$)\n",
    "\n",
    "  **Note**: Internally OpenFst stores weights in a special object that you will need to convert to a float, using the `float()` function, before adding your negative log probabilities.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sample_random_path_prob(f):\n",
    "    '''Given an FST, randomly sample a path through it and compute the negative log probability.\n",
    "    \n",
    "        Args:\n",
    "            f (fst.Fst()): an FST\n",
    "        \n",
    "        Returns:\n",
    "            input_label_seq (list(str)): the list of input labels from the arcs that were sampled\n",
    "            output_label_seq (list(str)): the list of output labels from the arcs that were sampled\n",
    "            neg_log_prob (float): negative log probability of the sampled path\n",
    "        '''\n",
    "    curr_state = f.start() # start from beginning\n",
    "    weight_type = f.weight_type() # type of weights used in the fst\n",
    "    input_label_seq = []\n",
    "    output_label_seq = []\n",
    "    neg_log_prob = 0.0 # log(1) = 0\n",
    "\n",
    "    while f.final(curr_state) == fst.Weight(weight_type, 'inf'): # the .final method returns the probability of a state being final\n",
    "                                                             # it's infinite when the state is NOT final\n",
    "            \n",
    "            # your code here\n",
    "            \n",
    "    return input_label_seq, output_label_seq, neg_log_prob\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Here is the code to help you debug your solutions\n",
    "input_label_seq, output_label_seq, neg_log_prob = sample_random_path_prob(f)\n",
    "\n",
    "print('\\n'.join(['{} {}'.format(input_label_seq[i], output_label_seq[i]) for i in range(len(input_label_seq))]))\n",
    "print(neg_log_prob)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get a similar result as below after running the above cell (although remember that the paths are randomly sampled)\n",
    "```\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_3 <eps>\n",
    "p_3 p\n",
    "eh_1 <eps>\n",
    "eh_2 <eps>\n",
    "eh_3 eh\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_3 p\n",
    "er_1 <eps>\n",
    "er_1 <eps>\n",
    "er_2 <eps>\n",
    "er_3 <eps>\n",
    "er_3 <eps>\n",
    "er_3 er\n",
    "z_1 <eps>\n",
    "z_1 <eps>\n",
    "z_2 <eps>\n",
    "z_3 <eps>\n",
    "z_3 z\n",
    "22.303673804999995\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "4. You are now given a set of observations, ($x_1, \\dotsc, x_t, \\dotsc$).  Can you use your WFST for the word \"*peppers*\" to compute $p(X,Q)$ for a randomly sampled path $Q$ through the HMM?  For now, we won't use real samples $x_t$, and will instead assume that you already have a function `observation_probability(state, t)` that computes $b_j(t) = p(x_t|q_t=j)$, provided here:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def observation_probability(hmm_label, t):\n",
    "    \"\"\" Computes b_j(t) where j is the current state\n",
    "    \n",
    "    This is just a dummy version!  In later labs we'll generate \n",
    "    probabilities for real speech frames.\n",
    "    \n",
    "    You don't need to look at this function in detail.\n",
    "    \n",
    "    Args: hmm_label (str): the HMM state label, j.  We'll use string form: \"p_1\", \"p_2\", \"eh_1\" etc  \n",
    "          t (int) : current time step, starting at 1\n",
    "          \n",
    "    Returns: \n",
    "          p (float): the observation probability p(x_t | q_t = hmm_label)\n",
    "    \"\"\"\n",
    "    \n",
    "    p = {} # dictionary of probabilities\n",
    "    \n",
    "    assert(t>0)\n",
    "    \n",
    "    # this is just a simulation!\n",
    "    if t < 4:\n",
    "        p = {'p_1': 1.0, 'p_2':1.0, 'p_3': 1.0, 'eh_1':0.2}\n",
    "    elif t < 9:\n",
    "        p = {'p_3': 0.5, 'eh_1':1.0, 'eh_2': 1.0, 'eh_3': 1.0}\n",
    "    elif t < 13:\n",
    "        p = {'eh_3': 1.0, 'p_1': 1.0, 'p_2': 1.0, 'p_3':1.0, 'er_1':0.5}\n",
    "    elif t < 18:\n",
    "        p = {'p_3': 1.0, 'er_1': 1.0, 'er_2': 1.0, 'er_3':0.7}\n",
    "    elif t < 25:\n",
    "        p = {'er_3': 1.0, 'z_1': 1.0, 'z_2': 1.0, 'z_3':1.0}\n",
    "    else:\n",
    "        p = {'z_2': 0.5, 'z_3': 1.0}\n",
    "        \n",
    "    for label in ['p_1', 'p_2', 'p_3', 'eh_1', 'eh_2', 'eh_3', 'er_1', 'er_2', 'er_3', 'z_1', 'z_2', 'z_3']:        \n",
    "        if label not in p:\n",
    "            p[label] = 0.01  # give all other states a small probability to avoid zero probability\n",
    "            \n",
    "    # normalise the probabilities:\n",
    "    scale = sum(p.values())\n",
    "    for k in p:\n",
    "        p[k] = p[k]/scale\n",
    "        \n",
    "    return p[hmm_label]\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Enter your code below.  You might want to convert the observation probabilities into negative log probabilities.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sample_random_path_obs_prob(f):\n",
    "    '''Given an FST and observation probabilities, randomly sample a path\n",
    "        through it and compute the negative log probability.\n",
    "    \n",
    "        Args:\n",
    "            f (fst.Fst()): an FST\n",
    "        \n",
    "        Returns:\n",
    "            input_label_seq (list(str)): the list of input labels from the arcs that were sampled\n",
    "            output_label_seq (list(str)): the list of output labels from the arcs that were sampled\n",
    "            neg_log_prob (float): negative log probability of the sampled path\n",
    "        '''\n",
    "    t = 1\n",
    "    curr_state = f.start() # start from beginning\n",
    "    weight_type = f.weight_type() # type of weights used in the fst\n",
    "    input_label_seq = []\n",
    "    output_label_seq = []\n",
    "    neg_log_prob = 0.0 # log(1) = 0\n",
    "\n",
    "    while f.final(curr_state) == fst.Weight(weight_type, 'inf'):\n",
    "        \n",
    "    # your code here\n",
    "    \n",
    "    \n",
    "    return input_label_seq, output_label_seq, neg_log_prob\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Here is the code to help you debug your solutions\n",
    "input_label_seq, output_label_seq, neg_log_prob = sample_random_path_obs_prob(f)\n",
    "print('\\n'.join(['{} {}'.format(input_label_seq[i], output_label_seq[i]) for i in range(len(input_label_seq))]))\n",
    "print(neg_log_prob)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get a similar result as below after running the above cell,\n",
    "```\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_3 <eps>\n",
    "p_3 <eps>\n",
    "p_3 <eps>\n",
    "p_3 <eps>\n",
    "p_3 p\n",
    "eh_1 <eps>\n",
    "eh_1 <eps>\n",
    "eh_1 <eps>\n",
    "eh_2 <eps>\n",
    "eh_2 <eps>\n",
    "eh_2 <eps>\n",
    "eh_3 eh\n",
    "p_1 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_2 <eps>\n",
    "p_3 <eps>\n",
    "p_3 p\n",
    "er_1 <eps>\n",
    "er_1 <eps>\n",
    "er_1 <eps>\n",
    "er_2 <eps>\n",
    "er_3 er\n",
    "z_1 <eps>\n",
    "z_2 <eps>\n",
    "z_2 <eps>\n",
    "z_3 z\n",
    "170.48377790878732\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You might have noticed that the dummy observation probability function above effectively allows the observation sequence $x_t$ to be arbitrarily long.  This is simply to allow it to match the length of your sampled path $Q$.  In real use, the observation sequence will have a fixed length $T$, and any matching path through the HMM will have to have the same length.  We'll explore this more when writing the Viterbi decoder in the next lab.\n",
    "\n",
    "## If you have more time\n",
    "\n",
    "You might like to start thinking about how to implement the Viterbi algorithm over HMMs in WFST form.  Try working with the \"*peppers*\" example above.  You'll need to write functions to compute and store the probabilities $V_j(t)$, giving the probability up to time step $t$ of the observation sequence $(x_1, \\dotsc, x_t)$ along the most likely path $(q_1, \\dotsc, q_t)$."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}