make_wfirst_Lyotstop.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Make a Lyot stop based on the WFIRST Cycle 5 telescope pupil"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "import skimage.transform\n",
    "import astropy.io.fits as fits\n",
    "import matplotlib.pyplot as plt\n",
    "import os\n",
    "import glob\n",
    "import PIL.ImageDraw\n",
    "import PIL.Image\n",
    "%pylab inline --no-import-all\n",
    "matplotlib.rcParams['image.origin'] = 'lower'\n",
    "matplotlib.rcParams['image.interpolation'] = 'nearest' \n",
    "matplotlib.rcParams['image.cmap'] = 'gray'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Set basic parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 458,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "overwrite = False\n",
    "secobs_key = \"wfirstCycle5\"\n",
    "flip = True\n",
    "iD = 40\n",
    "oD = 90\n",
    "D = 1000 # original pupil diameter\n",
    "N = 128 # quadrant width after binning\n",
    "symm = 'half' # set to either 'quart' or 'half'\n",
    "pad = 19 # feature padding, units of thousandths of pupil diameter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 459,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Destination /astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS already exists\n"
     ]
    }
   ],
   "source": [
    "#LS_dir = os.path.abspath(\"../InputMasks/LS\")\n",
    "#LS_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks/LS\")\n",
    "#LS_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks_v2/LS\")\n",
    "#LS_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks_v3/LS\")\n",
    "#LS_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks_v4/LS\")\n",
    "#LS_dir = os.path.normpath(\"/Users/neil/Box Sync/scda/InputMasks/LS\")\n",
    "#LS_dir = os.path.expanduser(\"~/Box Sync/scda/InputMasks_wfirst/LS\")\n",
    "LS_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS\")\n",
    "telap_dir = os.path.normpath(\"/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/TelAp\")\n",
    "\n",
    "if not os.path.exists(LS_dir):\n",
    "    os.mkdir(LS_dir)\n",
    "    print(\"created {:s} for LS arrays\".format(LS_dir))\n",
    "else:\n",
    "    print(\"Destination {:s} already exists\".format(LS_dir))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 460,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x19f458690>"
      ]
     },
     "execution_count": 460,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAFkCAYAAAAg+UxRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAHVxJREFUeJzt3X/MZNV93/H3lxBwaLSiqsOuUERNSrJO5AjKYhzkhtjF\nqmu7xnEsxayjOLLluk5AslAlYxqstqC2CKuAScBCjdUoJF0L4zjGVgTBuLEwYJAXl6g2QaWGULRm\nYxxaVmB7DZz+ce/A7Ozz7M7Mc3+cc+77JT3SPjN3nznn3HPOZ879MRMpJSRJOmbsAkiS8mAgSJIA\nA0GS1DIQJEmAgSBJahkIkiTAQJAktQwESRJgIEiSWgaCJAkwECRJLQNBkgQYCJKkloEgSQLg2LEL\n0JWI+AfAm4HHgB+MWxpJ6sQrgFcBt6eUvtf3i1UTCDRh8CdjF0KSevAbwH/r+0VqCoTHxi6AxrV3\n796xi9CbXbt2jV0EjeuxIV6kpkDwMNGEpJQ4//zzufXWW8cuyiDmv9lwVu+IGLFEGtgg81tNgaCK\n+VWvh1tsEwNCW2UgKEsGwOoMCG2VgaBsGALdmm9Pw0HLMBA0mq0GwO7duzsqSVnWqberBy0janlX\nFhFnAvVeZlKJWvpbTQyHIuxKKT3Q94u4QlDvDIG8eWhJMwaCemMQlGe2zwyGafKzjNQbJ5Uyud+m\ny0CQJAEGgiSpZSBoLSklzxFMlPu+XgaCVrLqZODx6LKssr8Mhvp4lZGOykGvI/Gy1XoYCNqUQaBV\nedlq2TxkpMOscijA0JiWVfqFfaM8BoIO0ccg9t1iGfrYT4ZCWTxkJMCBq/54GKkcBsLEGQQaisGQ\nPw8ZTVSXx3gNlWnosr/YZ/JkIEyQg1E5sB/mx0CYkDHfmXmYIG9j7R9XC3kxECbCQaec2T/z4Enl\nyjnQVApPOo/PQKiUQaBSGQzj8ZBRhcYIAwOobvapaTAQKuMgUk3sz8MyECpRwtUaHgLIU+77pYS+\nXQsDoQIOFk2B/bx/BkLhchokOZVF3clpv+ZUlhp5lVGhHBiaKq9C6o8rhAKVHAYO4ryUvD9KHge5\nMhAK4yCQXuZ46JaBUJASOn8JZdTyStifJZSxFAZCIez00uYcH90wEDJX4zXYJR+3rklt+6HGsTI0\nAyFjdm5pdY6b9RkImSq5U5dcdr2s5P1YctnHZCBkyM4sbZ3jaHUGQmam0olrO35dmqm0/1TGU1cM\nhIzYeaXuOa6WZyBkorZOW1t9pqa2/VdbffpiIGTAzir1z3F2dAbCyKbcSadyHDs3U273KY+3ZRgI\nI7JzSsNz3G3OQBiJnVIaj+NvYwbCCOyM0vgch4czEAY2pU44pbrWZEr7bUp1XYaBMCA7n5Qfx+XL\nDISB2Ok2NuUrXsZge2/M8dkwEAZgZ5Py5zg1ENQzB1lZ3F/TZiD0zAEmlWPq49VA6NHUO9eyPK49\nDNt5OVMetwZCT6bcqaTSTXX8dh4IEXFMRFwREd+OiOci4pGIuGyD7S6PiH3tNndExGkLzx8fEddH\nxFMRcSAibomIk7oubx+m2pmkmkxxHPexQvgo8K+A3wFeDXwE+EhEXDTbICIuAS4CPgicDTwL3B4R\nx839nWuBtwHvAs4FTgY+20N5OzXFTnQ0tkkZ3E+Hm1qbHNvD3zwH+HxK6bb298cj4j00E//Mh4Er\nUkpfBIiI9wL7gV8Fbo6IbcD7gQtSSl9pt3kf8FBEnJ1Sur+HcmtEETG5wTckzx9oGX2sEO4BzouI\nnwWIiNOB1wN/3v5+KrADuHP2H1JKzwD30YQJwFk0YTW/zcPA43PbZMcJTarPlMZ1HyuEK4FtwF9H\nxAs0ofO7KaVPt8/vABLNimDe/vY5gO3AwTYoNtsmK1PqNNLUpJQmscrqIxDeDbwHuAD4FnAG8ImI\n2JdSuqmH1xudYXB0UxlQpbIPH90U+nAfgXAV8J9SSp9pf/9mRLwKuBS4CXgSCJpVwPwqYTvwjfbf\nTwLHRcS2hVXC9vY5SVLH+jiHcALwwsJjL85eK6X0KM2kft7syfYk8utozj8A7AWeX9hmJ3AKcG8P\nZV6b76y6U/u7r7HYrt2pfbz3sUL4AnBZRDwBfBM4E7gY+IO5ba5tt3kEeAy4AngC+Dw0J5kj4lPA\n1RHxNHAAuA64O6crjGrvHJIOV/Ohoz4C4SKaCf564CRgH/DJ9jEAUkpXRcQJwI3AicBdwFtSSgfn\n/s7FNCuNW4DjgduAC3so71oMA2m6ag2FqGVii4gzaQ41DaKWdhvSMgPIdu2e7d6PgQNhV0rpgb5f\nxM8yWoODZz3LtFuN77rGZBj0p8Z2MxAkSYCBsLIa3xVIWk9t84GBsILadr6kratpXjAQNKiaBk8N\n3B+aZyAsyYEjaTO1zA8GwhJq2dml8EqjbtiOw6phnjAQJEmAgXBUNaR+bmzTPLgfuld6mxoIkiTA\nQDii0tO+ZB7/3hrbbzwlzxsGgiQJMBA2VXLKl8D2HZft369S29dA2ECpO1NSPkqcRwwEZcvj4Oux\n3bQuA2FBiakuKU+lzScGgiQJMBAOUVqal872HoftPqyS2ttAUNY8Hr4a20tbYSBIkgAD4SUlLesk\nlaWU+cVAkCQBBgJQTnrXyLYflu09nhLa3kCQJAEGQhGpPXVeObMc2yl/uc83kw8ESVJj0oGQe1pP\nhfthGLZzHnLeD5MOBEnSywwEFcHj40dm+6gLkw2EnJdtkuqW6/wz2UCQJB1qkoGQazpPmfukX7Zv\nfnLcJ5MMBJXJ4+Qbs13UlckFQo6pLGmacpuPJhcIkqSNGQjKRm7vlmphu2pZkwoEB4ak3OQ0L00q\nEFQ+T6AeyvZQlwwESRIAx45dgKHktCzT5lJKWbzr7aIMOfS5HMqgo8ul308mEKTN9DUQN/q7TtDK\nmYGg4kTElifWsd6Nzb9uqXVQvQwETUpOk+isLK4alItJnFR2wCkisgqDeTmXTcPJYZ5yhaCqlTTR\numLQ2CaxQlBZupgQS37X3VXZDRatqvpAcFDUabMJs+QgWHSkutRSRx1q7Pmq+kDQdNQ6SdZaL+XH\nQFAVap80a6+f8mAgKEvLLp1rOkR0NKvUdexDDypT1YHgoKjbVIJg0VTrPRVjzlu9BEJEnBwRN0XE\nUxHxXEQ8GBFnLmxzeUTsa5+/IyJOW3j++Ii4vv0bByLilog4qY/ySpJ6CISIOBG4G/gh8Gbg54F/\nDTw9t80lwEXAB4GzgWeB2yPiuLk/dS3wNuBdwLnAycBnuy6vJKkRXS9PIuJK4JyU0q8cYZt9wMdT\nSte0v28D9gO/lVK6uf39u8AFKaXPtdvsBB4CfimldP8Gf/NMYO/8Yx4yKpuHRtZn3y/bBn1/V0rp\ngb5ft49DRm8Hvh4RN0fE/oh4ICI+MHsyIk4FdgB3zh5LKT0D3Aec0z50Fs1d1PPbPAw8PrfNETkg\nyuc+XI/tVr6x9mEfgfAzwG8DDwP/DPgkcF1E/Gb7/A4g0awI5u1vnwPYDhxsg2KzbSRJHerjs4yO\nAe5PKX2s/f3BiHgN8CHgph5eb0Pnn3/+S//evXs3u3fvHuqlJWlle/bsYc+ePaOWoY9A+A7Nsf55\nDwG/1v77SSBoVgHzq4TtwDfmtjkuIrYtrBK2t88d1a233rpisSVpPItvXMc4h9bHIaO7gZ0Lj+0E\n/gYgpfQozaR+3uzJ9iTy64B72of2As8vbLMTOAW4t4cyS9Lk9bFCuAa4OyIuBW6mmeg/APzLuW2u\nBS6LiEeAx4ArgCeAz0NzkjkiPgVcHRFPAweA64C7N7rCaJEn1ergVUbr6eIb5TS+Mb5nufNASCl9\nPSLeCVwJfAx4FPhwSunTc9tcFREnADcCJwJ3AW9JKR2c+1MXAy8AtwDHA7cBF3ZdXuXJMNgaQ0Hr\n6Pw+hLHM34dQS52mzEDYOsdB+ebGQbH3IUhbYhh0w3bUqgwESRJgICgzvqvtlu2pVRgIkiTAQFBG\nfDfbD9tVy+rjPgRpJU5Y/Zu1sVce6UiqWyHs3bv36BtpVLOvgpzS11/mwrYvy9DzmSsEDcLJJ0/z\n+8XVgwwE9cYQKIvhIANBnTMIyuc5h2kyENQJQ6BOrhqmxUDQWgyA6Vnc5wZEfQwELc0Q0DxXD/Ux\nEHREhoCWYTjUwUDQYQwBbYXhUC4DQYMHwPwkYfgMb8j297xDWQyEiRozBMYqgxrz36Y2dDi7esib\ngTAhOYSA8mU4yEConCGgdRgO02QgVGiMQzEO4nrN9u1Q/cpwGI+BUAEDQENY3OdDrxw2KoO6ZSAU\nykNBGtvYAWGf7J6BUBBXAsqZh5bKZyBkzhBQaca4z8Rw6IaBkCFDQLUwHMpiIGTCcwKqnecc8mcg\njMgQ0JR5r0N+DISBGQLS4QyHPBgIAzAEpOUZDuMxEHpgAEjdGPq8w9RvhDMQOmIISP1z9dAvA2EL\nDAFpPIZD9wyEFRkCUn4Mh24YCEswBPqRUvJLckZQe/8yHNZnIGzAAJDq4Enp1RgILT8uQqqfH8B3\nZJMOBFcC0jT5MRobm1wguBKQtMiVQ6P6QDAA8uaJ5WHZN49s7JXDRmUY0jGjvXKPIuKln6GklF76\nkVSHMcb1GPPXS69dywQWEWcCe4d8zVrabmyuEIZjn+3GCH12V0rpgb5fpPpDRl1zQHXPw0bDsO92\nZ4wv/hmCgbAEB5KkzdQUDgbCBgyA4blK6Jd9ehhjnJTukoHQcsBI6lppq4dJB4IhkBdXCf2wn+eh\nhHCYXCA4OCSNLddwmEQgGALlcJXQLft+/nIKhyoDwUEgqURjn5Su8k5llc1A74btqFUZCJIkYIBA\niIiPRsSLEXH1wuOXR8S+iHguIu6IiNMWnj8+Iq6PiKci4kBE3BIRJ/VdXuXBd7dbY/tpHb0GQkS8\nFvgg8ODC45cAF7XPnQ08C9weEcfNbXYt8DbgXcC5wMnAZ5d83S2XXeNzUluP7VaHMeax3gIhIn4S\n+GPgA8D/XXj6w8AVKaUvppT+J/Bemgn/V9v/uw14P3BxSukrKaVvAO8DXh8RZ/dVZuXHyW01tpe2\nos8VwvXAF1JKX55/MCJOBXYAd84eSyk9A9wHnNM+dBbNFVDz2zwMPD63jSSpQ71cdhoRFwBn0Ezs\ni3YACdi/8Pj+9jmA7cDBNig220YT4b0Jy3F1oK3qPBAi4qdpjv+/KaX0o67//rLOP//8l/69e/du\ndu/ePVZR1AFD4cgMg/Lt2bOHPXv2jFqGzr8gJyLeAfwp8AIwG8E/RrMqeAF4NfAIcEZK6a/m/t9f\nAt9IKV0cEW8EvgT8/flVQkQ8BlyTUvrEBq972BfkOEjqYygczn5enw36+SBfkNPHOYQvAb9Ic8jo\n9Pbn6zQnmE9PKX0beBI4b/Yf2pPIrwPuaR/aCzy/sM1O4BTg3h7KrEI4+R3K9lCXOj9klFJ6FvjW\n/GMR8SzwvZTSQ+1D1wKXRcQjwGPAFcATwOfbv/FMRHwKuDoingYOANcBd6eU7u+6zCqLh48ahoG6\nNtRnGR3Sc1NKV0XECcCNwInAXcBbUkoH5za7mOYQ0y3A8cBtwIXDFFe5m3ooGAbqQ+fnEMay0TkE\ncODUboqhYJ+u2yZ9epBzCFV+2qmmYzY5TiEYDAL1zQ+3UxVqnyxrr5/yUH0gTOGdoxq1Tpq11kuH\nG3u+8pCRqlLTISSDQEMzEFSlkoPBINBYDARVraRgMAg0turPIUAZk4H6lVLKdsLNuWwaTg7zlCsE\nTUpOKwZDQLmZxApBmpdDGEA+5ZBmDARNSm6TcG7l0bRN5pBRRLhEn7h1Jt91+syqr2PfVC5vDFwh\naBJWHXBbOdG7zv/NZULQtBkIqt46YdAFQ0GlmVQgOOCmp7R9Xlp5tXU57fNJBYKmZcjDRF3+zZwm\nCE2LgaAqjXWYqKu/byhoDJMLBAda/XILg3Vfx75av9z28eQCQXXLNQzWfb3cJgzVbZKB4CCrT0Rk\nHwbrvu46dVP+ctynkwwE1WWoG866NMQNb9KqJhsIDq46lBgGM4bCdOW6HycbCCpfyWEwYygoJwaC\nilRDGMwYCsrFpAPBQVWmHG4465o3sE1Hzvtt0oGg8pRyJdG6DAWNafKB4IAqR+1hMGMo1Cv3fTX5\nQFAZphIGM4aCxmAg4GDK3dTCYMZQqEsJ+8dAULZKuvu4L97VrCEZCC0HUV5quqx0q7wstXyl7A8D\nQdkxDA5nKGgIBoKyYhhszlBQ3wyEOQ6ecdV4w1nXvIGtPCW1v4GgLEz95PGqDAX1wUBY4MAZnmGw\nHkMhf6W1uYGgURkGW2MoqEsGwgYcNMMwDLphKOSpxHY2EDZR4s4siWHQLUMhL6W2r4GgQXn3cX+8\nq1lbZSAcgYOlW95j0D/vVRhfye1pIChbhsF6bDet69ixC5C7iHCAbZGHiIY3a8Nl2362nW2/NSWv\nDsAVgnpW+gCZGvfXtBkIS3CQrMePohifH3UxnBrazUBYUg07e0geJsqLodCvWtrLQFDnDIM8GQo6\nGgNhBQ6QozMM8mYodK+mNjIQVlTTzu+aYVAGQ6E7tbWNgaBOGAZlMRS0EQNhDQ6OQxkGZTIUtqbG\n9jAQtCWGQdkMBc3rPBAi4tKIuD8inomI/RHxuYj4uQ22uzwi9kXEcxFxR0SctvD88RFxfUQ8FREH\nIuKWiDip6/Kua+oDww+pq4cfire6Wuvfxwrhl4HfA14HvAn4ceAvIuInZhtExCXARcAHgbOBZ4Hb\nI+K4ub9zLfA24F3AucDJwGd7KO/aau0UR+OH1NXHD8VbXs31jr4HakS8Evhb4NyU0lfbx/YBH08p\nXdP+vg3YD/xWSunm9vfvAheklD7XbrMTeAj4pZTS/Ru8zpnA3l4rs4kpTXauCurnPt7ciGGwK6X0\nQN8vMsQ5hBOBBPwdQEScCuwA7pxtkFJ6BrgPOKd96CyaD96b3+Zh4PG5bTQwJ4pp8LzCdPUaCNH0\nlGuBr6aUvtU+vIMmIPYvbL6/fQ5gO3CwDYrNtsnGFAaEYTAthsLhplDHvj/++gbgF4DX9/w6o6v5\nY7INg2lKKa207x0D5etthRARvw+8FXhDSuk7c089CQTNKmDe9va52TbHtecSNtsmOzV2GsNg2lwp\n1FmnzfQSCG0YvAN4Y0rp8fnnUkqP0kzq581tv43mqqR72of2As8vbLMTOAW4t48yd6WmzmMYCKYd\nCjXVZRmdHzKKiBuA3cD5wLMRMVsJ/L+U0g/af18LXBYRjwCPAVcATwCfh+Ykc0R8Crg6Ip4GDgDX\nAXdvdIWRumcYaJ6Hj6ahj3MIH6I5afyXC4+/D/gjgJTSVRFxAnAjzVVIdwFvSSkdnNv+YuAF4Bbg\neOA24MIeytu5kgeD9xhoM+uEwuz/lWhqqwMY4D6EoYx5H8JmSmtbw0DLmEI/yTAMqrkPYbIy7FSb\nmsIgVzdqv6u5pLJ2zUDoWQmdyzDQqmoNhRLK2CcDYeIMA62r1lCYMgNhALkOgnWuJDIMNG+dPlHL\neKiRgTCQ3Dqbl5WqS6WHQm7lGYuBMKBcOp1hoD6UGgq5lCMHBsLAxu58hoH6VFoojP36uTEQRjBW\nJzQMNIRSQsEwOJyBMJKhO6NhoCHlHgqGwcYMhBEN0Sn97mONJdfvajYMNmcgjKzPzuk9Bhpbbvcq\nGAZHZiBkoI9OahgoF7mEgmFwdAZCJrrsrN5wptyMfQObYbAcAyEjXXRazxcoZ2OEgmGwPAMhM1vp\nvIaBSjBkKBgGqzEQMrROJzYMVJIhQsEwWJ2BkKl1vplqWYaBctBnKBgG6zEQMrZMpzYMVLI+QsEw\nWJ+BkLkj3axjGKgGXYXCUDe21cxAKMRiRzcMVJOthoJB0A0DoSCzTm8YqEbrhoJh0J1jxy6AVrNK\n5zcIVJpZn122nxsG3XKFUCnDQCWz/47DQKiQg0k1sB8Pz0CojININbE/D8tAqIiDRzWyXw/HQKiI\nJ9hUI/v1cAyEyjh4VBP787AMhAo5iFQD+/HwvA+hUrPB5PFXlcYgGI+BUDmDQaUwCMbnIaOJcLAp\nZ/bPPLhCmBBXC8qNQZAXVwgT5CBUDuyH+XGFMFGuFjQWgyBfBsLEGQwaikGQPwNBgMGg/hgE5fAc\ngg7h4FWX7E9lcYWgw7ha0FYZBGUyELQpg0GrMgjKZiDoqOYHueGgRYZAPQwErcRVg2YMgvoYCFqL\nwTBdBkG9DARticEwHQZB/QwEdcLzDHUyBKbFQFDnXDWUzyCYJgNBvTEYymMQTJuBoN55OClvhoBm\nDAQNynDIgyGgjRgIGs3ipGRA9McA0DL8cDtlIyJe+lnGnj17ei5Rnpat96rtKWUfCBFxYUQ8GhHf\nj4ivRcRrxy6T+jc/mW02oRkIh1qmzaQjyToQIuLdwH8G/i3wj4EHgdsj4pWjFkyDc7I7nG2irkXO\nx20j4mvAfSmlD7e/B/B/gOtSSlctbHsmsHf4UioXOfflrXLCn7xdKaUH+n6RbFcIEfHjwC7gztlj\nqRnxXwLOGatcytfiO+aS3kGXXHbVI+erjF4J/Biwf+Hx/cDODbZ/Re8lUtGWnVj37u1uoblr167O\n/pYmbZD5LedAWNWrxi6A6uAkrgy9Crin7xfJORCeAl4Ati88vh14coPtbwd+A3gM+EGvJZOkYbyC\nJgxuH+LFSjyp/DjNSeWPj1o4SapMzisEgKuBP4yIvcD9wMXACcAfjlkoSapR1oGQUrq5vefgcppD\nRf8DeHNK6bvjlkyS6pP1ISNJ0nCyvQ9BkjQsA0GSBFQUCDV9CF5EXBoR90fEMxGxPyI+FxE/t8F2\nl0fEvoh4LiLuiIjTFp4/PiKuj4inIuJARNwSEScNV5OtiYiPRsSLEXH1wuPV1TsiTo6Im9oyPxcR\nD7YfxzK/TVX1johjIuKKiPh2W6dHIuKyDbarqt5ZSykV/wO8m+beg/cCrwZuBP4OeOXYZVuzPn8O\n/Cbw88AvAl+kub/iJ+a2uaSt478AXgP8GfC/gePmtvlk+/9+hebDAe8B7hq7fku2wWuBbwPfAK6u\nud7AicCjwB/QfFzLPwTeBJxaeb3/DfC3wD8HTgF+DXgGuKjmeuf8M3oBOqkEfA34xNzvATwBfGTs\nsnVUv1cCLwL/ZO6xfcDFc79vA74P/Prc7z8E3jm3zc7275w9dp2OUt+fBB4G/inw3xcCobp6A1cC\nXznKNjXW+wvAf1l47Bbgj2qud84/xR8ymsiH4J0IJJp3SkTEqcAODq3zM8B9vFzns2guK57f5mGa\nG/tyb5frgS+klL48/2DF9X478PWIuLk9RPhARHxg9mTF9b4HOC8ifhYgIk4HXk+zQq653tnK+j6E\nJa36IXhFae/Ovhb4akrpW+3DO2gCYqM672j/vR042A6gzbbJTkRcAJxBM9AX1VrvnwF+m+a7P/4D\ncDZwXUT8MKV0E/XW+0qad/h/HREv0JzT/N2U0qfb52utd7ZqCITa3QD8As07p6pFxE/ThN+bUko/\nGrs8AzoGuD+l9LH29wcj4jXAh4CbxitW794NvAe4APgWzRuBT0TEvjYINbDiDxmx+ofgFSMifh94\nK/CGlNJ35p56kuY8yZHq/CRwXERsO8I2udkF/BTwQET8KCJ+RHOi8MMRcZDmXV+N9f4O8NDCYw/R\nnGiFevf3VcCVKaXPpJS+mVL6E+Aa4NL2+Vrrna3iA6F9J7kXOG/2WHuY5TwG+LjYvrRh8A7gjSml\nx+efSyk9StPZ5+u8DXgdL9d5L/D8wjY7aSaZe3st/Pq+RHNV1RnA6e3P14E/Bk5PKX2bOut9N4cf\n3twJ/A1Uvb9PoHkzN+9F2nmp4nrna+yz2l38AL8OPMehl51+D/ipscu2Zn1uAJ4Gfpnmnc7s5xVz\n23ykrePbaSbRPwP+F4dejncDzeWMb6B59303hV2Ox+FXGVVXb5rzJT+keWf8j2gOoxwALqi83v+V\n5uTvW2kutX0nzWWo/7Hmeuf8M3oBOqsI/A7Ntcjfp3lncNbYZdpCXV6keee0+PPehe3+Hc1lec/R\nfF76aQvPHw/8Hs1htQPAZ4CTxq7fim3x5flAqLXe7aT4V22dvgm8f4Ntqqo38PdoPtH4UeDZdqL/\n98CxNdc75x8/3E6SBFRwDkGS1A0DQZIEGAiSpJaBIEkCDARJUstAkCQBBoIkqWUgSJIAA0GS1DIQ\nJEmAgSBJav1/x1BoBIMpvkAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x19f36d510>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "telap_fname_tail = \"TelAp_half_wfirstCycle5_N{:04d}.dat\".format(D/2)\n",
    "telap_fname = os.path.join(telap_dir, telap_fname_tail)\n",
    "telap_rhs = np.loadtxt(telap_fname)\n",
    "wfirst_telap = np.concatenate([telap_rhs[:,::-1], telap_rhs], axis=1)\n",
    "plt.imshow(wfirst_telap)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define annular stop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 461,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "xs = np.linspace(-D/2 + 0.5, D/2 - 0.5, D)\n",
    "[Xs, Ys] = np.meshgrid(xs, xs)\n",
    "inside_ann = np.less_equal(Xs**2 + Ys**2, (iD/100.*D/2)**2)\n",
    "outside_ann = np.greater_equal(Xs**2 + Ys**2, (oD/100.*D/2)**2)\n",
    "annular_stop = np.ones((D,D))\n",
    "annular_stop[inside_ann] = 0.\n",
    "annular_stop[outside_ann] = 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 462,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x1a085e210>"
      ]
     },
     "execution_count": 462,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAFkCAYAAAAg+UxRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAGI1JREFUeJzt3X+sJWV9x/H3l+KuUrMhKbIbYihY6qrRYAFBQkUtpBat\ntZZEqEaNRqkKCSFNQCqmLaQtwRRWFAyppEa0MbgGCcZAEVuDC7JhsTQFJKWAhCy7gtLuhlWWH9/+\nMXPlcPbe3ftj5szMM+9XcpM9M8/eM8/M8zyfeWbmnBuZiSRJ+3W9AZKkfjAQJEmAgSBJqhkIkiTA\nQJAk1QwESRJgIEiSagaCJAkwECRJNQNBkgQYCJKkmoEgSQIMBElSzUCQJAGwf9cb0JSI+C3gHcDD\nwK+63RpJasRLgcOAmzLz522/WTGBQBUGX+96IySpBR8A/qXtNynpktHDXW+AJLXk4Vm8SUmB4GUi\nSaWayfhWUiBIklbAQJAkAQaCJKlmIEiSAANBklQzECRJgIEgSaoZCJIkwECQJNUMBEkSYCBIkmoG\ngiQJKOvrr6VFycyuN+HXIqLrTZB+zUBQsfo08C9koW00KNQFA0FFGMLgvxTz1ceQUNsMBA1SaQGw\nGNN1NiDUNANBgzDGANgXA0JNMxDUSwbA0hkQWikDQb1hCDRrcn8aDloMA0GdMgRmw3DQYhgImjlD\noFuGgxZiIGhmDIL+mTsmBoPAQFDLDIFhcNYgMBDUEoNguJw1jJeBoEYZBOUwGMbHQFAjDIJyGQzj\nYSBo2QyBcfE+Q/n8ewhaFsNg3Dz+ZXKGoCVxINAcLyWVx0DQohgEWojBUA4vGWmfDAMthu1k+Jwh\naF52bi2HN56HzRmC9mAYqAm2o+FxhqBfswOrad5fGBZnCAIMA7XL9jUMzhBGzo6qWXG20H/OEEbM\nMFAXbHf95QxhhOyQ6pqzhX5yhjAyhoH6xPbYLwbCiNj51Ee2y/7wktEI2OHUd15C6gdnCIUzDDQk\nttduGQgFs3NpiGy33TEQCmWn0pDZfrvhPYTC2JFUCu8rzJ4zhIIYBiqR7Xp2DIRC2GlUMtv3bBgI\nBbCzaAxs5+0zEAbOTqIxsb23y5vKA2XH0Fh5s7k9zhAGyDCQ7AdtMBAkSYCBMDieFUkvsD80y0AY\nEBu/tCf7RXMaD4SI2C8iLoqIByNiV0Q8EBEXzFPuwojYWpe5OSKOmFq/OiKuiIgnImJnRGyMiIOb\n3t4hyEwbvbQX9pFmtDFD+DTwF8CngNcA5wLnRsRZcwUi4jzgLOAM4FjgKeCmiFg18Xs2AO8CTgVO\nBA4BvtXC9vaajVxaPPvLykTTOzAibgC2ZebHJ5ZtBHZl5ofq11uBz2XmZfXrNcB24MOZeW39+nHg\n9My8ri6zHrgPeHNmbp7nfY8CtjRamR6wgUtLU+jjqEdn5l1tv0kbM4TbgJMi4ncBIuJI4ATgu/Xr\nw4F1wC1z/yEzdwB3AMfXi46h+ozEZJn7gUcmyhTPMJCWzn6zfG18MO1iYA3wk4h4jip0PpOZ36jX\nrwOSakYwaXu9DmAtsLsOioXKFM1GLS1fZpY6U2hVG4FwGvB+4HTgXuCNwOcjYmtmXtPC+xXHMJBW\nzlBYujYC4RLgHzLzm/XreyLiMOB84BpgGxBUs4DJWcJa4Mf1v7cBqyJizdQsYW29rliGgdQcQ2Fp\n2riHcADw3NSy5+feKzMfohrUT5pbWd9EPo7q/gNUN4efnSqzHjgUuL2Fbe4Fw0Bqnv1q8dqYIdwA\nXBARjwL3AEcB5wBfniizoS7zAPAwcBHwKHA9VDeZI+Jq4NKIeBLYCVwObJrvCSNJ0sq1EQhnUQ3w\nVwAHA1uBL9XLAMjMSyLiAOAq4EDgVuCUzNw98XvOoZppbARWAzcCZ7awvb3gWYzUHi8dLU7jn0Po\nypA/h1DKMZD6bsChMNjPIWgJDANpduxve2cgdMjGKc2e/W5hBoIkCTAQOuNZitQd+9/8DIQO2Bil\n7tkP92QgzJiNUOoP++OLGQiSJMBAmCnPRqT+sV++wECYERud1F/2z4qBIEkCDISZ8OxD6j/7qYHQ\nOhuZNBxj768GgiQJMBBaNfazDWmIxtxvDYSWjLlRSUM31v5rIEiSAANBklQzEFow1ummVJIx9mMD\noWFjbERSqcbWnw0ESRJgIDRqbGcT0hiMqV8bCJIkwEBozJjOIqSxGUv/NhAaMJbGIo3ZGPr5/l1v\ngDQEEbHPMmMYMFQ2A0GqLWbQX+7/Nyw0BAbCCtnRh2mlg38T72fbGZ7MnHnbmSUDQaPRt448uT2G\ng/rAQFgBO3H/9S0EFmI4DEfJswQDQUUacoed23aDQbPmY6fLZGftp4gYdBhMKqkupSm1/ztD0OCV\nPmh6OUmz4gxhGeyU/VF6GEwbW337rMRxwBmCBmnMA6P3GNQWZwganDGHwST3g5pmICyRZ2Xd8Sbr\nntwn3SptPPCSkXrPAW/fvIykJjhDWAI72+wZBkvj/pq9ksYFA0G95eC2PO43LZeBoF5yUFsZ95+W\nw0BYpJKmhX3nYNYM9+PslDI+GAjqFQexZrk/tRQGgnrDwasd7lctloGwCKVMB/vMQatd7t/2lTBO\nGAjqnIPVbLiftS8GgjrlIDVb7m/tjYGwDyVMA/vKwakb7vf2DH28MBDUCQelbrn/NR8DQTPnYNQP\nHgdNMxA0Uw5C/eLx0CQDYS+Gfj1Q0uwNedxoJRAi4pCIuCYinoiIXRFxd0QcNVXmwojYWq+/OSKO\nmFq/OiKuqH/HzojYGBEHt7G9mg3PRvvJ46I5jQdCRBwIbAKeBt4BvBb4S+DJiTLnAWcBZwDHAk8B\nN0XEqolftQF4F3AqcCJwCPCtprdXklSJpqc3EXExcHxmvnUvZbYCn8vMy+rXa4DtwIcz89r69ePA\n6Zl5XV1mPXAf8ObM3DzP7zwK2NJkXYY89esbz0L7z/benBba+9GZeVfTv3RaG5eM3g3cGRHXRsT2\niLgrIj42tzIiDgfWAbfMLcvMHcAdwPH1omOo/prbZJn7gUcmyrTKztEcw2AYPE7NGer40UYgvAr4\nJHA/8IfAl4DLI+KD9fp1QFLNCCZtr9cBrAV210GxUBlJUoPa+JvK+wGbM/Oz9eu7I+L1wCeAa1p4\nP/WYZ53DEhGDPbvVyrUxQ3iM6lr/pPuAQ+t/bwOCahYwaW29bq7MqvpewkJl1HOGwTB53MarjUDY\nBKyfWrYe+ClAZj5ENaifNLeyHviPA26rF20Bnp0qs54qVG5vYZtfxDMkSSs1xHGkjUtGlwGbIuJ8\n4Fqqgf5jwMcnymwALoiIB4CHgYuAR4HrobrJHBFXA5dGxJPATuByYNN8TxhJklau8UDIzDsj4r3A\nxcBngYeAszPzGxNlLomIA4CrgAOBW4FTMnP3xK86B3gO2AisBm4Ezmx6e9UOLzsMm/cSxqnxzyF0\npcnPIZSyT7pkIAyf/WDlGuwHg/0cgkbOMCiDx3F8DIQpnhVJasrQxhMDQY3yrLIsHs9xMRAkSYCB\nIEmqGQhqjJcXyuRxHQ8DQZIEGAgvMrQnAiT135DGFQNBjfCyQtk8vuNgIEiSAANBklQzECRJgIEg\nSaoZCFoxbziOg8e5fAaCJAkwECRJNQNBkgQYCJKkmoEgSQIMBElSzUDQivgo4rh4vMtmIEiSAANB\nklQzECRJgIEgSaoZCJIkwECQJNUMBEkSYCBIkmoGglYkM7veBM2Qx7tsBoIkCTAQJEk1A0GSBBgI\nkqSagSBJAgwESVLNQNCK+SjiOHicy2cgSJIAA0GSVDMQJEmAgSBJqhkIaoQ3HMvm8R0HA2FCRHS9\nCZIKM6RxxUCQJAEGghrkZYUyeVzHw0CQJAEGgiSpZiCoUV5eKIvHc1wMhClDeiJAUr8NbTwxENQ4\nzyrL4HEcHwNBkgTMIBAi4tMR8XxEXDq1/MKI2BoRuyLi5og4Ymr96oi4IiKeiIidEbExIg5ue3vV\nDM8uh83jN06tBkJEvAk4A7h7avl5wFn1umOBp4CbImLVRLENwLuAU4ETgUOAb7W5vZI0Zq0FQkS8\nHPga8DHgf6dWnw1clJnfycz/Aj5ENeD/af1/1wAfBc7JzB9k5o+BjwAnRMSxbW3zxLa3/RaSCjfE\ncaTNGcIVwA2Z+f3JhRFxOLAOuGVuWWbuAO4Ajq8XHQPsP1XmfuCRiTLqOS87DJPHbbz2b+OXRsTp\nwBupBvZp64AEtk8t316vA1gL7K6DYqEyGoDMHOSZ0lgZBuPWeCBExCuprv+fnJnPNP37JUntaOOS\n0dHAK4C7IuKZiHgGeCtwdkTspjrLD6pZwKS1wLb639uAVfW9hIXKtMqz2uZ41jkMHqfmDHX8aCMQ\nvge8geqS0ZH1z51UN5iPzMwHqQb1k+b+Qz3wHwfcVi/aAjw7VWY9cChwewvbrJY52PSbx0fQwiWj\nzHwKuHdyWUQ8Bfw8M++rF20ALoiIB4CHgYuAR4Hr69+xIyKuBi6NiCeBncDlwKbM3Nz0NkuSWrqp\nPI8XnX5k5iURcQBwFXAgcCtwSmbunih2DvAcsBFYDdwInDmbzVUbvMHcT84ONCdKaQwRcRTVpaZG\nlbJ/+sJA6B/beLNaauNHZ+ZdbfziSX6XkWbKwadfPB6aZCBo5hyE+sHjoGkGgjrhYNQt97/mYyDs\ng9e82+Og1A33e3uGPl4YCOqUg9Nsub+1NwaCOucgNRvuZ+2LgbAIQ58GDoGDVbvcv+0rYZwwENQb\nDlrtcL9qsQwE9YqDV7Pcn1oKA2GRSpgODoWDWDPcj7NTyvhgIKiXHMxWxv2n5TAQ1FsOasvjftNy\nGQhLUMq0cEgc3JbG/TV7JY0Ls/r6a2nZ5ga5kjpe0wwCNcEZwhI5KHUnMx34prhPulXaeGAgaHAc\nACvuBzXNS0YapDFfRjII1BZnCMswxkGor8Y2OI6tvn1W4jjgDEGDNzlIlthJDQHNijOEZSpx4ClB\nSTdZS6pLaUrt/84QVKQh32MwBNQVA2EFIsLO23NDuZxkOxqOPrejlTIQNBp9CwdDQH1jIKyQs4Rh\nmu+YtRkStpEy9OFEok0GglTb26C9mIHAQV9DZyBIi+BgrzHwsdMGlD6NlDSOfm4gNGQMjUUaq7H0\nbwNBkgQYCI0ay1mENCZj6tcGgiQJMBAaN6azCal0Y+vPBkILxtaIpBKNsR8bCJIkwECQJNUMhJaM\ncboplWKs/ddAaNFYG5U0ZGPutwaCJAkwEFo35rMNaWjG3l8NhBkYeyOThsB+aiBIkmoGwox49iH1\nl/2zYiDMkI1O6h/75QsMBEkSYCDMnGcjUn/YH1/MQOiAjVDqnv1wTwZCR2yMUnfsf/MzECRJgIHQ\nKc9SpNmz3y3MQOiYjVOaHfvb3hkIPWAjldpnP9u3xgMhIs6PiM0RsSMitkfEdRHx6nnKXRgRWyNi\nV0TcHBFHTK1fHRFXRMQTEbEzIjZGxMFNb29f2Fil9ti/FqeNGcJbgC8AxwEnAy8B/jUiXjZXICLO\nA84CzgCOBZ4CboqIVRO/ZwPwLuBU4ETgEOBbLWyvJAmIzGz3DSIOAn4GnJiZP6yXbQU+l5mX1a/X\nANuBD2fmtfXrx4HTM/O6usx64D7gzZm5eZ73OQrY0mplZqDt4yGNTSGzg6Mz866232QW9xAOBBL4\nBUBEHA6sA26ZK5CZO4A7gOPrRccA+0+VuR94ZKJMkQppvFIv2J+WptVAiOpobAB+mJn31ovXUQXE\n9qni2+t1AGuB3XVQLFSmWDZiaeXsR0u3f8u//0rgdcAJLb9PcSLCy0fSMhkGy9PaDCEivgi8E3hb\nZj42sWobEFSzgElr63VzZVbV9xIWKlM8G7W0dPab5WslEOoweA/w9sx8ZHJdZj5ENaifNFF+DdVT\nSbfVi7YAz06VWQ8cCtzexjZL0tg1fskoIq4E/hz4E+CpiJibCfxfZv6q/vcG4IKIeAB4GLgIeBS4\nHqqbzBFxNXBpRDwJ7AQuBzbN94RRybx0JC2es4OVaeMewieobhr/+9TyjwBfBcjMSyLiAOAqqqeQ\nbgVOyczdE+XPAZ4DNgKrgRuBM1vY3t6ba+QGgzQ/g6AZrX8OYVZK+RzCvpRyvKSmjCQMivkcgho0\nksYvLYr9oVkGgiQJMBAGybMiyX7QhrY/mKaWeKNZY2UQtMcZwsDZOTQmtvd2GQgFsJNoDGzn7TMQ\nCmFnUcls37NhIBTETqMS2a5nx5vKhfFms0phEMyeM4RC2Zk0ZLbfbhgIBbNTaYhst90xEApn59KQ\n2F675T2EEfC+gvrOIOgHZwgjYqdTH9ku+8NAGBk7n/rE9tgvXjIaIS8hqWsGQT85QxgxO6W6YLvr\nL2cII+dsQbNiEPSfMwQBdla1y/Y1DM4Q9GvOFtQ0g2BYnCFoD3ZiNcF2NDzOEDSvyc7sjEGLZQgM\nmzME7ZOdXIthOxk+ZwhaFO8vaCEGQTkMBC2JwaA5BkF5vGSkZXEwGDePf5mcIWjZvPE8LoZA+QwE\nNcJLSeUyCMbDQFCjDIZyGATjYyCoFQbDcBkE42UgqFXeZxgGQ0BgIGiGnDX0j0GgSQaCZs5ZQ7cM\nAS3EQFCnDIfZMAS0GAaCesNwaJYhoKUyENRL04OZAbFvBoBWykDQIBgQezIA1DQDQYM0xoAwANQ2\nA0FFmG+wHHJIOPirCwaCirXQoNqnoHDgV58YCBodB2Fpfv49BEkSYCBIkmoGgiQJMBAkSTUDQZIE\nGAiSpJqBIEkCDARJUs1AkCQBBoIkqdb7QIiIMyPioYj4ZUT8KCLe1PU2SVKJeh0IEXEa8I/AXwO/\nB9wN3BQRB3W6YZJUoF4HAnAOcFVmfjUzfwJ8AtgFfLTbzZKk8vQ2ECLiJcDRwC1zy7L63uLvAcd3\ntV2SVKreBgJwEPAbwPap5duBdfOUf2nrWyRJ3ZjJ+NbnQFiqw7reAElqyWGzeJM+/4GcJ4DngLVT\ny9cC2+YpfxPwAeBh4FetbpkkzcZLqcLgplm8WfTpzwlOi4gfAXdk5tn16wAeAS7PzM91unGSVJg+\nzxAALgW+EhFbgM1UTx0dAHyly42SpBL1OhAy89r6MwcXUl0q+g/gHZn5eLdbJknl6fUlI0nS7JT0\nlJEkaQUMBEkSUFAglPQleBFxfkRsjogdEbE9Iq6LiFfPU+7CiNgaEbsi4uaIOGJq/eqIuCIinoiI\nnRGxMSIOnl1NViYiPh0Rz0fEpVPLi6t3RBwSEdfU27wrIu6OiKOmyhRV74jYLyIuiogH6zo9EBEX\nzFOuqHr3WmYO/gc4jeqzBx8CXgNcBfwCOKjrbVtmfb4LfBB4LfAG4DtUn6942USZ8+o6/jHweuDb\nwP8AqybKfKn+f2+l+nLA24Bbu67fIvfBm4AHgR8Dl5Zcb+BA4CHgy1Rf1/LbwMnA4YXX+6+AnwF/\nBBwK/BmwAzir5Hr3+afzDWikEvAj4PMTrwN4FDi3621rqH4HAc8Dvz+xbCtwzsTrNcAvgfdNvH4a\neO9EmfX17zm26zrto74vB+4H/gD4t6lAKK7ewMXAD/ZRpsR63wD809SyjcBXS653n38Gf8loJF+C\ndyCQVGdKRMThVN/nNFnnHcAdvFDnY6geK54scz/VB/v6vl+uAG7IzO9PLiy43u8G7oyIa+tLhHdF\nxMfmVhZc79uAkyLidwEi4kjgBKoZcsn17q1efw5hkfb2JXjrZ785zao/nb0B+GFm3lsvXkcVEHv7\n4r+1wO66Ay1Upnci4nTgjVQdfVqp9X4V8Emqv/3xd8CxwOUR8XRmXkO59b6Y6gz/JxHxHNU9zc9k\n5jfq9aXWu7dKCITSXQm8jurMqWgR8Uqq8Ds5M5/pentmaD9gc2Z+tn59d0S8nurvf1zT3Wa17jTg\n/cDpwL1UJwKfj4itdRBqxgZ/yYilfwneYETEF4F3Am/LzMcmVm2juk+ytzpvA1ZFxJq9lOmbo4FX\nAHdFxDMR8QzVjcKzI2I31VlfifV+DLhvatl9VDdaodzjfQlwcWZ+MzPvycyvA5cB59frS613bw0+\nEOozyS3ASXPL6sssJ1FdoxykOgzeA7w9Mx+ZXJeZD1E19sk6rwGO44U6bwGenSqznmqQub3VjV++\n71E9VfVG4Mj6507ga8CRmfkgZdZ7E3te3lwP/BSKPt4HUJ3MTXqeelwquN791fVd7SZ+gPdR/WnN\nycdOfw68outtW2Z9rgSeBN5CdaYz9/PSiTLn1nV8N9Ug+m3gv3nx43hXUj3O+Daqs+9NDOxxPPZ8\nyqi4elPdL3ma6sz4d6guo+wETi+83v9MdfP3nVSP2r6X6jHUvy+53n3+6XwDGqsIfIrqWeRfUp0Z\nHNP1Nq2gLs9TnTlN/3xoqtzfUD2Wt4vq+9KPmFq/GvgC1WW1ncA3gYO7rt8S98X3JwOh1HrXg+J/\n1nW6B/joPGWKqjfwm1TfaPwQ8FQ90P8tsH/J9e7zj19uJ0kCCriHIElqhoEgSQIMBElSzUCQJAEG\ngiSpZiBIkgADQZJUMxAkSYCBIEmqGQiSJMBAkCTV/h+DFG9+QJKw6QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x19d6ae150>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(annular_stop)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Pad the telescope aperture"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 463,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "max_shift = int(round(D*float(pad)/1000))\n",
    "shift_range = range(-max_shift,max_shift+1,1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 464,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(39, 39)"
      ]
     },
     "execution_count": 464,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "[Xshifts, Yshifts] = np.meshgrid(shift_range, shift_range)\n",
    "allowed_shifts = np.less_equal(Xshifts**2 + Yshifts**2, max_shift**2)\n",
    "allowed_shifts.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 465,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1129"
      ]
     },
     "execution_count": 465,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "XYshifts_allowed = zip(Xshifts[allowed_shifts], Yshifts[allowed_shifts])\n",
    "len(XYshifts_allowed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 466,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "padded_wfirst_mask = np.ones(wfirst_telap.shape)\n",
    "for (xshift,yshift) in XYshifts_allowed:\n",
    "    wfirst_telap_shifted = np.roll(np.roll(wfirst_telap, yshift, 0), xshift, 1) \n",
    "    padded_wfirst_mask *= wfirst_telap_shifted"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 467,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x1a18b3790>"
      ]
     },
     "execution_count": 467,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1YAAAM1CAYAAABzPGCfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X/IfdldH/r3SsfEhjIErs0MoYgW22lLS2yi0SCNLZEr\nVvoL/zC1VFE0tU0g5BbU0Egh4XJDxGSaNkqggpj2ImmCBEtJqiktGjUhExu4GgPBpLlhnKmxuU2Y\nVCfGdf8453HO95nnec55nr332Wut/XpBYL7nOfl+97PPWp+13nutvU+ptQYAAIC7e8baBwAAANA7\nwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoAAGAiwQoA\nAGAiwQoAAGAiwQoAAGAiwQoAAGCi+9Y+gLmUUv63JN+a5JNJfn/dowEAAFb05Um+Ksl7a62/d45/\ncJhglV2o+rdrHwQAANCMf5Dk/z7HPzTSVsBPrn0AAABAUz55rn9opGBl+x8AAHDobBlhpGAFAACw\nCsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEKAABgIsEK\nAABgIsEKAABgovvWPgAAxlNrXfsQVlVKWfsQADgzwQqAW9l6aDrFKedI+AIYi2AFwJUEqGVdd34F\nLoA+CVYAGyEo9eG2n5MgBtAGwQpgQELUdlz1WQtbAOcnWAEMQJDi0OX2IGgBLE+wAuiQIMVtCFoA\nyxOsABolPLGUm9qW0AVwN4IVQEOEKdZ22AaFLIDTCVYAKxKkaJkthACnE6wAzkCAYgS+ewvgeoIV\nwEKEKbbC9kEAwQpgNoIU2D4IbJdgBXAHQhScxhcYA1shWAHcgkAF0130IwELGIlgBXCEMAXLcG8W\nMBLBCuAKwhScl5AF9E6wAoggBS3xAAygR4IVsFnCFPTBahbQg2esfQAAaxCqoE/6LtAqwQoAAGAi\nwQoAAGAiwQoAAGAiwQoAAGAiwQrYJE8Wgz7pu0CrBCsAAICJBCsAAICJBCtgs2wpgr7os0DLBCsA\nAICJBCsAAICJBCsAAICJBCtg09yzAX3QV4HWCVYAAAATCVYAAAATCVbA5tliBG3TR4EeCFYAAAAT\nCVYAAAATCVYAAAATCVYAAAATCVYAAAAT3bf2AQB9qbX+8X97UhfA7aihMC7BCjjqcCIAwDyELBiL\nrYDAtWqtmwlVJjXQpq30zS3VWxiVYAVcyQAPcH5qL/RLsALucZurpiYAAKe7TW1VX6E/7rECDOB7\npRTnAhqylW2A13EPFvTFihVsnCAB0D61GtonWMGGzTFQG+wBjlNvYXy2AsIGGZwB+nRRv20NhPZY\nsYINcUP0cSYr0AZ98WbqObTHihVsgMEXYExWsKAdVqxgcEIVwPjUelifFSsYlEEWYFusXsG6BCsY\njEAFsG0CFqzDVkAYyFqhSpgDuJ7aDNtgxQoGYPAE4CpWr+B8rFhBxzxudxkmILAufXB+xgtYnmAF\nAAAwka2A0JGWrzbWWl1lBriktbp9+XjUbZiPFSsAAICJBCvogL3x5+cqLqxD3zsv4wvMR7CCxvU0\n4PV0rABL66km9nSs0Cr3WEGjDHIAnJNHs8M0VqygMbZltMPkAs5Ln2uDcQjuRrACAACYyFZAaMBI\nVwY9dh1gjLp++Duo63CcFStY2QiDLwBjM1bBcYIVrMhA1T5XaeE89LX2GbPgZoIVrMQABUBvjF1w\nPfdYwZkZlADomceyw9UEKzgTgQqAkQhYcC9bAQEAACayYgUL2uoqlUeuA1u2tdrvseywY8UKAABg\nIsEKFlBr3dwVy5G5AgvL0sfGYfxjywQrmJkBZcd5ALZI7dtxHtgiwQpmZCABgB1jIlsjWMEMbH0Y\nn61KsAx9a2zGR7ZEsIKJDBgAcDNjJVsgWMEduQp3nPMDbImadzPjJqMTrOAODAwAcDfGUEYlWMEt\nGRC2y70gMC99aruMpYxIsAIAAJjovrUPAHrgyhoAzOtwbLV6yQisWAEAAEwkWMERVqsAYFnGWkZg\nKyBcQYGfT63VFg9geMaN6WwNpHdWrAAAACYSrOASVx25iauoMA99iZsYi+mRYAUHFPJlOK/AyNS4\nZTiv9Eawgj0FHADaYmymJ4IVm1drVbi5FVuYYBp9iNswTtMLwYpNU6gBoA/GbFonWLFZCvR5Od/A\niNS283K+aZlgxSYpzADQJ2M4rRKs2BwFmTm4RwTuRt9hDsZyWiRYsSkK8bqcf2Akatq6nH9ac9/a\nBwDnoPgCwHguxncrobTAihXDE6pYioEcbkefYSnGelogWDE0hRYAtsGYz9oEK4alwALAthj7WZNg\nxZAUVgDYJnMA1iJYMRwFFQC2zVyANQhWAAAAEwlWDKPW6gpVB3xGwAjUsvaZF3BughVDUDhZi8dH\nw2n0FdZijsC5CFZ0T8EEAG5irsA5CFZ0yxJ/v3xuQM/UsD6ZN7A0wYouKYwAwF2YQ7AUwQpgIveO\nwM30EWALBCu640rTGHyOQI/UrjH4HFmCYEVXFEIAYA7mFMxNsKIbCiAts9UJrqZv0DJzC+YkWAEA\nAEx039oHAMe4mgQALOVinmF1lamsWNE0oQoAOAdzDqYSrGiWAgcAnJO5B1MIVjRJYdsGnzPQEzVr\nG3zO3JVgRXMUNABgTeYi3IVgRVMUMnrmxme4lz5Bz8xJuC3BCgAAYCLBima4MrRNPnegB2rVNvnc\nuQ3BiiYoXABAi8xROJVgxeoULEbinhLY0RcYibkKpxCsWJVCRaIdAG1To0i0A44TrFiNAgUA9MTc\nhZsIVqxCYWJktkCxdfoAIzOH4TqCFQAAwESCFWfnSg9X0S6AFqlNXEW74CqCFWelEAEAIzCn4TLB\nirNRgNgS95iwVdo+W2Juw6H71j4AxqfoAACjupjnuKiAFSsAAICJBCsAAICJBCsWZRsgALAF5jzM\nHqxKKc8opby+lPLbpZQvlFI+Xkp57RXve10p5dH9e36hlPI1l37+rFLKW0spnymlfL6U8s5SynPn\nPl6WUWtVYLg1bQZoiZrEbZn/bNsSK1Y/kuQfJfknSf5Ckh9K8kOllFdevKGU8sNJXpnk5UlelOSJ\nJO8tpTzz4O95OMm3J/mOJC9J8rwk71rgeAEAACZZ4qmAL07y7lrre/Z//lQp5buyC1AXXpXk9bXW\nf58kpZTvTvJ4kr+b5B2llPuTfF+Sl9Va/8v+Pd+b5KOllBfVWj+4wHEDzKqU4solm+KpaMCWLbFi\n9StJXlpK+XNJUkp5fpJvSvIf9n/+6iQPJnnfxf+h1vq5JB/ILpQlyddlF/oO3/OxJJ86eA8NsgQO\nAGyd+dA2LbFi9YYk9yf5rVLKl7ILb/+s1vqz+58/mKRmt0J16PH9z5LkgSRP7gPXde+hMQoIc6i1\nuuoNrM6YxhyMaduyRLD6ziTfleRlSX4zydcm+RellEdrrW9f4N+jAQYgAICnE662Y4lg9cYk/1et\n9d/t//wbpZSvSvKaJG9P8liSkt2q1OGq1QNJfn3/348leWYp5f5Lq1YP7H8G0AX3WbEVJo7A1i1x\nj9Wzk3zp0mt/dPFv1Vo/kV04eunFD/cPq/iG7O7PSpJHkvzhpfc8lOQrk/zqAsfMBCaNzE2bAtak\nBjE3bWobllix+vkkry2lfDrJbyR5QZJXJ/nXB+95eP+ejyf5ZJLXJ/l0kncnu4dZlFJ+KsmbSimf\nTfL5JG9J8n5PBGyLQgEAcJwtgeNbIli9Mrug9NYkz03yaJKf3L+WJKm1vrGU8uwkb0vynCS/lOTb\naq1PHvw9r85u5eudSZ6V5D1JXrHA8QIsynZARmeyCJCUUQb7UsoLsttCyJmM0nZo02gTNf2Fkemv\ncLrR+ksHXlhr/fA5/qEl7rFiAww6AAC3Zw41LsGKW1MQAADuzlxqTIIVAADARIIVJ6u1usLC2Whr\nwBrUHs7FvGo8ghUAAMBEghXAGXgKFKPStgF2lvgeKwZjmRoAYBkX8ywXKfpnxQpollAPnJOaA0wh\nWHEjgwwAwPLMufonWAGciW0ejEabBniKYMW1XDmhBdohcA5qDS3QDvsmWHElHRsA4PzMwfolWPE0\nOjQsx9YpRqEtw3LMxfokWHEPHRkAYH3mZP0RrIDmGVyAJakxwBwEKwAAgIkEK5Lsrta5Ygfn4d4U\neqcNw3mYn/VFsAIAAJhIsMKVEACAhpmr9UGwAgAAmEiw2jhXQAAA2mfO1j7BasN0UHqivQJLUFvo\nifbaNsEKAABgIsFqgzy6E9bncdX0StuFdZnHtUuwArphIAHmpKYAcxKsAAAAJhKsNsbVOQCA/pnT\ntUew2hAdENriXhV6o81CW8zt2iJYAV0xiABzUEuAuQlWG2EAAQAYjzleOwQrgBXZWkUvtFWAmwlW\nG+BKBgDAuMz12iBYDU5HAwAYnznf+gQrAACAiQSrgblyAQCwHeZ+6xKsgO4YOIAp1BBgCYLVoAwa\nAADbYw64HsFqQDoU9MVjrGmdNgp9MRdch2AFdMmgAdyF2gEsRbAajAEDAABzwvMTrAAAACYSrAbi\nygT0yz0stErbhH6ZG56XYAV0y4AB3IaaASxJsBqEwQIAgMvMEc9HsBqADgNjsOWK1miTMAZzxfMQ\nrICuGSyAU6gVwNIEq84ZKAAAOMaccXmCFQAAwESCVcdceYDxuKeFVmiLMB5zx2UJVgAAABMJVp1y\nxQEAgNsyh1yOYNUhHQIAgLsyl1yGYAV0zwAB3ESNAM5BsOqMwQEAgKnMKecnWAEAAEwkWAE0xmOu\nWZs2CHB7glVHLNnC9fQP4CpqA1xP/5iXYAUAADCRYNUJVxRgW2zFYi3aHmyLOeZ8BKsOaPBwGn0F\nOKQmwGn0lXkIVgAAABMJVgAAABMJVo2zNAvb5V4Xzk2bg+0y55xOsAIAAJhIsGqYKwcAAJyLuec0\nglWjNGwAAM7NHPTuBCtgKAYEIFELgPMTrBpkMAAAYC3moncjWAEAAEwkWAE0zOOvORdtDWAawaox\nll5hOv0Itk0NgOn0o9sTrAAAACYSrBriygBwFVu0WJo2BlzF3PR2BCsAAICJBKtGuCIA89KnYJv0\nfZiXPnU6wQoAAGAiwQqgA+6BYSnaFsA8BKsGWGKFZehbsC36PCxD3zqNYAUAADCRYLUyVwCAU9my\nxdy0KeBU5qzHCVYAAAATCVYrkvwBAOiFuevNBKuVaJgAAPTGHPZ6ghUAAMBEghUwNFfWYBv0dWBt\nghUAAMBEgtUKXFUD7srjsZmLtgTclbns1QQrAACAiQQrYHiurMHY9HGgBYLVmSn+AAD0zpz26QSr\nM9IAgTm4N4aptCFgDua29xKsgE1Q/GFM+jbQivvWPgAAuK2pKy4m4wDMTbAC6FApZZhwsMa2tLv8\nm843ADcRrM5klAEZYIqeJ/VXHbvaDmxdrbXr2j4nwQqARY084B7+bkIWwLYJVmdgsAVGN3J4OtVN\n58A4AIzMqtWOYAXAnRhET2dlC2B8Hre+MAMoMJpSilA1gfMHjMicV7ACNkTRn0YgmJfzOZ0+DbTE\nVsAFKfjAkpZ+5LpJ/3lcPs8+U6BXW7/XSrAC4I9teUBshfuxAPpkKyCwKSaqV7MtrU0+l+vpy0Br\nrFgBbJQJez+sYgG0z4oVQMfuEo6sgvTtrp+fzxxgWYLVQlxRBFpkcj0OnyXQoi3PgW0FBDZna08t\n2tLvujVb3SK4pd8V6IcVqwUo+MA5XRecbPnblps+b+0AOKetzoUFK4ABmUhvl88eYB22AgIMxKSa\n5Kl2sNWrxgBrEKxmZhAD1iBQcRXtAljL1u5nTmwFBAAAmEywAgAAmEiwAjbJtl3ok74LtEqwmpFi\nDwAAO1ubGwtWAAAAEwlWAAAAEwlWM9naUieMQL+Fvuiz0J8t9VvBCgAAYCLBCgAAYCLBagZbWuKE\n0ei/0Ad9Ffq1lf4rWAEAAEwkWAEAAEwkWE20laVNAAC4qy3MmQUrYPO2UOyhZ/oo0APBCgAAYCLB\nCgAAYCLBCgAAYCLBagJ7vgEA4DSjz50FKwAAgIkEKwAAgIkEqzsafSkTtkafhjbpmzCWkfu0YAUA\nADCRYAUAADCRYHUHIy9hwpbp29AWfRLGNGrfFqwAAAAmEqwAAAAmEqwAAAAmEqwADoy67xt6oy8C\nvRGsAAAAJhKsAAAAJhKsbsnWBAAAmGbEObVgBQAAMJFgBQAAMNEiwaqU8rxSyttLKZ8ppXyhlPKR\nUsoLLr3ndaWUR/c//4VSytdc+vmzSilv3f8dny+lvLOU8twljvdUIy5ZAgDAGkabW88erEopz0ny\n/iR/kORbk/zFJP80yWcP3vPDSV6Z5OVJXpTkiSTvLaU88+CvejjJtyf5jiQvSfK8JO+a+3gBLhut\n0ENv9EGgR2Xu4lVKeUOSF9dav/mG9zya5MdqrW/e//n+JI8n+Z5a6zv2f/7dJC+rtf7c/j0PJflo\nkm+stX7wir/zBUkemfWXuUShh+0opax9CLBZxlvYjjOMty+stX546X8kWWYr4N9K8qFSyjtKKY+X\nUj5cSvn+ix+WUr46yYNJ3nfxWq31c0k+kOTF+5e+Lsl9l97zsSSfOnjPWSnyAAAwr5Hm2EsEqz+b\n5B8n+ViS/z3JTyZ5SynlH+5//mCSmt0K1aHH9z9LkgeSPLkPXNe9B2AxIxV66Im+B/TqvgX+zmck\n+WCt9Uf3f/5IKeUvJ/nBJG9f4N8DAABY1RIrVr+T3b1Qhz6a5Cv3//1YkpLdqtShB/Y/u3jPM/f3\nWl33HgAAgCYsEazen+ShS689lOS/JUmt9RPZhaOXXvxwH6C+Icmv7F96JMkfXnrPQ9mFs19d4JgB\nAADubImtgG9O8v5SymuSvCO7wPT9SX7g4D0PJ3ltKeXjST6Z5PVJPp3k3cnuYRallJ9K8qZSymeT\nfD7JW5K8/6onAgLMyRMBYT0X/c+9VkBvZg9WtdYPlVL+XpI3JPnRJJ9I8qpa688evOeNpZRnJ3lb\nkuck+aUk31ZrffLgr3p1ki8leWeSZyV5T5JXzH28AAAAU83+PVZrWfJ7rEY5R8BprFjB+oy9sC0L\njr1df48VQLeEKmiDvgj0RrAC2DORg7bok0BPBCsAAICJBKsj7PGGbXBlHNqkb8I2jDDnFqwAAAAm\nEqyAzXNFHNqmjwI9EKwAAAAmEqxuMMJeT+BmroRDH/RVGF/vc2/BCtgsEzXoiz4LtEywAgAAmOi+\ntQ8A4Fxc7Yb+Xe7HvW8dAsYhWAHDE6hgXBf9W8AC1mYr4DUUaOhfKUWogo3Q32EMPc/BrVgBQzGx\ngm07rAE9T9CA/ghWQLeEKOAmV9UIYQtYimAFdEegAu7KPVnAUgQroAvCFDAnWwaBuQlWQJMEKeBc\nPMIdmINgBTRDmAJaYDULuAvBClidQAW0yj1ZwKkEK2AVwhTQE6tYwDGCFbA4IQoYice4A1cRrIDF\nCFTAVtgyCAhWwOwEKmCrBCzYLsHqCoohnE6IAng62wXh7mqtXc4vBCvg1nosdgBr8wAMGJtgBZxM\noAKYhy2DMB7BCriRMAWwHKtYMA7BCriHIAWwjsv1V9CCvghWlyhibJEwBdAeq1lsWY8PsBCsYMN6\nK1gAW+WeLGifYAUbIEABjOG6ei5wwfoEKxiYQAWwDVa0YH2CFQxIoALYJgEL1iNYwSCEKQAuePAF\nnJ9gBR0SogA41VVjhrAF8xOsoCMCFQBzsGUQ5idYQeOEKQCWYssgzEewOqCg0AJBCoA1XB5/zItY\nW29fEixYQQN6KhoAbIPVLLgdwQpWIkwB0AshC44TrODMBCoAeubBF3A1wQoWJEQBMCqPcYd7CVaw\nAIEKgC2ymsWWCVYwE2EKAHbck8UWCVZwB0IUAJzGlkG2QrCCWxCoAGA6WwYZkWAFJxCoAGB+AhYj\nEazggAAFAOd33fgrcNETwQoiUAFAi6xo0RPBis0SpgCgD54ySA8EKzZHoAKAflnFolWCFUMTogBg\nTB7jTmsEK4YkUAHA9ljNYk2CFcMQpgCAxD1ZrEOwoktCFABwClsGORfBiq4IVADAVLYMsgTBii4I\nVADA3AQs5iRY0SxhCgA4B/dkMQfBimYIUgDA2i7PRwQtTiVYsSphCgBomdUsTiVYcXbCFADQIyGL\nmwhWLE6QAgBGY8sglwlWLEKYAgC2xGoWghWzEqgAgK3zGPdtEqyYTJgCAHg6q1jbIlhxK0IUAMDt\nXTWHErbGIlhxEoEKAGBetgyO5RlrHwBtK6UIVTTLQATchppBq8y3xmDFiiRWpGifCREwh6tqiTGQ\nVlzXFo2BfRCsNs5gQssMJMA5HNYa4yItsmWwD4LVBhk0aJ2BA1jLRf0xVtIiTxlsm2C1IQYJWmeQ\nAFohYNE6q1jtEawGZSCgZQYBoBfX1SvjLK3wGPd2CFYDUeRpnUIPjMKKFi2zZXAdgtUAFHVapqAD\nI/PgC1pny+D5CFadUbRpncINbNXl+mfMpiW2DC5PsOqE4kzLFGaAp7OaReusZs1LsGqYIkzLFGGA\n0wlZtMw9WfMQrBqk4NIyBRdgGg++oGVWse5OsFqZokrLeiiqtVb9CDiqxXp21TGpZ7TCPVm3J1it\nROGkZQonwDqsZtEyq1k3E6zOSJGkdQolQBsELFrmnqyrCVYLUQhpmSII0Ifr6rV5Bq2wZfApgtXM\nFDpattVCBzAaK1q0bKtbBgWrGShqtG5rhQ1gKwQsWra1LYOC1QSKGC3bQgEDYMf3ZNG6LaxiCVYn\nUqRo3ciFCoDTXR4PzGFoycj3ZAlWRyhGtGyUQgTAcqxm0bpRVrMEq2soPLSq96IDwHqELFrWe8AS\nrA4oMLSo1+ICQNtsGaRVvbZFwQoaJEwBcG5Ws2CaZ6x9AMBTaq1C1R04Z8BN1IjbMx7B7VmxggYY\nvABoke/JgtMJVnBGAhQAPbpu/BK44CmCFZyBQAXAiKxowVMEK1iQQAXAFghYIFjBLASo9dVaDejA\n06jP52XLIFsmWMEEBmwAOM6KFlsgWMEtCVMAcDe+K4uRCVZwhCDVD9sBgUPqd9sufz7qN70TrOAK\nBmMAOC+rWfROsII9YQoA2iBk0SPBis0TqACgXR58QS8EKzZFiAKAPl01hgtbtESwYhMEKgAYj9Us\nWiJYMTSBCgDGJ2DRAsGKIQhQXPDIdSAxLmzVdZ+7cYFzEKzomoETADjGihbnIFjRJYEKALgtAYsl\nCVZ0Q5jiVLYDwrYZLzjG92SxBMGKZhkYAYClXZ5vCFrclWBFU4Qp5mLVCrbJOMJUVrO4K8GK1RkE\nAYAWCVnchmDF2QlSnItVK9gW4wtLsmWQYwQrzsJgBwCMxGoWlwlWLEqgAgBG5zHuJIIVCxGoAICt\nEbC2TbBiEgGK1rnPCrbBeERLrmuPxqOxCVbciQEMAOB2rGiNTbDiZMIUAMB0HnwxJsGKawlSjMJ2\nQBib8YqeeYz7OAQr7mFwYlTCFYzJuMVorGb1S7DCoAQA0CAhqy+C1YYJVAAAffDgi/YJVhshRAEA\n9O+qOZ2w1YZnrH0ALKvWKlTBnr4AY9GnYcd8rw1WrAakYwEAbI97stYlWHVOiILb8XRAGIPxD25m\ny+D52QrYKUu+cHf6DvRNH4a7MX9clhWrzugMAABM4QmDyxCsOiBMAQAwN/dkzUuwaowQBefhXivo\nk3ESluGerOkEq0YYKAAAaIktg7cjWK1ImIJ1WbWCvhg3YR22DJ5GsDojAwK0R7iCPhhDoQ2X+6Ix\n9CmC1YIMAgAAjEzQeorvsVqA7wiAvuiv0DZ9FPqx5XmwFasZbbURAQDAoS0++EKwuiMhCsbiXito\nk/EW+ralx7gLVrekwAMAwN2NuprlHqsTbXm/KGyFPg5t0SdhbKPNrwWrE4z0gQM309+hDfoibMco\n/V2wOsFoy5QAANCKUebaiwerUsqPlFL+qJTypkuvv66U8mgp5QullF8opXzNpZ8/q5Ty1lLKZ0op\nny+lvLOU8tyljxcAAOC2Fg1WpZSvT/LyJB+59PoPJ3nl/mcvSvJEkveWUp558LaHk3x7ku9I8pIk\nz0vyriWPFyAZb8839ET/A3q1WLAqpfypJP8myfcn+f8u/fhVSV5fa/33tdb/J8l3Zxec/u7+/3t/\nku9L8upa63+ptf56ku9N8k2llBctdcwAAAB3seSK1VuT/Hyt9T8dvlhK+eokDyZ538VrtdbPJflA\nkhfvX/q67B4Ff/iejyX51MF7ABblqjmclz4H9GyR77EqpbwsyddmF5AuezBJTfL4pdcf3/8sSR5I\n8uQ+cF33HoDF+eJgOA+hCujd7MGqlPJnsrs/6ltqrV+c++8HAABozRJbAV+Y5E8n+XAp5YullC8m\n+eYkryqlPJndqlPJblXq0ANJHtv/92NJnrm/1+q695yVK9awXa6kw7L0MdiukebYSwSrX0zyV7Lb\nCvj8/f8+lN2DLJ5fa/3t7MLRSy/+D/sA9Q1JfmX/0iNJ/vDSex5K8pVJfnWBYwa4kYkfLEPfAkYx\n+1bAWusTSX7z8LVSyhNJfq/W+tH9Sw8neW0p5eNJPpnk9Uk+neTd+7/jc6WUn0ryplLKZ5N8Pslb\nkry/1vrBuY8Z4BTut4J5CVXASBZ5eMUV7qmctdY3llKeneRtSZ6T5JeSfFut9cmDt706yZeSvDPJ\ns5K8J8krznO4VyulGAQAAGAGo12sLKMEhVLKC7LbQrioUc4XcHejDQSwBuMpcKbx9IW11g+f4x9a\n8nusAIZkQgjT6EPAiASrW3KlGkh2E0OTQ7gd/Qa4MOKcWrACAACYSLACmMDVdziNvgKMTrACmMiE\nEW6mjwBbIFgBzMDEEa6mbwBbIVgBAABMJFgBzMSVebiXPgFsiWB1ByM+HhKYh8dJg34A3GzUubRg\nBbAAk0q2StsHtkqwAgAAmEiwuqNRlzCB+dgOxZZo78ApRp5DC1YACzPZZHTaOIBgBQAAMJlgNcHI\nS5nAvFzRZ1TaNnCq0efO9619AABbcTEBHX1gYRsEKoB7WbECODMTUnqnDQM8nWAFsAITU3ql7QJc\nzVZAgJWo+0m2AAAaiElEQVTYGkhPBCqAm1mxmsiECAAAbraFObNgBbAyKwG0ThsFOE6wAmiAiSut\n0jYBTiNYzWALS5vA8mqtJrE0Q3sE5rKVubJgBdAYk1nWpg0C3J5gBdAgE1vWou0B3I1gNZOtLHEC\n52MrFuekvQFL2NIcWbACaJzJLkvTxgCm8wXBAB3wZcIsQaACmI8VqxmZ8AAAwM7W5sZWrAA6crjC\nsLUBi3lYpQJYhhUrgE6ZIHNb2gzAcqxYAXTMvVecQqACWJ4Vq5mZ3ABr8KhsrqJdAGvZ4pxYsAIY\nxBYHMW6mTQCcj2AFMAATaK6jbQCch2C1AIMYcE5qDsdoI8A5bbXmeHgFQKe2OnBxNxftxT1XAMuw\nYrUQEx4AALZmy3NgwQqgQ1seuJhG2wFYhq2AAB0xKWYOtgUCzM+KFQAAwESCFUAnrFYxN20KYD62\nAi6olGKbBTCJiS9LO2xjxixgiq2PWVasAAAAJhKsFrb15A7cnfrBuWlzwF2pH4IVQJMMUKxF2wO4\nG8HqDAxSwG2oGaxNGwRuQ83YEawAGmJwohXaIsDtCFZnYoACjlEnaI02CRyjTjzF49YBVmZQOs3F\no8DnOl9z/32jujg/HsUOcDPBCmBFJvWnOZzUzz3Br7X6HE7guxkBbmYrIMBKTOZPc47JvMBwGm0W\n4HqC1RkZkIAL6sFpzhl4hKvTaLvABfXgXoLVmWmAgDpwXK11laCz1r/bG20YUAeeTrACOCMDEaPQ\nlgHuJVgBnImJ6GlaWDFq4Rh6oE0DPEWwWoGBCLZHvz+utW14rR1Pq7Rt2B79/mqCFcDCDECMThsH\nEKwAFmXCeVzrK0OtH18rtHVg6wQrgIWYaB7XU2Dp6VjXos0DWyZYrcTgA+MqpejjJ+gxqPR4zOem\n/cPY9O/rCVYr0jBhPPr1aXoOKD0f+znpCzAe/fpmghXATAw4pxkhmIzwO5yDPgFsiWC1MoMOjEFf\nPs1IgWSk32VJ+gaMQV8+TrACmMhgc5oRg8iIv9MS9BFgCwSrBhhwoF/673GjP6589N9vLvoK9Ev/\nPY1gBQAAMNF9ax8AQG9cuTvN1lZxLn5f7eN6h+dma+0DGJ8Vq0YYiKEP+upptjxp3vLvfhv6EvRB\nXz2dYAXArAQL5wBgiwSrhrgiAO0qpeijJxAonuJcHKdfQdv0z9sRrACOMLCcRpB4OufkNPoYMALB\nqjEGF6BHAsT1nBugR+aktydYAdzAwALnoa8BvfO4dYArmOQdZyXmdIfnStu63sW50baAHlmxapBB\nF9alDx5n4nt3zt1x+iCsSx+8G8EK4IDBBNqgLwK9sRWwUaUUVzXhjEziTqMuzePiPGp3N7M1EM5P\nXbo7K1YN07CBlpjczs85BVpi7jmNYAVsnoEE2qaPAj0QrBpnMIFl6WOnsbKyHOf2NPoqLEsfm849\nVsAmGUBOY9J/Hu65Oo17roCWWbEC4Eomr+fnnAP0S7DqgCuYMC996jgT/PU498fpwzAvfWoetgJ2\nwuPXYToDx3HqTBtsDTzOtkCYhzozHytWwCYYOGBM+jbQCsGqIwYPuBt95zSu/LfHZ3IafRzuRt+Z\nl2AFDM2gcVyt1QS+YT6f0+jrwNoEK2BYJlrHmbD3w2d1nD4PrEmw6oxBA06jr8A26ftwGn1lfoIV\nMByDxXG2l/XJ53YaNQBYg2DVIQMGXK2Uon+cwMS8fz7D49QDuJ6+sQzBqlM6BNxLnziNCfk4fJan\nURvgXvrEcgQroHsGidOYiI/HZ3oaNQI4B8GqYwYK0A9O4b6csfl8T6NWgH6wNMEK6JYBArgNNQNY\nkmDVOYMEW6Xtn8ZKxnb4rE+jdrBV2v7yBCugOwaH05hob4/P/DRqCLAEwWoABgi2RHs/zj032+bz\nP41awpZo7+chWA1Ch2F0vpPmNCbUXNAWjlNX2AJt/HwEK4BBmEhzmTYBcD6C1UBckWBEriifxgSa\n62gbx6kzjEq7Pi/BCmiWAeE0Js4co42cRs0BphCsBmNQYBTa8mlMmDmVtnIatYdRaMvnJ1gBdMpE\nmdvSZgCWI1gNyBUKeuZeh9OYIHNX2s5x6hC9037XIVgNSoeiR9rtaUyMmUobOo2aRI+02/UIVgAd\nMSFmLtoSwLzuW/sAWE4pxcBJ81xZO42+zBIu2pV+eL3Dc6Mf0jp9eV1WrIDVGABOYzLH0rSx06hZ\nwE0Eq8EZBKBvJryci7YGfTPnW59gtQE6Gi3SLoEeqV20SLtsg3usgLNS/E9j9YA1uOfqNBfnRz8F\nDlmx2giDJPTDZI21aYPQD3O8dghWwNko/sBI1DTgkK2AG+Lx66zF5OM4fZPWHLZJffh6tgWyJn2z\nLVasNkYH5Ny0ueNMyGidNnqcWse5aXPtEayAxSj6wJaoebBtghWwCBOM01gJoBfa6mnUPtgu91ht\nkP3gLM3E4jj9jx55HPtp3NPMkvS/dlmxAmal4AOohbBFgtWGKfrMTZs6jSvZ9E4bPo2ayNy0qbbZ\nCrhxtiswB4X+OP2M0Xgc+2lsv2cu+ln7rFgBAABMJFjhCgh3VkrRfk7gSjWj08aPUy+ZQtvpg2AF\n3IkifxoTTrZCWz+N2gnjEqxI4koat6OtnMZEk63R5k+jhnIq87O+CFbArSjwpzHBZKu0/dOopTAe\nwQo4mYnAaUws2Tp94DRqKozF49a5h8fCch0TgOP0G3jKRX9QO27ma0+4in7TJytWwFEK/HEmRnA1\nfeM4NRbGIFhxJUWeC9rCcSaOcDN95Di1lgvaQr8EK66lY6MNHGfCCKfRV45Tc9EG+iZYAU/j8a6n\nMVGE29FnjlN/oV+CFTdS3LfHZ34aE0S4G33nNGrx9vjM+ydYcZSODvcyMYRp9CG4l7nWGAQrTqLD\nj8/2k9OYEMI89KXj1OVt8BmPQ7ACFPUTmQjCvPSp06jR0AfBCjbOgH0aE0BYhr51GrUa2nff2gdA\nPy6KukFwHAbq47R3WN5FP1OTblZKUZMGor2Px4oVt6YQjMHneJwJDJyXPnec2j0Gn+OYBCvYIAX9\nOBM8WIe+d5waDm0SrLgTRZ2RmdjBuvRBRmYONS73WHFn9nr3RSE/jTYNbai1qls3ODw36lY/tOmx\nzb5iVUp5TSnlg6WUz5VSHi+l/Fwp5c9f8b7XlVIeLaV8oZTyC6WUr7n082eVUt5aSvlMKeXzpZR3\nllKeO/fxwhYo5KcxOYG26JOnUeOhDUtsBfxrSf5lkm9I8i1JvizJfyyl/MmLN5RSfjjJK5O8PMmL\nkjyR5L2llGce/D0PJ/n2JN+R5CVJnpfkXQscLxMo5ozCBA7apG8yCnOm8ZWlC1Yp5SuS/PckL6m1\n/vL+tUeT/Fit9c37P9+f5PEk31Nrfcf+z7+b5GW11p/bv+ehJB9N8o211g9e8e+8IMkji/4yXMvA\n1yZF/DhtF/qhph2nprVJ213VC2utHz7HP3SOh1c8J0lN8j+SpJTy1UkeTPK+izfUWj+X5ANJXrx/\n6euyu//r8D0fS/Kpg/fQEAWjPT6T40xAoC/67HFqf3t8JtuxaLAqu5b0cJJfrrX+5v7lB7MLWo9f\nevvj+58lyQNJntwHruveA3BnJmjQJ30XaNXSTwX8iSR/Kck3Lfzv0ABPCWyHq2MA22U8bofxeFsW\nW7EqpfyrJH8zyV+vtf7OwY8eS1KyW5U69MD+Zxfveeb+Xqvr3kODSimKyEouzr3zf7NaqwkHdE4/\nPs6YsC7nfpsWCVb7UPV3kvyNWuunDn9Wa/1EduHopQfvvz+7pwj+yv6lR5L84aX3PJTkK5P86hLH\nDIzPRAzGok8DLZl9K2Ap5SeS/P0kfzvJE6WUi5Wp/1lr/f39fz+c5LWllI8n+WSS1yf5dJJ3J7uH\nWZRSfirJm0opn03y+SRvSfL+q54ICFvnqhgA17E1EM5jiXusfjC7h1P850uvf2+Sn0mSWusbSynP\nTvK27J4a+EtJvq3W+uTB+1+d5EtJ3pnkWUnek+QVCxwvC7iY6CvkyxKoTqMdwrgu+rd6eDPj8nlo\nh9u2+PdYnYvvsWrXKG2sRQr4cdofbIeaeJyauBztr1lDfY8VsAAF/DgTCNgWff44YwcsR7CCDhkY\njzPBgm3S948zhsAylv4eK7Cve0YGw+O0M8B9V8cZm+ejnXHBihVno/BM4/wBMDdjyzTOH4cEK85K\nAbob5+00rrwCh9SE0xhj7sZ54zLBirNTiG7H+Tqu1moCBVxJfTiNseZ2nC+uIlixCgXpNM4TAOdi\nzDmN88R1BCtolMJ9GleigVOoFacx9sDdCVasRvG+nnNznO09wG2pG6cxBl3PueEmghWrUqCezjk5\nzsQImEINOc5Y9HTOCccIVqxOoXqKc3GcCREwB7XkOGPSU5wLTiFY0QQFyzk4hYkQMCc15Thjk3PA\n6QQrmrHVwlVK2ezvfhsmQMAS1JbjtjxObfX35m4EK1iRgn0aEx9gSWrMaYxZcDPBiqZsqWhv6Xed\nwoQHOAe15jRbGru29LsyD8GK5myhkG3hd5yDiQ5wTmrOabYwhm3hd2R+ghVNGrmgjfy7zcV3zQBr\nUX9OM/JYNvLvxrIEK5o1YmEb8XeamwkN0AK16LgRx7QRfyfOR7CiaSMVuJF+l6WYyAAtUZOOG2ls\nG+l3YR2CFc0b4TGvvR//OZjAAC1Sm47rfYwbYZ5BGwQrWJBifRoTF6BlatRxxjsQrOiIgj0mExag\nB2rVmMwtmJNgRVd6KoA9HetaTFSAnqhZx/U09vV0rPThvrUPAG6rlNL04KZQH9fy5wdwk4v6pdZf\n7+LctFzrfX4swYoVXWq1ILZ6XC1peaAFOJVadlyrY2Krx0X/BCu6pTACALdh7sCSbAWka61sC1So\nj2vhcwKYk22Bx7W0LdDnxNKsWNG9tQvl2v9+D1oYUAGWosYdt/ZYufa/zzYIVgxBwWyXCQewBWpd\nu8wROBdbARnGObcbKNKnMdEAtqTWany4weG5MVYzIitWcEsK9WmEKmCL1L7TGEsZkWDFcBTr9ZlY\nAFumBq7PXIA1CFYMaamCqlAfZ0IBoBaewljNaNxjxbDmfBS7In2cSQTAvTyO/bi57492rlmTFSuG\npsCeh1AFcD018jyM+axNsGJ4UwutQg0AyzJWMwJbAdmEu2w1UKSPcxUW4DS2BR5nrKZ3VqzYlFML\nsEJ9nFAFcHtq53HGanolWLE5xwqxQg0A6zJW0yPBik26riAr1KdxxRXg7tTQ0xir6Y1gxWZdLswK\n9XG1VhMCgBmop6cxVtMTD69g0xTo05kAAMyv1mosOsL5oRdWrICjhCqA5aixMAYrVsC1DPYA5+Fx\n7NA/K1YAAAATCVbAlaxWAZyf2gv9EqyApzGwA6xHDYY+uccK+GMGc4A2uOcK+mPFCkgiVAG0SG2G\nfghWgIEboGFqNPRBsIKNM2ADtE+thvYJVrBhBmqAfqjZ0DbBCjbKAA3QH7Ub2iVYwQYZmAH6pYZD\nmwQr2BgDMkD/1HJoj2AFG2IgBhiHmg5t8QXBsAEGX4Ax+SJhaIcVKxicUAUwPrUe1idYwcAMtADb\noebDugQrGJQBFmB71H5Yj2AFAAAwkYdXwGBcrQTYNg+0gHVYsYKBCFUAXDAmwHkJVjAIAygAlxkb\n4HxsBYTOGTQBuImtgXAeVqygY0IVAKcyZsCyBCsAAICJBCvolCuPANyWsQOW4x4r6IxBEYAp3HMF\ny7BiBR0RqgCYizEF5iVYAQAATGQrIDTOFUUAlnI4xtgaCNNYsYLGGegAWJqxBqazYgUduBjwrF4B\nMCeBCuZjxQoAAGAiK1bQkcMri1avALgLq1SwDCtW0CkDIwC3ZeyA5Vixgo659wqAUwhUsDzBCgYg\nYAFwFYEKzsdWQBiIARSAC8YEOC8rVjAYq1cA2yZQwToEKxiUgAWwLQIVrMtWQBicgRZgfGo9rM+K\nFWyA1SuAMQlU0A7BCjZEwAIYg0AF7RGsYIMELIA+CVTQLvdYwYYZoAH6oWZD26xYwcZZvQJom0AF\nfRCsgCQCFkBrBCroi2AF3EPAAliXQAV9co8VcCUDO8D5qb3QLytWwLWsXgGch0AF/ROsgKMOB3wh\nC2AewhSMRbACbsUqFsA0AhWMyT1WwJ2YGADcntoJ47JiBdyZLYIAxwlTsA2CFTALWwQB7iVQwbbY\nCgjMykQCQC2ELbJiBczO6hWwVQIVbJdgBSzGPVjAFghTQCJYAWdiFQsYjUAFHBKsgLMSsIDeCVTA\nVQQrYBW2CQI9EaaAYwQrYHVWsYBWCVTAqQQroBkCFtAKgQq4Ld9jBTTHhAZYkxoE3IUVK6BJlyc2\nVrGApQhSwBwEK6ALHnYBzEmYAuZmKyDQHRMiYAo1BFiCFSugS7YKAqcSpIBzEKyAIXiiIHCZQAWc\nk2AFDEXAAgQqYA2CFTAkAQu2R6AC1iRYAUO7aqIlbEH/hCigNYIVsDlWs6BfAhXQKsEK2CzfjQV9\nEKaAHghWAPH4dmiJIAX0SLACuILVLDgvYQronWAFcISQBcsQpoCRCFYAtyBkwTTCFDAqwQrgjjzK\nHW4mRAFbIlgBzMiKFlsnTAFbJVgBLETIYiuEKQDBCuAsbpp4Cl30QHgCuJlgBbAyK1u0SpgCOJ1g\nBdAQIYu1CVMAdyNYATTK9kGWIjwBzE+wAujQ5YmxoMVNBCmA5QlWAAMQtDgkSAGcn2AFMCBfXrwd\nQhRAGwQrgI24zQRcCFuXsATQH8EKgKe5bmIvcM1LgAIYh2AFwMlOCQLC147QBLAtghUAsxIoANii\nZ6x9AAAAAL0TrAAAACYSrAAAACYSrAAAACYSrAAAACYSrAAAACYSrAAAACYSrAAAACYSrAAAACYS\nrAAAACZqPliVUl5RSvlEKeV/lVJ+rZTy9WsfEwAAwKGmg1Up5TuT/HiSf57kryb5SJL3llK+YtUD\nAwAAONB0sEry6iRvq7X+TK31t5L8YJIvJPm+dQ8LAADgKc0Gq1LKlyV5YZL3XbxWa61JfjHJi9c6\nLgAAgMuaDVZJviLJn0jy+KXXH0/y4BXv//LFjwgAAOjJ2TJCy8Hqtr5q7QMAAACa8lXn+ofuO9c/\ndAefSfKlJA9cev2BJI9d8f73JvkHST6Z5PcXPTIAAKBlX55dqHrvuf7BsrttqU2llF9L8oFa66v2\nfy5JPpXkLbXWH1v14AAAAPZaXrFKkjcl+elSyiNJPpjdUwKfneSn1zwoAACAQ00Hq1rrO/bfWfW6\n7LYA/tck31pr/d11jwwAAOApTW8FBAAA6MFITwUEAABYhWAFAAAw0TDBqpTyilLKJ0op/6uU8mul\nlK9f+5joRynlNaWUD5ZSPldKebyU8nOllD9/xfteV0p5tJTyhVLKL5RSvubSz59VSnlrKeUzpZTP\nl1LeWUp57vl+E3pTSvmRUsoflVLedOl1bY3JSinPK6W8fd9OvlBK+Ugp5QWX3qOtMUkp5RmllNeX\nUn57344+Xkp57RXv09YY2hDBqpTynUl+PMk/T/JXk3wkyXv3D76AU/y1JP8yyTck+ZYkX5bkP5ZS\n/uTFG0opP5zklUlenuRFSZ7Irp098+DveTjJtyf5jiQvSfK8JO86xy9Af/YXgF6eXc06fF1bY7JS\nynOSvD/JHyT51iR/Mck/TfLZg/doa8zhR5L8oyT/JMlfSPJDSX6olPLKizdoa2zBEA+vuOb7rv7f\n7L7v6o2rHhxd2ofy/57kJbXWX96/9miSH6u1vnn/5/uTPJ7ke/ZPsLw/ye8meVmt9ef273koyUeT\nfGOt9YMr/Co0qpTyp5I8kuQfJ/nRJL9ea/0/9j/T1pislPKGJC+utX7zDe/R1pislPLzSR6rtf7A\nwWvvTPKFWut37/+srTG87lesSilfluSFSd538VrdpcVfTPLitY6L7j0nSU3yP5KklPLVSR7Mve3s\nc0k+kKfa2ddl9xUGh+/5WHZfaq0tctlbk/x8rfU/Hb6orTGjv5XkQ6WUd+y3OH+4lPL9Fz/U1pjR\nryR5aSnlzyVJKeX5Sb4pyX/Y/1lbYxOa/h6rE31Fkj+R3VWPQ48neej8h0Pv9iueDyf55Vrrb+5f\nfjC7oHVVO3tw/98PJHlyP1hc9x5IKeVlSb42u4nEZdoac/mz2a2I/niS/zO77VdvKaX8Qa317dHW\nmM8bktyf5LdKKV/K7sL9P6u1/uz+59oamzBCsIK5/USSv5Td1TaYVSnlz2QX3L+l1vrFtY+HoT0j\nyQdrrT+6//NHSil/OckPJnn7eofFgL4zyXcleVmS38zuwtG/KKU8ug/xsAndbwVM8pkkX8ruSseh\nB5I8dv7DoWellH+V5G8m+eu11t85+NFjSUpubmePJXnmfp/4de+BFyb500k+XEr5Yinli0m+Ocmr\nSilPZnd1VltjDr+T3f0phz6a5Cv3/62uMZc3JnlDrfXf1Vp/o9b6b5O8Oclr9j/X1tiE7oPV/orv\nI0leevHafivXS7Pb8wsn2Yeqv5Pkb9RaP3X4s1rrJ7Ir7Ift7P7sniJ40c4eSfKHl97zUHaTmF9d\n9ODpyS8m+SvZXdF9/v5/H0ryb5I8v9b629HWmMf78/Qt8Q8l+W+Jusasnp3dRe5Df5T9PFNbYytG\n2Qr4piQ/XUp5JMkHk7w6u07+02seFP0opfxEkr+f5G8neaKUcnFV7X/WWn9//98PJ3ltKeXjST6Z\n5PVJPp3k3cnuRtxSyk8leVMp5bNJPp/kLUne72lGXKi1PpHdVpk/Vkp5Isnv1VovVhe0Nebw5iTv\nL6W8Jsk7spvEfn+SHzh4j7bGHH4+u3b06SS/keQF2c3F/vXBe7Q1hjdEsNo/pvMrkrwuuyXj/5rk\nW2utv7vukdGRH8zuxtr/fOn1703yM0lSa31jKeXZSd6W3VMDfynJt9Vanzx4/6uzu2r3ziTPSvKe\nJK9Y9MgZwT3fe6GtMYda64dKKX8vuwcL/GiSTyR51cEDBbQ15vLK7ILSW5M8N8mjSX5y/1oSbY1t\nGOJ7rAAAANbU/T1WAAAAaxOsAAAAJhKsAAAAJhKsAAAAJhKsAAAAJhKsAAAAJhKsAAAAJhKsAAAA\nJhKsAAAAJhKsAAAAJhKsAAAAJvr/AdD/fZuf22VZAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x19f41f3d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,10))\n",
    "plt.imshow(padded_wfirst_mask)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 468,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "wfirst_LS = padded_wfirst_mask*annular_stop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 469,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "if flip:\n",
    "    Lyot_stop = np.copy(wfirst_LS[::-1,:])\n",
    "else:\n",
    "    Lyot_stop = np.copy(wfirst_LS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 470,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x1a3172650>"
      ]
     },
     "execution_count": 470,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA/EAAAPRCAYAAAC1bSijAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3W2obdt5F/D/uBwTDRICam4IEqxUr4pSSWrSUEyVFKUW\n3+gHo6JFaWs1gRIE22L8knywpJBeo1EK9otVkXhLKRVprBWlTUtDEynYxkBoYwjXXFuNJiTobdPp\nh713ss8++2W9zDHnePn9IJCz97rnrL3WGM/z/NeYa+2yLEsAAACA9j219x0AAAAADiPEAwAAQCeE\neAAAAOiEEA8AAACdEOIBAACgE0I8AAAAdEKIBwAAgE4I8QAAANAJIR4AAAA6IcQDAABAJ4R4AAAA\n6IQQDwAAAJ0Q4gEAAKATQjwAAAB0QogHAACATjza+w6spZTy25L8ySSfSPJ/9703AAAATOA3J/ld\nST6wLMv/3OIfHCbE5yLA/4u97wQAAADT+ctJ/uUW/9BIl9N/Yu87AAAAwJQ+sdU/NFKIdwk9AAAA\ne9gsj44U4gEAAGBoQjwAAAB0QogHAACATgjxAAAA0AkhHgAAADohxAMAAEAnhHgAAADohBAPAAAA\nnRDiAQAAoBNCPAAAAHRCiAcAAIBOCPEAAADQCSEeAAAAOiHEAwAAQCeEeAAAAOiEEA8AAACdEOIB\nAACgE0I8AAAAdEKIBwAAgE4I8QAAANAJIR4AAAA6IcQDAABAJ4R4AAAA6IQQDwAAAJ0Q4gEAAKAT\nQjwAAAB0QogHAACATgjxAAAA0AkhHgAAADohxAMAAEAnhHgAAADohBAPAAAAnRDiAQAAoBNCPAAA\nAHRCiAcAAIBOCPEAAADQCSEeAAAAOiHEAwAAQCeEeAAAAOiEEA8AAACdEOIBAACgE0I8AAAAdEKI\nBwAAgE482vsOAMDMlmXZ+y40p5Sy910AgGYJ8QBQiYB+mkMeN0EfgFkJ8QBwAgF9X4I+ALMS4gHg\nHsJ6v+567oR7AHomxAMwPUF9Lvc93wI+AK0T4gGYisDOfW5bH4I9AC0R4gEYkrDOWlyWD0BLhHgA\nhiC0s7Wba06oB2ALQjwA3RHYaZFL8QHYghAPQNMEdnom2AOwNiEegKYI7YzOZfgAnEOIB2B3gjsz\nu77+BXoAHiLEA7A5oR1u55QegIcI8QBUJ7TDaYR6AG4S4gFYhaAO9d23zwR8gDkI8QCcTHCHdnhv\nPcAchHgADia0Qx9chg8wLiEegAcJ79C3qz0szAP0T4gH4FaCO4zHJfcA/RPiAfgSwR3mIdAD9EmI\nB5iY0A4k3kMP0BMhHmAygjvwEKf0AO16au87AMB2BHjgWOoGQFucxAMMzgAOnMvJPEA7hHiAAQnu\nQC0CPcC+hHiAQQjuwNYEeoDtCfEAHRPcgVYI9ADbEOIBOiK0Az3wK+sA6hHiATogvAM9u6phwjzA\n+YR4gEYJ7sBoXHIPcD4hHqAxwjswA6fzAKcR4gEaILgDs3I6D3Ccp/a+AwCzE+ABLqiHAA9zEg+w\nA4MqwO2czAPcT4gH2JDwDnA475sHeJIQD7AB4R3gdMI8wJcJ8QAVCe8A6xHmAYR4gNUJ7gB1ed88\nMDOfTg+wIgEeYFvqLjAbJ/EAKzBEAuzHZfbATIR4gDMI7wDtEOaBGQjxAEcQ2gHad7NWC/XASIR4\ngAMI7wD9ckIPjESIB7iH8A4wDmEeGIFPpwe4gwAPMCb1HeiZk3iAGwx3AONzKg/0SogHuCS8A8xH\nmAd6I8QD0xPeARDmgV4I8cC0hHcAbhLmgdYJ8cB0hHcAHiLMA60S4oFpCO8AHEuYB1rjV8wBUxDg\nATiHPgK0QogHhmfwAmAN+gnQApfTA8MybAGwNpfXA3sT4oFhCO0AbOVmzxHqga24nB4YggAPwJ70\nIWArTuKBrhmaAGiFS+2BLQjxQJeEdwBaJcwDNbmcHuiOAA9AD/QroAYhHuiKgQiAnuhbwNpcTg90\nwRAEQK9cXg+syUk80DwBHoAR6GfAGpzEA80y7AAwGqfywLmEeKA5wjsAoxPmgVO5nB5oigAPwEz0\nPeBYQjzQDIMMADPS/4BjCPFAM1xSCMCM9D/gGEI80BSDDAAz0feAY/lgO6A5VwONywsBGJXwDpzK\nSTzQLAMOACPS34BzCPFA0ww6AIxEXwPOJcQDzTPwADAC/QxYgxAPdMHgA0DP9DFgLUI80A0DEAA9\n0r+ANQnxAAAA0Am/Yg7oil8/B0AvnMADNTiJB7pkMAKgZfoUUIsQD3TLgARAi/QnoCYhHuiaQQmA\nluhLQG1CPNA9AxMALdCPgC0I8cAQDE4A7EkfArbi0+mBYfjkegC2JrwDW3MSDwzHQAXAFvQbYA9C\nPDAkgxUANekzwF5cTg8My+X1AKxNeAf25iQeAAAAOiHEA8NzagLAGvQToAVCPDAFgxcA59BHgFYI\n8cA0DGAAnEL/AFoixANTMYgBcAx9A2iNEA9Mx0AGwCH0C6BFQjwwJYMZAPfRJ4BWCfHAtAxoANxG\nfwBaJsQDUzOoAXCdvgC07tHedwBgb1cD27IsO98TAPYivAO9EOIBLgnzAPMR3oHeuJwe4AYDHcAc\n1HugR0I8wC0MdgBjU+eBXgnxAAAA0AkhHuAOTmkAxqS+Az0T4gHuYdADGIu6DvROiIdB+ET1egx8\nAGNQz+sxh8B2hHgYwFXj1EDrMfgB9E0dr8ccAtsS4qFjy7I80TA10HoMgAB9Ur/ruW0OMYtAXUI8\ndOq+Bql51lNKMQwCdELNrsssAvsQ4qFDGuP+DIUAbVOn92degTqEeOjMoQ1R46zPgAjQJvW5PvMI\n7EeIh44c2wg1zvoMigBtUZfrM4/Avh7tfQeAh53T/K7/twabOq4eV0MKwH70uLrO7XFX/73nCc7n\nJB4at2YwFDLrMpgA7EP9rcssAm0R4mEymicAcChzA7RHiIeG1WqcGnI9ToMAtqXu1mMOgTYJ8dCg\nZVmqNzgNtB4DJcA21Nt6tphDzCJwGiEeGrNlQ9M86ymlGC4BKlFj6zKLQNuEeGjIHo1M86zLkAmw\nLnW1LrMItE+Ih0bs2cA0z7oMnADrUE/rMotAH4R4aEALjauF+zAygyfAedTRulqYA1q4D9ADIR52\npmHNwwAKcBr1cx7mIniYEA87aq1R+aTY+gyiAMdRN+tqsfe3dn+gNUI87KTlBtXyfRuBgRTgMOpl\nXS33+5bvG+xNiIcd9NCYeriPPTOYAtxPnayrhz7fw32EPQjxsLGeGlJP97VHBlSA26mPdfXU33u6\nr7CVR3vfAZhJj43o6j4bqOq4elx7XBsAa9Nr6uq11yzLYm3ANU7iYSO9Ns4rvd//1hlOgNmpg3X1\n3sd7v/+wJiEeNjBK4xnl52iVARaYlfpX1yj9e5SfA84lxENlozWc0X6e1hhkgdmoe3WN1rdH+3ng\nFEI8VDRqoxn152qFgRaYhXpX16j9etSfCw4lxEMlozeY0X++vRlsgdGpc3WN3qdH//ngPkI8VDBL\nY1mWZZqfdQ8GXGBU6ls9M/XmWX5OuEmIh5VpKKzJoAuMRl1jTeYuZiTEw4pmbSSz/txbMfACo1DP\n6pq1H8/6czMvIR5WMnsDmf3nr83gC/ROHatr9j48+8/PXIR4WIHGccHjUJcBGOiV+lWX/nvB48As\nhHhgVRpoXaUUwzDQDTWrPn0X5rN6iC+lPFVKeVcp5ZdKKV8opXy8lPKOW273zlLK85e3+fFSylfe\n+P5LSynvK6X8ainlc6WU50opr1z7/sK5NM8neUzqMxQDrVOn6tNvn+QxYQY1TuK/K8nfSPK3kvy+\nJH8nyd8ppbzt6gallO9M8rYk35bk9Uk+n+QDpZSXXPt7nk3yjUm+Kcmbkrw6yQ9VuL9wkpl+hcsp\nPD71GZCBVqlPdemx9/P4MLqy9gIvpfxokk8vy/Kt1772XJIvLMvyVy///HyS712W5fsu//zyJC8k\n+eZlWd5/+edfSfKWZVl++PI2zyT5aJKvWZblQ7f8u69N8uFVfxi4g8ZwHMNcXdYj0BI1vy41/zjW\nIxt63bIsH9niH6pxEv/TSd5cSvk9SVJK+aokX5vk317++SuSvCrJT1z9B8uyfDbJzyZ54+WXvjrJ\noxu3+ViST167DdAJr4jXZUABWqEe1aOXAlceVfg7vyfJy5P811LKF3PxQsHfXZblX11+/1VJllyc\nvF/3wuX3kuTpJC9ehvu7bgO70EBPtyyLAa+SUoq1CexKfa9HfT+d2YMR1QjxfyHJX0ryliS/mOQP\nJ/kHpZTnl2X5wQr/HmxGEz2fZlqPIA/sRV2vR10/n9mD0dQI8e9O8veXZfnXl3/+hVLK70ry3Ul+\nMMmnk5RcnLZfP41/Osl/vvz/n07yklLKy2+cxj99+T3YnCa6Hs20HkEe2Jp6Xo96vh6zByOp8Z74\nlyX54o2v/cbVv7Usyy/nIoi/+eqblx9k94ZcvJ8+ufiAul+/cZtnkrwmyc9UuM9wL02UnhhSgK2o\nN/TEPMcoapzE/2iSd5RSPpXkF5K8Nsnbk/zTa7d59vI2H0/yiSTvSvKpJD+SXHzQXSnlB5K8p5Ty\nmSSfS/LeJB+87ZPpoSYFvw6viNd19dhav0AN6nd96ncd5g9GUCPEvy0Xofx9SV6Z5Pkk/+Tya0mS\nZVneXUp5WZLvT/KKJD+Z5BuWZXnx2t/z9lyc6D+X5KVJfizJWyvcX2AnGml9Lq8H1qZu16duA/dZ\n/ffE78XviaeGUfZH6wyE9VnLwBrU6/rU621Yy1TQ9e+JhyFootvxWNdnWAHOpY7Upx9ux2NNz4R4\nuIXCvj2PeX0GcOBU6kd9+uD2POb0SoiHGxT0/SzL4vGvzCAOHEvdqEvv25fHnh4J8XCNQt4GzwMA\nM9Dv2uB5oDc1Pp0euqSAt+Xq+XACtD6/fg44hPpbj/rbHr8xh544iYdopi3z3NRjWAHuoj7Uo6+1\ny3NDL4R4pqdgt89zVE8pxbAOfImaUJd+1j7PET0Q4oEuaKp1GdoBdaAufQxYixDP1DTUvni+6jLA\nw7zs/7r0r754vmidEA/AlxjkYT72PUBfhHim5VXWPnne6jPQwzzs9/r0rT553miZEM+UFOa+ef7q\nM9jD+Ozz+vSrvnn+aJUQz3QU5DF4Husz4MO47O/69KkxeB5pkRDPVBTisXg+6zPow3js6/r0p7F4\nPmmNEA90bVkWzbUyAz+Mw36uS08CtiDEMw1NdWyeXwD2pA+NzfNLS4R4pqDwzsHzXI/TO+iffVyP\n/jMHzzOteLT3HYCaFNv5XD3nhtX1XX9M7S3og1pYl1o4H3MGLXASDwzJYFWX4QXaZ5/Wpc8AexHi\nGZbmijVQl4AA7bI/69JfsAbYkxDPkBRWrlgLdQkK0B77si59hSvWAnsR4hmOggrbEhigHfYjbMvc\nyR6EeGB4Gmx9ggPszz6sTz8BWiDEMxTNlbtYG/UJELAf+68+fYS7WBtsTYhnGAooD1mWxTqpTJCA\n7dl3dekdHMIaYUtCPDAdjbYugQK2Y7/VpV8ALRLigSkZzOoqpQgXUJE9Vp8+AbRKiGcIGi2nsG4A\nuI3+wCmsG7YixNM9BZNzWD91OSmE9dlXdekLnMP6YQtCPF1TKFmDdVSXwAHrsZ/q0g9Yg3VEbUI8\nQDTc2gQPOJ99VJc+APTi0d53AE6h0VLD1boyKNdx9bjav3AcNakuNYkazBTU5CQe4AYDXV0GGjic\n/VKXeg/0SIinOxouW7DO6hJM4GH2SV3qPFuwzqhBiKcrCiFbst7qElDgbvZHXeo7W7LeWJsQD3AP\njbeuUoqwAtfYE/Wp60DvhHi6oemyF2uvPqEF7IMtqOfsxdpjTUI8XVD42Js1WJ8Aw8ys//rUcfZm\nDbIWIR7gQJpvfYIMM7Lu61O/gZEI8TRP46Ul1iNAX9RtWmI9sgYhnqYpdLTIuqzLqSQzsd7rUq9p\nkXXJuYR4gBNowHUJNszAOq9LnQZGJcTTLM2X1lmjdQk4jMz6rkt9pnXWKOcQ4mmSwkYvrNW6BB1G\nZF3XpS7TC2uVUz3a+w4A9O56Ezacr+/6Y2rgoVdqQ11qAzATJ/E0RyMG7iII0SPrFriLuZdTCPEA\nK9KM6xOI6In1Wp+6C8xGiKcpGjEjsI7rE4zogXVan3rLCKxjjiXE0wwFjJFYz/UJSLTM+qxPnWUk\n1jPHEOIBKtGQ6xOUaJF1WZ/6CsxMiKcJmjGjsrbrE5hoifVYn7rKqKxtDiXEszsFi9FZ4/UJTrTA\nOqxPPWV01jiHEOIBNqApA5xHHQW4IMSzKw2ZmVjvdTkFZU/WX13qJzOx3nmIEA+woWVZNOeKBCn2\nYN3Vo2YCPEmIZzeaMjOz/usRqNiS9VaPOsnMrH/uI8SzC4UJ7IOaSinCFVVZY3Wpj2AfcDchHmBH\nGnRdQhY1WFd1qYsA9xPi2ZzmDI+zJ+pyYsparKX61EN4nD3BbYR4NqUQwe3sjfqEL85h/dSnDsLt\n7A1uEuLZjAIE97NH6hPEOIV1U5/6B/ezR7hOiAdoiCZdn0DGMayX+tQ9gOMI8WxCg4bD2S/1CWYc\nwjqpT72Dw9kvXBHiARq0LItmXZmAxn2sj7rUOIDTCfFUp0nD6ewfYDTqGpzO/iER4qlMoYHz2Uf1\n+JVh3GRN1KWewfnsI4R4gA5o2HUJbSTWQW3qGMA6hHiAThiA6xLg5ub5r0v9AliPEA/QEYNwXYLc\nnDzvdalbAOsS4qlG04Y67K26BLq5eL7rUq+gDntrbkI8QIc077oEuzl4nutSpwDqEOKpQuOG+uyz\nunxK+bg8t/WpT1CffTYvIZ7VKSiwnWVZ7LnKhL2xeD7rUpNgW/bbnIR4VqWQwD7svbqc3PbPc1if\nOgT7sPfmI8QDDEITr08I7JPnrT71B2A7Qjyr0cBhf/ZhfQJhXzxf9ak7sD/7cC5CPMBgNHJgK+oN\nwPaEeFahiUNb7Mm6nO72wfNUlzoDbbEn5yHEAwxKM69LQGyb56cu9QVgP0I8Z9PIoV32Z12CYps8\nL3WpK9Au+3MOQjzA4DT0ugTGtng+6lJPAPYnxHMWzRz6YK/WJTi2wfNQlzoCfbBXxyfEA0xCU69L\ngNyXx78u9QOgHUI8J9PQoT/2bV2C5D487nWpG9Af+3ZsQjzAZDT2ugTKbXm861IvANrzaO87QJ80\ndejb1R4WgOq4/riql+uzbuuzbqF/y7Kol4NyEs/RNHYYh/1cnwFqXR7P+tQFGIf9PCYn8QCT80p9\nfaWUTQepLZ/PUX+uWRn4AdonxHMUzR3GJMj3o7Xn6ZD7o3f0wfMEY9LjxyPEA5BEk6/tmNP40Z6H\n+36eWR+T1gjwAP0Q4jmYBg/jE+Trui3Iz/543/bze4y2pb/D+PT3sQjxADxGo6/LY/swj9F2BHiA\n/vh0egCeYLCH8dnnAH1yEg/ArfwueRiT8A7QNyfxHETDh3nZ/zAO+xnmZf+PQ4gHAACATgjxPMir\ndoA6AP2zjwF1YAxCPAAH0fihX/YvwDiEeO6l6QPXLcuiLkBH7FngJjWhf0I8d7LBgbuoD9A++xS4\ni/rQNyEegJMYAKBd9ifAuIR4AE4mKEB77EuAsQnx3MoAABzKe26hDfYicAz1ol9CPACrMAzAfuw/\ngHkI8TzBIACcSv2A7dl3wKnUjz4J8QCsykAA27HfAOYjxPMYwwCwBrUE6rPPgDWoJf0R4gGowlAA\n9dhfAPMS4vkSAwEAAMxHDuiLEA9ANYYCWJ99BTA3IR6AqgQOWI/9BIAQTxJDAVCXGgPns4+AmtSY\nfjza+w6wPxsW2ML1WlNK2fGeQD/0aGBLy7Lo0R1wEg/A5gQTeJh9AsBthHgAdiGgwN3sDwDuIsRP\nzpAA7EkNgifZF8Ce1KD2eU88ALvwnju43dXeMEgDcBsn8RMzHAB7EeDhYfYJsBc5oW1O4gHYjFAC\nx3EqD8BNTuIB2IQAD6ezfwC4IsRPyiv6wFZKKQIIrMBeArYkL7RLiAegGoED1mdfAcxNiJ+QV9WA\nLQgaUI/9BWxBbmiTEA/A6gQMqM8+A5iTEA/AqgQL2I79BjAfv2IOgLMJErCf6/vPpa8A43MSPxnN\nHVibAA/tsB+BtckP7XESD8BJhAVo09XeNHgDjMlJPABHE+ChffYpwJiE+Il4RR5Yg2AA/bBfgTXI\nEW0R4gE4mEAA/bFvAcYixANwEEEA+mX/AoxDiJ+ES2CAcwgA0D/7GDiHPNEOIR6Aexn8YRz2M0D/\nhHgA7mTgh/HY1wB9E+In4NIX4BQGfRiX/Q2cQq5ogxAPwBMM+DA++xygT0I8AI8x2MM87HeA/gjx\ng3PJC3AMAz3Mx74HjiFf7E+IH5gNBhyqlGKQh4mpAcAx5Ix9CfEAkzO4A1fUA4D2CfEAEzOwAzep\nCwBtE+IH5RIX4CEGdeAu6gPwEHljP0I8wIQM6MBD1AmANgnxAJMxmAOHUi8A2iPEAwAAQCeE+AF5\nfwpwF6dqwLHUDeAucsc+hHiASRjEgVOpHwDtEOIBJmAAB86ljgC0QYgfjEtagJsM3sBa1BPgJvlj\ne0I8wMAM3MDa1BWAfQnxAIMyaAO1qC8A+xHiB+JSFuCKARuoTZ0Brsgh2xLiAQZjsAa2ot4AbE+I\nH4RXv4DEQA1sT90BEnlkS0I8wCAM0sBe1B+A7QjxAAMwQAN7U4cAtiHEA3TO4Ay0Qj0CqE+IH4D3\nnwAAAHuTS7YhxAN0zKkX0Bp1CaAuIR6gUwZloFXqE0A9QnznXLICAAC0Qj6pT4gHAACATgjxAB1y\nqSrQOnUKoA4hHqAzBmOgF+oVwPqEeICOGIiB3qhbAOsS4gE6YRAGeqV+AaxHiAfogAEY6J06BrAO\nIb5jfn0DAADQGjmlLiEeoHFOr4BRqGcA5xPiARpm4AVGo64BnEeIBwAAgE4I8Z3yPhMYn9MqYFTq\nG4xPXqlHiAdokAEXGJ06B3AaIR4AAAA6USXEl1JeXUr5wVLKr5ZSvlBK+flSymtv3OadpZTnL7//\n46WUr7zx/ZeWUt53+Xd8rpTyXCnllTXuL0BLnE4Bs1DvAI63eogvpbwiyQeT/L8kfzLJ70/yt5N8\n5tptvjPJ25J8W5LXJ/l8kg+UUl5y7a96Nsk3JvmmJG9K8uokP7T2/e2R95fAuAy0wGzUPRiX3FJH\nWfuBLaV8T5I3Lsvydffc5vkk37ssy/dd/vnlSV5I8s3Lsrz/8s+/kuQty7L88OVtnkny0SRfsyzL\nh275O1+b5MOr/jCNshlgXIZZYEZmGxjXRLPN65Zl+cgW/1CNy+n/dJKfK6W8v5TyQinlI6WUb7n6\nZinlK5K8KslPXH1tWZbPJvnZJG+8/NJXJ3l04zYfS/LJa7cBGMpETQ7gMeofwOFqhPjfneRvJvlY\nkj+R5J8keW8p5a9cfv9VSZZcnLxf98Ll95Lk6SQvXob7u24zJa9Uw5gMsMDs1EEYk/yyvkcV/s6n\nknxoWZa/d/nnny+l/MEk357kByv8ewAAADCFGifx/z0X712/7qNJXnP5/z+dpOTitP26py+/d3Wb\nl1y+N/6u2wAAAMBUaoT4DyZ55sbXnkny35JkWZZfzkUQf/PVNy/D+huS/PTllz6c5Ndv3OaZXLwQ\n8DMV7jPAblxCCnBBPQR4WI3L6b8vyQdLKd+d5P25COffkuRbr93m2STvKKV8PMknkrwryaeS/Ehy\n8UF3pZQfSPKeUspnknwuyXuTfPC2T6afhfeTAAAAvVmWxYt0K1o9xC/L8nOllD+f5HuS/L0kv5zk\nO5Zl+VfXbvPuUsrLknx/klck+ckk37Asy4vX/qq3J/likueSvDTJjyV569r3F2BPGhrA40opDi4A\n7rH674nfywy/J36U5wr4MiEe4ElmHhjPBDNP178nHoADTNDMAE6iPgLcTYgHAACATgjxnXBZGYzF\nKRPA/dRJGIs8sx4hHgAAADohxANszOkSwGHUS4AnCfEAAADQCSEeYENOlQCOo24CPE6I74APgYAx\nGEQBTqN+whjkmnUI8QAAANAJIR5gA06RAM6jjgJcEOIBAACgE0J847xvBPrn9AhgHeop9E++Od+j\nve8A97verCx4AACgR16EW48Q3xGBHvqjYQGsq5RiDoJOmIPqEOI7dbUhNDEAAKAlwntdQnznbm4Q\noR4AANiS0L4tIX4wLrmHdmhoAHW4pB72Z87ZjxA/MIEeAABYi+DeBiF+EgI9bEuTA6jLaTxsw0zT\nHiF+Qj4UDwAAuI/w3i4hfmL3bUwBH06n6QFsw2k8nMfM0ichnlu5/B4AAMYjuPfvqb3vAO0rpdjs\ncCB7BWBb6i4cxkw/DifxHMzpPAAA9ENoH5MQz0kEegAAaI/gPj4hnrPdLBRCPQAAbENon48Qz+qc\n0gMAQD2C+9yEeKryO+kBAGAdwjuJEM9GnM4zA40VYB9+XzwjM19wkxDP5gR6AAC4m+DOfYR4diXQ\nAwCA4M7hhHiaIdDTM40XYF8uqadH5gdOIcTTJIEeAIARCe6cS4ineQI9AAA9E9xZkxBPV+4qgMI9\ne9KYAdrgknr2ZiZgC0I8Q/D76AEA2IvwzpaEeIbi0nsAALYguLMXIZ5h3SysQj01aOAAbXFJPbXo\n+bRCiGcaTukBADiG4E6LhHimJNADAHAbwZ3WPbX3HYC9lVIUa05i3QC0SX3mFGZCeuEkHi45nQcA\nmIvQTo+EeLiFQA8AMCbBnd4J8fAAgR4AoG+COyPxnng4gvdKccU6AGibOk1idmNMTuLhBH4HPQBA\newR2ZiDEwwpccg8AsA/BndkI8bCyuxqJcA8AcDphHS54TzxsxHuyxuF5BOiDej0GMxQ8zkk8bMyl\n9wAA9xPa4W5CPOzIB+QBAAjtcAwhHhrilB4AmIXgDqcR4qFRV41NmAcARiK8w3mEeGic0/m2GDwA\n+lJK0T8boH/CeoR46IhADwD0QnCHOoR46JRADwC0RnCH+oR4GIBPuQcA9iC0w/aEeBiQD8UDAGoS\n3mE/QjwMzCX36zKwAPTJh9utQx+ENgjxMAmBHgA4luAO7RHiYUICPQBwF8Ed2ibEw+TuatTCPQCM\nTViHPgl5W5a+AAAgAElEQVTxwK2c1gPAeAR36J8QDzxIoAeAfgnuMBYhHjiKQA8A7RPcYVxP7X0H\ngH6VUqYZEmb5OQFGNUsdn6k3w6ycxANnuzksOKEHgG0I7DAfIR5YnUvuAaAewR3mJsQDVQn0AHA+\nwR24IsQDmxHoAeBwgjtwGyEe2MXVYCLMA8DjhHfgPkI8sCun8wAguAOHE+KBZgj0AMxEcAdOIcQD\nTfJr6wAYjdAOrEGIb8TNgKLIw+Nu2xNbBXv7EWAMpRS9AxoiA51GiG/U9QVtMcPtXH4PQEvMbPAw\nM9v5hPgOCPTwMIEegD2YzeBhZrN1CfGdEejhYQI9ADWZweBhZrB6hPiOXW0MjQTu5vfRA7AWMxc8\nzMxVnxA/AB8IAQ9zOg/AKcxVcD9z1faE+AG55B7u59fXAXAXsxM8zOy0LyF+cAI9PMwpPcDczEjw\nMDNSO4T4iXgPPTxMoAeYg3kIDmMeao8QPyGn83AY+wNgPGo7PExwb5sQPzmBHgAAENz7IcTzJT7l\nHgAA5iC090uI505O6QEAYByC+xiEeA7iQ/EAAKBPwvtYhHiOclcBEO4BAGBfwvochHhW4dJ7AADY\nnuA+HyGe1Qn0AABQj+A+NyGeqgR6AAA4n+DOFSGezQj0AABwOMGd2wjx7EKgBwCAJwnuPESIZ3cC\nPQAAMxPcOYYQT1MEegAAZiC4cyohnmYJ9AAAjERwZw1CPF0Q6AEA6JHgztqEeLoj0AMA0DLBnZqE\neLom0AMA0ALBna0I8QzjqnAK8wAAbEV4Z2tCPMO5r5AK+AAAHEtQpyVCPFNx+T0AAIcQ3GmVEM+0\nXH4PAMBNwjutE+KZntN5AIC5Ce70RIiHawR6AIA5CO70SoiHOwj0AABjEdwZgRAPBxDoAQD6JLgz\nGiEejiTQAwC0TXBnZEI8nEGgBwBog+DOLIR4WIlADwCwLcGdGQnxUIFADwBQh+DO7IR4qEygBwA4\nj+AOXybEw4auGpAwDwDwMOEdniTEww6czgMA3E5wh/sJ8bCz+xqVgA8AjEhQh9MJ8dAwl98DACMR\n3uF8Qjx0wOX3AECvBHdYlxAPnRHoAYDWCe5QjxAPHRPoAYBWCO6wDSEeBnGzcQr1AEBNQjvsQ4iH\nQflQPACgBuEd9iXEN6KUoiBShUvuAYBzmVOpzZx6OCEeJiLQAwCHEtyhTUI8TEqgBwBuEtyhfUI8\nINADwMQEd+iLEA88RqAHgPEJ7tAvIR64k0APAOMQ3GEMQjxwkLsav3APAG0R1mFsQjxwFr+PHgDa\nILzDHIR4YBU3BwehHgDqEtphTkJ8Q0opijHD8H56AFifWZERmRWPI8QD1Qn0AHA6wR24TogHNuU9\n9ABwGOEduI0QD+zC6TwAPElwBx4ixAO7E+gBmJngDhxDiAea4lPuARid0A6cQ4gHmuaUHoARCO7A\nWoR4oBs+FA+A3gjvwNqEeKA7TucBaJngDtQkxDemlKLwwxEEegBaYH6D05jfjifEA8O4a4DSHABY\ni7AO7E2IB4bntB6AcwjuQEuEeGAqAj0AhxDcgVYJ8cC0BHoArhPcgR4I8QAR6AFmJbgDvXlq7zvA\nkwQI2NeyLIY6gMGp9bA/uec0TuIB7uB0HmAsQjswAiEe4AACPUCfBHdgNEI8wJFuDoRCPUA7hHZg\ndEI8wJmc0gPsS3AHZuKD7QAAAKATQjwAAAB0QohvlEtyoU8u6QTYlroLfZJ3TifEAwAAQCeEeICV\nORUC2IZ6C8xIiAcAAIBOCPEAFTgdAqhLnQVmJcQ3zIc9AAAAo5FzziPEA1TilAigDvUVmJkQDwAA\nAJ0Q4gEAAKATQjxARS75BFiXugrMTogHAACATgjxjfPJjdA/p0YA61BPoX/yzfmEeAAAAOiEEA+w\nAadHAOdRRwEuCPEAAADQCSG+A943AmNwigRwGvUTxiDXrEOIB9iQQRTgOOomwOOEeAAAAOiEEA+w\nMadKAIdRLwGeJMQDAABAJ4T4TvgQCBiL0yWA+6mTMBZ5Zj1CPAAAAHSieogvpXxXKeU3SinvufH1\nd5ZSni+lfKGU8uOllK+88f2XllLeV0r51VLK50opz5VSXln7/gJsxSkTwO3UR4C7VQ3xpZQ/kuTb\nkvz8ja9/Z5K3XX7v9Uk+n+QDpZSXXLvZs0m+Mck3JXlTklcn+aGa9xcAAABaVi3El1J+a5J/nuRb\nkvzvG9/+jiTvWpbl3yzL8l+S/NVchPQ/d/nfvjzJX0/y9mVZ/tOyLP85yV9L8rWllNfXus8AW3Pa\nBPA4dRHgfjVP4t+X5EeXZfkP179YSvmKJK9K8hNXX1uW5bNJfjbJGy+/9NVJHt24zceSfPLababj\nwyAAAIDeyDHrelTjLy2lvCXJH85FGL/pVUmWJC/c+PoLl99LkqeTvHgZ7u+6DcAQlmXR3ADiFB7g\nEKuH+FLK78zF+9m/flmWX1v77wcAAIBZ1bic/nVJfkeSj5RSfq2U8mtJvi7Jd5RSXszFaXrJxWn7\ndU8n+fTl//90kpdcvjf+rtsAAADAVGqE+H+f5A/l4nL6r7r838/l4kPuvmpZll/KRRB/89V/cBnW\n35Dkpy+/9OEkv37jNs8keU2Sn6lwn7vhklsYk0tIgdmpgzAm+WV9q19OvyzL55P84vWvlVI+n+R/\nLsvy0csvPZvkHaWUjyf5RJJ3JflUkh+5/Ds+W0r5gSTvKaV8Jsnnkrw3yQeXZfnQ2vcZoAXeGw/M\nSoAHOFyVD7a7xWOVeVmWd5dSXpbk+5O8IslPJvmGZVlevHaztyf5YpLnkrw0yY8lees2dxcAAADa\nU0Z55bOU8tpcXIY/hVGeN+BJTuOBmZhpYFyTzTSvW5blI1v8QzV/TzwAJzDQArNQ7wCOJ8QDAABA\nJ4R4gAY5nQJGp84BnEaI79Rk7y+BKRlwgVGpbzA+eaUeIR4AAAA6IcQDNMxpFTAadQ3gPEI8QOMM\nvMAo1DOA8wnxHfM+EwAAoDVySl1CPEAHnF4BvVPHANYhxAN0wgAM9Er9AliPEA/QEYMw0Bt1C2Bd\nQjxAZwzEQC/UK4D1CfEAHTIYA61TpwDqEOIBAACgE0J85/z6BgAAoBXySX1CPECnXKoKtEp9AqhH\niAfomEEZaI26BFCXED8Al6wAAAB7k0u2IcQDdM6pF9AK9QigPiEeYAAGZ2Bv6hDANoR4gEEYoIG9\nqD8A2xHiB+H9J0BikAa2p+4AiTyyJSEeYDAGamAr6g3A9oT4gXj1C7hisAZqU2eAK3LItoR4AAAA\n6IQQDwAAAJ0Q4gfjUhYAAGAr8sf2hHgAAADohBAPAAAAnRDiB+SSFkAdAGpTZwB1YB9CPAAAAHRC\niAcAAIBOCPEAAADQCSF+UN6fAvOy/4GtqDcwL/t/P0I8AAAAdEKIBwAAgE4I8QNziQsAALA2OWNf\nQvzgbDCYiz0PbE3dgbnY8/sT4gEAAKATQjwAAAB0QoifgEteAACAc8kVbRDiAQAAoBNCPMAgvDoO\n7EX9AdiOED8JzRUAADiVPNEOIR4AAAA6IcQDAABAJ4T4ibgEBsZlfwN7U4dgXPZ3W4R4AAAA6IQQ\nDwAAAJ0Q4ifjUhgAAOBQ8kN7hHiAzmmuQCvUI4D6hHgAAADohBAPAAAAnRDiJ+RSNwAA4CFyQ5uE\neICOaa5Aa9QlgLqE+ElpsAAAwF3khXYJ8QAAANAJIR4AAAA6IcRPzCUyAADATXJC24R4gE5psECr\n1CeAeoT4yWmyAADAFfmgfUI8AAAAdEKIBwAAgE4I8bhkBgAAkAs6IcSTxIYFAICZyQP9EOIBAACg\nE0I8AN1YliXLsux9N+Bg1iwAaxPi+RKX0EA/Zt+vghGts0bVKeiJ/dqXR3vfAQA41VVIMnzQitmD\nOwD1OYnnMQZhoFX3hSPBiRZYo0CPzP/9cRIPwBCcyrMXAR2ALTmJ5wkGYKBn3ovMVqw1oHfm/j4J\n8QCd0XAPI2BRi7V1OPUKYH0up+dWpRQDCtCMc+rR9f9WoOBUa/TEZVmsQaAZ6lG/hHgApuF98xzL\nC9oAtEaIB2A6wjwPEd4BaJX3xHMnwy20x75cl/c2c5M1sT51C9pjX/bNSTz38t54YG9b1KC7/g1D\nzrj26m3eFw/sTQ3qnxAPAHfwoXhj8aI0ACMQ4nmQ03gAgb5X+hfAl+lfYxDiATqh8bZDoG+b4N4e\nBwIA6xHiOYjmC+yhh7oj0Lehh7VyxfvigT2oO+MQ4gFgJTeDpIGpnp5COwCsSYgHgEqE+vUI7QBw\nQYgH6IDwN4bbgqjn9kkC+5i8NQ9gHUI8B9N8AdY3e7DXVwDqm6mvzECIB6BJM4e7+372HgexmZ/L\nu/hwOwBOJcRzFKfxAPs6pAZvGQ71BIC2ecFwPEI8AAxGsAaAcT219x2gP17Ng23Zc8Ao1DPYlj03\nJiGekygIQE1OkpmBdQ7UZF4flxAPAAAAnRDiOZlX9wAAoD3m9LEJ8QAN04SB0ahrAOcR4jmLRgys\nzfuEmYn1DqzNfD4+IR4AAAA6IcRzNq/2AQDA/szlcxDiARqlEQOjUt8ATifEswrNGAAA9mMen4cQ\nDwAAAJ0Q4lmNV/8AAGB75vC5CPEAAADQCSGeVXkVENZhLwGjU+dgHfbSfIR4VqeQAKdalmXvuwCb\ns+6BU5m75yTEU4WCAgAA9Zi35yXEAzRGUwZmod4BHE+IpxqNGTiGS4qZmfUPHMOcPTchHgAAADoh\nxAMAAEAnhHgAAADohBBPVd6vA8exZ4DZqHtwHHsGIZ7qFBrgIT7UC+wD4GHmahIhHgAAALohxLMJ\nrxoCAMDpzNNcEeIBGqE5A7NS/wAOJ8SzGQ0aAACOZ47mOiGeTSlAwE0+zAu+zH4AbjI/c5MQz+YU\nIgAAeJi5mdsI8QAN0KSB2amDAIcR4tmFRg0kLh2G29gXQGJe5m5CPLtRmAAA4EnmZO4jxAMAAEAn\nhHh25VVGAAD4MvMxDxHiAQAAoBNCPLvzaiPMyYd3wd3sD5iTuZhDCPE0QcECAGBm5mEOJcQDAABA\nJ4R4muHVR2Zl7QNcUA+ZlbXPMYR4mqKAwRy83xceZp/AHMy/HEuIBwAAgE4I8TTHq5EAAMzA3Msp\nhHiAHWneAI9TFwHuJ8TTJA0cxuV9vnA4+wXGZd7lVEI8zVLYAAAYkTmXcwjxAAAA0AkhnqZ5lZKR\nWd8At1MfGZn1zbmEeJqn0ME4vL8XjmffwDjMtaxBiAcAAIBOCPF0wauWAAD0zDzLWoR4uqHwMRLr\nGeB+6iQjsZ5ZkxAPAAAAnRDi6YpXMaFfPpwLTmf/QL/Mr6xNiKc7CiEAAD0wt1KDEA8AAACdEOLp\nUinFK5t0y9oFOIx6Sa/MqtQkxANQnffzwvnsIwASIZ7OeYUTAICWmE+pTYinewolPbFeAY6jbtIT\n65UtCPEMQcEEAGBP5lG2IsQDUJX38cJ67CcAhHgAAADohBDPMFzCBADAHsyhbEmIZygKKC2zPgFO\no37SMuuTrQnxAAAA0AkhnuF4NRTa4UO4YH32FbTD3MkehHiGpKACAFCTeZO9CPEMS2GlJdYjwHnU\nUVpiPbInIR4AAAA68WjvOwA1Xb1K6v2DAACcywk8LXASzxQUXAAAzmGepBVCPNNQeNmLtQewDvWU\nvVh7tESIBwAAgE4I8UzFq6iwHZ9FAfXYX7Ad8yOtEeKZjkIMAMAhzI20SIhnSgoyW7HWANalrrIV\na41WCfFMS2EGAOA25kRaJsQDAABAJ4R4puZVVqjDh25BffYZ1GE+pHVCPAAAAHRCiGd6Xm2lFmsL\noA71lVqsLXogxEMUbACA2ZkH6YUQD5cUbliH9+nCduw3WIc5kJ4I8XCNAs5arCWAutRZ1mIt0Rsh\nHm5QyAEA5mDuo0dCPNxCQQcAGJt5j16tHuJLKd9dSvlQKeWzpZQXSik/XEr5vbfc7p2llOdLKV8o\npfx4KeUrb3z/paWU95VSfrWU8rlSynOllFeufX/hLgo7HM/7c2F79h0cz5xHz2qcxP/RJP8wyRuS\nfH2S35Tk35VSfsvVDUop35nkbUm+Lcnrk3w+yQdKKS+59vc8m+Qbk3xTkjcleXWSH6pwf+FOCjwA\nwFjMd/Su1H71tpTy25P8jyRvWpblpy6/9nyS712W5fsu//zyJC8k+eZlWd5/+edfSfKWZVl++PI2\nzyT5aJKvWZblQ7f8O69N8uGqPwxTcsLBsWYeDuwX2Ie6A4ebeb9Q1euWZfnIFv/QFu+Jf0WSJcn/\nSpJSylckeVWSn7i6wbIsn03ys0neePmlr07y6MZtPpbkk9duA5tQ6AEAxmCuYwRVQ3y52CXPJvmp\nZVl+8fLLr8pFqH/hxs1fuPxekjyd5MXLcH/XbWAzCj48zGkY7Mf+g4eZ5xjFo8p//z9O8geSfG3l\nfweqK6UYkgAAOiTAM5JqJ/GllH+U5E8l+WPLsvz3a9/6dJKSi9P2656+/N7VbV5y+d74u24Dm9MA\neIg1ArAtdZeHWCOMpkqIvwzwfzbJH1+W5ZPXv7csyy/nIoi/+drtX56LT7P/6csvfTjJr9+4zTNJ\nXpPkZ2rcZziURgAA0AdzGyNa/XL6Uso/TvIXk/yZJJ8vpVyduP+fZVn+7+X/fzbJO0opH0/yiSTv\nSvKpJD+SXHzQXSnlB5K8p5TymSSfS/LeJB+87ZPpAQAAYAY13hP/7bn44Lr/eOPrfy3JP0uSZVne\nXUp5WZLvz8Wn1/9kkm9YluXFa7d/e5IvJnkuyUuT/FiSt1a4v3A074+Hx9kPsL9lWZw6wjX2A6Oq\n/nvit+L3xLOXUfYQ65h1YLAPoA1qEMy7D9jdUL8nHoamUQAAtMFcxgyEeADO4gQM2mE/AoxPiIcV\neNUXAGBf5jFmIcTDSjQOrAGAfai/WAPMRIiHFWkgAADbMn8xGyEeVqaRAABsw9zFjIR4qEBDYRY+\nRAvaY18yC/MWsxLioRKNBQCgDnMWMxPioSINZh6ea4B9qcPz8FwzOyEeKtNoAADWYa4CIR42oeEA\nAJzHPAUXhHjYiMbDaHx4FrTL/mQ05ij4MiEeNqQBjcnzCtAG9XhMnld4nBAPG9OIAAAOY26CJwnx\nsAMNid65VBfaZ5/SO/MS3E6Ih51oTAAAtzMnwd2EeNiRBtU/zyFAW9Tl/nkO4X5CPOxMowIAuGAu\ngocJ8dAADQsAmJ15CA4jxEMjNC564cOyoB/2K70wB8HhhHhoiAYGAMzG/APHEeKhMRpZPzxXAG1S\nn/vhuYLjPdr7DgBPumpoLoMEAEYkvMPpnMRDwzQ4WuOFJeiPfUtrzDdwHiEeAAAAOiHEQ+O8Wt0m\nzwtA29TpNnle4HzeEw8d8B55AKBnwjusx0k8dEQDZE9eRIJ+2b/syfwC6xLioTMaIQDQC3MLrE+I\nhw5piPvy+AP0Qb3el8cf6hDioVMaIwDQKnMK1CPEQ8dKKZokm/B+2v/f3v2HblvddQB/f5xua4Ts\nD5syIlZUrljYZqsk2goHUlERRbqCotiWpSAWrEmTgSMQwx9ZLgYFMgvCjFFCZLkImuVE1/xjP4Ix\n1xruedIabeg2zZ3+uO6v3X59vs/z/XXf13Wu6/UC4bmv+3r8nvv7nPuc877Ouc4F/fM9ZhuMTWDz\nhHiYAZ0lADA24xHYDiEeZkLHCQCMxTgEtkeIhxnRgW6e3zFAX7Tbm+d3DNslxMPM6EgBgG0x7oDt\nE+JhhnSoAMCmGW/AOIR4mCkdKwCwKcYZMB4hHmZMB3u8/D4B+qT9Pl5+nzAuIR5mzvNaOSrPlob5\n8H3mKIwpYBrOHrsAwHbsdLoGcADAQQjuMC1m4mFhdMQAwH4ZN8D0CPGwQDrkg/M7A+ibdvzg/M5g\nmoR4WCgdMwCwF+MEmC4hHhbMBjWciT0UYH58rzkdYwOYPiEe0FkDAMYD0AkhHkii4waAJTMOgH4I\n8QBnYGADMA/ac2AOhHjgee6DY537ZmG+fL/Zoe+H/gjxwIvozAFg/vT30CchHjglHTsAzJd+Hvol\nxAN7ssTOIAdgbpberuvboX9CPHBGOnsA6J/+HObh7LELAPRhveO3IdL8+TeG+WutCXUL4N8Y5sdM\nPHBgBgQAMH36a5gnIR44lCUMDJbwGQGWaAnt+xI+IyyVEA8cms1x5slSelgO3/f50TfD/AnxAAAA\n0Akb2wFHtnPF34wOAIzD7Dssh5l44NgYQADA9ul/YVnMxAPHyqw8AGyH8A7LZCYe2AgDiz65+ALL\n43vfJ/0sLJeZeGBjzMoDwPES3gEz8cDGGXAAwNHpT4FEiAe2pLfn1vZUVgAOrqd2vrc+FNgsIR7Y\nKoOQ6XLbAyyX7/906TeB3dwTD2yde+UB4PSEd2AvZuKB0Ux1eeAUywTA8Ztiez/VvhGYDiEeGJ3B\nCgDoD4H9EeKBSTDzMC63NgDagfHoA4GDEOKBSTGIAWBJ9HvAQdnYDpgcG98BMHfCO3BYZuKByRpj\neaFBFcCyjNHP6GuAoxDigckz2NksKx6AHdqDzdKfAcfBcnqgC5bYA9Ar4R04Tmbiga4YCAHQE/0W\ncNyEeKA7mxoQGWgBLJN+BeiJ5fRAlyyvB2CqhHdgk8zEA12zy+/RuAgC7KZdODx9ErANQjwwCwZN\nAIxJPwRsi+X0wGzsHkAdZDbJ4Atg2apKvwF0wUw8MFuWNZ6eJbPAXrQPe9O3AGMT4oHZM9gC4Djo\nT4ApEOKBRTDwAuAo9CPAVAjxwGLstQTSwAyA5NT9geXzwNTY2A5YHM+YB+BMBHdgqoR4YLGWPEBz\nAQM4k9baItvJJX5moC+W0wMAAEAnhHgAAADohBAPAAAAnRDiARbG/fDAfmkvAKZHiAcAAIBOCPEA\nAADQCSEeAAAAOiHEAwAAQCeEeAAAAOiEEA8AAACdEOIBAACgE0I8wIJ45jNwUNoNgGkR4gEAAKAT\nQjwAAAB0QogHAACATgjxAAAA0AkhHmAhbE4FHJb2A2A6hHgAAADohBAPAAAAnRDiAQAAoBNCPMAC\nuJ8VOCrtCMA0CPEAAADQCSEeAAAAOiHEAwAAQCfOHrsAAGxeVSVxTytwODttCADjMxMPsCAG4sBB\naTcApsVMPMDCrA/IzcwDpyK4A0yXEA+wYLsH6kI9LJPQDtAPIR6A55mlh+UQ3AH6JMQDcEoCPcyP\n4A7QPyEegDOyuz30TXgHmA8hHoB9cw899EFoB5gvIR6AQ7PkHqZDcAdYBiEegGNxugAh4MPxENQB\nEOIB2DjL8OFwhHYAdhPiAdg6oR5OTWgH4EyEeABG5956lkxwB+AghHgAJsUsPXMntANwFEI8AJN2\nqsAj2NMLgR2A4ybEA9AdwZ4pEtgB2AYhHoBZsAyfbRPaARiDEA/ALO0VsIR7DkpYB2BKhHgAFsVS\nfE5HYAdg6oR4ABbvdMFNwJ8fQR2AngnxAHAaluX3S1gHYI6EeAA4hP0EREF/cwR0AJZKiAeADRH0\nD0dAB4C9CfEAMCKBFQA4iLPGLgAAAACwP0I8AAAAdEKIBwAAgE4I8QAAANAJIR4AAAA6IcQDAABA\nJ4R4AAAA6IQQDwAAAJ0Q4gEAAKATQjwAAAB0QogHAACATgjxAAAA0AkhHgAAADohxAMAAEAnhHgA\nAADohBAPAAAAnRDiAQAAoBNCPAAAAHRCiAcAAIBOCPEAAADQCSEeAAAAOiHEAwAAQCeEeAAAAOiE\nEA8AAACdEOIBAACgE0I8AAAAdEKIBwAAgE4I8QAAANAJIR4AAAA6IcQDAABAJ4R4AAAA6IQQDwAA\nAJ0Q4gEAAKATQjwAAAB0QogHAACATgjxAAAA0AkhHgAAADohxAMAAEAnhHgAAADohBAPAAAAnZh0\niK+qq6rqsar6SlU9WFVvHLtMAAAAMJbJhviqujzJzUnek+T1SR5Ncl9VnTdqwQAAAGAkkw3xSa5N\n8v7W2gdaa59KcmWSp5P86rjFAgAAgHFMMsRX1TlJLk7yoZ1jrbWW5P4kl4xVLgAAABjTJEN8kvOS\nvCTJyV3HTya5YPvFAQAAgPFNNcQfxsvHLgAAAACLtLU8OtUQ/2SS55Kcv+v4+UlO7PF3XrPJAgEA\nAMAeXrOtH3T2tn7QQbTWnq2qR5JcmuSvk6SqavX69j3+2n1JfjHJZ5N8dQvFBAAAYNleniHA37et\nH1jDfnHTU1U/n+TODLvSP5Rht/qfS/La1toTIxYNAAAARjHJmfgkaa3dvXom/A0ZltF/LMllAjwA\nAABLNdmZeAAAAOCFprqxHQAAALCLEA8AAACdmE2Ir6qrquqxqvpKVT1YVW8cu0ywH1V1XVU9VFVf\nqqqTVfXBqvrOU5x3Q1U9XlVPV9XfV9W373r/ZVV1R1U9WVVfrqp7qupV2/sksH9V9a6q+npV3bLr\nuHpO16rq1VV116qOPl1Vj1bVG3ado57Trao6q6reW1WfWdXhT1fVu09xnnoOGzKLEF9Vlye5Ocl7\nkrw+yaNJ7lttjAdT98NJ/iDJDyR5S5JzkvxdVX3DzglV9dtJrk7yjiTfn+SpDHX8pWv/n9uS/ESS\nn03ypiSvTvKX2/gAcBCri6zvyNBWrx9Xz+laVb0yyQNJvpbksiTfleS3knxx7Rz1nN69K8mvJfmN\nJK9N8s4k76yqq3dOUM9hs2axsV1VPZjkI621a1avK8l/JLm9tXbTqIWDA1pdfPrPJG9qrX14dezx\nJL/XWrt19frcJCeT/PLqSQ7nJnkiyRWttQ+uzrkwySeT/GBr7aERPgq8SFV9Y5JHkvx6kuuT/Gtr\n7TdX76nndK2qbkxySWvtzac5Rz2na1V1b5ITrbW3rx27J8nTrbVfWr1Wz2GDup+Jr6pzklyc5EM7\nx9pwZeL+JJeMVS44glcmaUn+O0mq6luTXJAX1vEvJflI/r+Of1+GR0aun/NvST4X3wOm5Y4k97bW\n/mH9oHrOTPxkkoer6u7V7VEfraq37bypnjMT/5zk0qr6jiSpqouS/FCSv1m9Vs9hwyb7nPgDOC/J\nS+AfjhkAAAMQSURBVDJc3Vt3MsmF2y8OHN5qFcltST7cWvvE6vAFGUL9qer4Bas/n5/kmVUnudc5\nMKqquiLJ92YYvO2mnjMH35ZhlcnNSX43wzLi26vqa621u6KeMw83Jjk3yaeq6rkMk4K/01r789X7\n6jls2BxCPMzJ+5J8d4Yr2jAbVfXNGS5QvaW19uzY5YENOSvJQ62161evH62q1yW5Msld4xULjtXl\nSX4hyRVJPpHh4uzvV9Xjq4tVwIZ1v5w+yZNJnstwRW/d+UlObL84cDhV9YdJfjzJj7TWvrD21okk\nldPX8RNJXrq6x2yvc2BMFyf5piQfrapnq+rZJG9Ock1VPZNh9kU9p3dfyHBP77pPJvmW1Z+158zB\nTUlubK39RWvt4621P0tya5LrVu+r57Bh3Yf41YzOI0ku3Tm2WpJ8aYZ7dmDyVgH+p5P8aGvtc+vv\ntdYey9ChrdfxczPsZr9Txx9J8r+7zrkww8DxXzZaeNif+5N8T4YZm4tW/z2c5E+TXNRa+0zUc/r3\nQF58K9+FSf490Z4zG6/IMIG27utZ5Qr1HDZvLsvpb0lyZ1U9kuShJNdmaGDuHLNQsB9V9b4kb03y\nU0meqqqdK9f/01r76urPtyV5d1V9Oslnk7w3yeeT/FUybBhTVX+S5Jaq+mKSLye5PckDdnhlClpr\nT2VYdvm8qnoqyX+11nZmLtVzendrkgeq6rokd2cILW9L8va1c9Rzendvhjr8+SQfT/KGDGPvP147\nRz2HDZpFiF89quK8JDdkWIbzsSSXtdaeGLdksC9XZtgA5h93Hf+VJB9IktbaTVX1iiTvz7B7/T8l\n+bHW2jNr51+b4cr4PUleluRvk1y10ZLD0bzgGafqOb1rrT1cVT+TYeOv65M8luSatQ2/1HPm4OoM\nofyOJK9K8niSP1odS6Kew6bN4jnxAAAAsATd3xMPAAAASyHEAwAAQCeEeAAAAOiEEA8AAACdEOIB\nAACgE0I8AAAAdEKIBwAAgE4I8QAAANAJIR4AAAA6IcQDAABAJ4R4AAAA6IQQDwAAAJ34Pwe2E1Yj\nrLioAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1a18772d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(12,12))\n",
    "plt.imshow(Lyot_stop)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bin to final array resolution, crop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 471,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "N_orig = D/2\n",
    "#scalefac = int(N_orig/N)\n",
    "#print(\"Binning the original LS array {0:d}x\".format(scalefac))\n",
    "#Lyot_stop_binned = np.reshape(Lyot_stop, (Lyot_stop.shape[0]/scalefac, scalefac,\n",
    "#                                          Lyot_stop.shape[1]/scalefac, scalefac)).mean(1).mean(2)\n",
    "\n",
    "Lyot_stop_binned = skimage.transform.resize(Lyot_stop, (2*N, 2*N), order=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 472,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(256, 128)\n"
     ]
    }
   ],
   "source": [
    "L_bin = Lyot_stop_binned.shape[0]\n",
    "\n",
    "Lyot_stop_binned_crop = Lyot_stop_binned[L_bin/2-N:L_bin/2+N,L_bin/2:L_bin/2+N]\n",
    "Lyot_stop_binned_crop_binary = np.round(Lyot_stop_binned_crop).astype(int)\n",
    "print(Lyot_stop_binned_crop.shape)\n",
    "\n",
    "Lyot_stop_full_binary = np.round(Lyot_stop_binned[L_bin/2-N:L_bin/2+N,L_bin/2-N:L_bin/2+N]).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 473,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x19f2cfa90>"
      ]
     },
     "execution_count": 473,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbcAAAMuCAYAAACaYU1JAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X+srHddJ/D3x4UWqSmNKLfLIitY7a5ZZbe3W2y0bN2S\nNUJEDYn2rgkrZMOyBeM2m4BkWUHIroqxrUibEMMq+CMGS1jQABURWahI48Wf/FK0tQulVyqsrVRa\noN/9Y+bW4fTee86cM3Nm5jOvVzJJz/M8Pfd5zjn3vs/nPd9npsYYAYBOvmzVJwAAiybcAGhHuAHQ\njnADoB3hBkA7wg2AdoQbAO08bNUnsBdV9egk35HktiSfW+3ZALBCj0jytUluGmP8zekO2ohwyyTY\nfnnVJwHA2viBJL9yup2bUkvetuoTAGCt3HamnZsSbqpIAGadMRc2JdwAYM+EGwDtCDcA2hFuALQj\n3ABoR7gB0I5wA6Ad4QZAO8INgHbmCreqenFV3VJVd1fViap6U1V9w45jfr6qHtjxeOuOY86uquur\n6q6quqeqbqyqxyziggBg3sntsiQ/m+TJSZ6a5OFJfrOqvnzHcW9LciTJ+dPHsR37r0vy9CTPTPKU\nJI9N8sY5zwUATmmudwUYYzxt9uOq+sEkf53kaJL3zuy6b4zxqVN9jqo6N8lzklw5xnj3dNuzk3y4\nqi4ZY9wyzzkBwE4Hfc7tvCQjyad3bL98Wlt+pKpuqKqvnNl3NJNQfefJDWOMjya5PcmlBzwfANj/\n+7lVVWVSL753jPGhmV1vy6RivDXJ1yX58SRvrapLxxgjk5ry/jHG3Ts+5YnpPgA4kIO8WekNSb4x\nybfObhxjvGHmww9W1Z8k+Ysklyd51wH+PADYk33VklX16iRPS3L5GOOTZzp2jHFrkruSXDDddGeS\ns6bPvc06Mt0HAAcyd7hNg+27k3z7GOP2PRz/uCSPTnIyBI8n+UKSK2aOuTDJ45O8b97zAYCd5qol\nq+qGTJb1PyPJZ6vqyHTX344xPldV5yR5aSbPud2ZybT2k0n+LMlNSTLGuLuqXpvkmqr6TJJ7krwq\nyc1WSgKwCPM+5/a8TFZH/s6O7c9O8vokX0zyzUmelclKyjsyCbUfHWN8fub4q6fH3pjk7CRvT/L8\nOc8FAE6pJgsY11tVXZRJnQkASXJ0jPGB0+302pIAtCPcAGhHuAHQjnADoB3hBkA7wg2AdoQbAO0I\nNwDaEW4AtCPcAGhHuAHQjnADoB3hBkA7wg2AdoQbAO1sVLgdP348Y4y1eACwvjYq3ABgL4QbAO08\nbNUnsKmWVU1W1VI+L8A2MbkB0I5wA6AdteSaWVTdqd4EtpnJDYB2hBsA7Qg3ANrxnFtT8z535zk6\noBOTGwDtCDcA2lFLkkSNCfRicgOgHeEGQDtqSfZlLzWm6hJYFZMbAO0INwDaUUuyNHtdgam+BBbN\n5AZAO8INgHbUkqyclZfAopncAGhHuAHQjlqSjXC66lJdCZyKyQ2AdoQbAO2oJdlo6krgVExuALQj\n3ABoRy1JS+pK2G4mNwDaEW4AtKOWZKuoK2E7mNwAaEe4AdCOWhKiroRuTG4AtCPcAGhHLQln4F3C\nYTOZ3ABoR7gB0I5aEg5otrpUUcJ6MLkB0I5wA6AdtSQskJvBYT2Y3ABoR7gB0I5aEg6BuhIOl8kN\ngHaEGwDtqCVhhdwADsthcgOgHeEGQDvCDYB2hBsA7Qg3ANqxWhJWyApJWA6TGwDtCDcA2lFLwiFT\nRcLymdwAaEe4AdCOWhIOgSoSDpfJDYB2hBsA7aglYUlUkbA6JjcA2hFuALSjloQFUkXCejC5AdCO\ncAOgHbUkHJAqEtaPyQ2AdoQbAO2oJWEfVJGw3kxuALQj3ABoRy0JZ6B+hM1kcgOgHeEGQDtqSdhB\nFQmbz+QGQDvCDYB21JIQVeRejDFOud3XjnVkcgOgHeEGQDtqSbaWOm13p6siT3eMrynrwuQGQDvC\nDYB21JJsFbXZ7vZSRe7l//W1ZpVMbgC0I9wAaEctSXvqsd0dpIrcy+f0PeCwmdwAaEe4AdCOWpKW\n1GC7W0YVuZc/y/eGw2ByA6Ad4QZAO2pJ2lB37e4wq8i9nIPvGcticgOgHeEGQDtqSTaaWmt361BF\nno6KkmUxuQHQjnADoB3hBkA7nnNj43huZnfr/Dzb6Xj+jUUyuQHQjnADoB21JGtLNTWfTawiT0dF\nyUGZ3ABoR7gB0I5akrWigppPpyrydFSU7IfJDYB2hBsA7aglWTlV05ltQ/W4VypK9srkBkA7wg2A\ndtSSHDp10u5UkbtTUXImJjcA2hFuALSjluRQqI12p4rcPxUlO5ncAGhHuAHQzlzhVlUvrqpbquru\nqjpRVW+qqm84xXEvr6o7qureqnpHVV2wY//ZVXV9Vd1VVfdU1Y1V9ZiDXgzrpaoefHBqY4wHHyyG\nrynJ/JPbZUl+NsmTkzw1ycOT/GZVffnJA6rqRUlekOS5SS5J8tkkN1XVWTOf57okT0/yzCRPSfLY\nJG/c5zUAwJeYa0HJGONpsx9X1Q8m+eskR5O8d7r5h5O8YozxG9NjnpXkRJLvSfKGqjo3yXOSXDnG\nePf0mGcn+XBVXTLGuGX/lwMAB3/O7bwkI8mnk6SqnpDk/CTvPHnAGOPuJO9Pcul008WZhOrsMR9N\ncvvMMWwoVeTu1GawfPsOt5r863VdkveOMT403Xx+JmF3YsfhJ6b7kuRIkvunoXe6YwBg3w5yn9sN\nSb4xybcu6Fx2dfXVV+dRj3rUl2w7duxYjh07dlinAMAG2Fe4VdWrkzwtyWVjjE/O7LozSWUync1O\nb0eS/MHMMWdV1bk7prcj032nde211+aiiy7azymzRCrI3akgV2Pn193P6vaYu5acBtt3J/n2Mcbt\ns/vGGLdmElBXzBx/biarK393uul4ki/sOObCJI9P8r55zwcAdpprcquqG5IcS/KMJJ+tqiPTXX87\nxvjc9L+vS/KSqvpYktuSvCLJx5O8OZksMKmq1ya5pqo+k+SeJK9KcrOVkgAswry15PMyWTDyOzu2\nPzvJ65NkjPHKqnpkktdkspryPUm+c4xx/8zxVyf5YpIbk5yd5O1Jnj/vyXO4VDrzUUWuH69BuT3m\nvc9tTzXmGONlSV52hv33Jfmh6QMAFsprSwLQjre84YxUN/NRRW4OFWVvJjcA2hFuALSjluQhVDTz\nUUVuPhVlPyY3ANoRbgC0o5YkiSpmXqpIWG8mNwDaEW4AtKOW3GKqyPmoIreDlZM9mNwAaEe4AdCO\nWnLLqFl2p37kJBXl5jK5AdCOcAOgHbXkFlCn7E4VCb2Y3ABoR7gB0I5asilV5O5UkczDysnNYnID\noB3hBkA7aslGVCW7U0WyCCrK9WdyA6Ad4QZAO2rJDacS2Z0qEraPyQ2AdoQbAO2oJTeQKnJ3qkgO\ni5WT68nkBkA7wg2AdtSSa0zFMR9VJKumolwfJjcA2hFuALSjllwzqoz5qCKBUzG5AdCOcAOgHbXk\nGlBFzkcVySawcnK1TG4AtCPcAGhHLbkiaor5qCKBeZjcAGhHuAHQjlryEKki56OKpAsrJw+fyQ2A\ndoQbAO2oJZdMBbE79SPbREV5OExuALQj3ABoRy25BKqGM1NDAstmcgOgHeEGQDtqyQVRRZ6ZKhI4\nTCY3ANoRbgC0o5bcJzXk7lSRcGZu6F4ekxsA7Qg3ANrZqFry6NGjS/vce6nQ1Aa7U0UC68DkBkA7\nwg2AdjaqllwmleP+qSLh4KycXCyTGwDtCDcA2lFLsi+qSGCdmdwAaEe4AdCOcAOgHc+5sWeeZ4PD\n4baAgzO5AdCOcAOgHbUkD6F+BDadyQ2AdoQbAO2oJUmiioR1ZeXk/pjcAGhHuAHQjlpyi6kiga5M\nbgC0I9wAaEctuWVUkbC5rJzcO5MbAO0INwDaUUtumdNVGepKoBOTGwDtCDcA2lFLkuRL60oVJaw/\nKyfPzOQGQDvCDYB21JI8hIoS2HQmNwDaEW4AtKOW5Izc9A1sIpMbAO0INwDaUUuyL1ZUwvpwQ/dD\nmdwAaEe4AdCOWpIDU1EC68bkBkA7wg2AdtSSLJSKElgHJjcA2hFuALSjlmRpVJRw+NzQPWFyA6Ad\n4QZAO2pJDoWKEjhMJjcA2hFuALSjluTQ7VzBpaYEFs3kBkA7wg2AdtSSrJyVlLAc23xDt8kNgHaE\nGwDtCDcA2hFuALQj3ABox2pJ1oqVk8AimNwAaEe4AdCOWpK1paIE9svkBkA7wg2AdtSSbAQVJRzM\ntr3OpMkNgHaEGwDtCDcA2hFuALQj3ABox2pJNo6Vk8BuTG4AtCPcAGhHuAHQjnADoB3hBkA7Vkuy\n0aychPltw+tMmtwAaEe4AdDO3OFWVZdV1Vuq6hNV9UBVPWPH/p+fbp99vHXHMWdX1fVVdVdV3VNV\nN1bVYw56MWy3qnrwAWy3/Uxu5yT5wyRXJTndkxxvS3IkyfnTx7Ed+69L8vQkz0zylCSPTfLGfZwL\nADzE3AtKxhhvT/L2JKnT/4p83xjjU6faUVXnJnlOkivHGO+ebnt2kg9X1SVjjFvmPScAmLWs59wu\nr6oTVfWRqrqhqr5yZt/RTEL1nSc3jDE+muT2JJcu6XwA2CLLuBXgbZlUjLcm+bokP57krVV16Zis\nPz0/yf1jjLt3/H8npvsA4EAWHm5jjDfMfPjBqvqTJH+R5PIk71r0nwcAOy39VoAxxq1J7kpywXTT\nnUnOmj73NuvIdB8AHMjSw62qHpfk0Uk+Od10PMkXklwxc8yFSR6f5H3LPh8A+pu7lqyqczKZwk6u\nlHxiVT0pyaenj5dm8pzbndPjfjLJnyW5KUnGGHdX1WuTXFNVn0lyT5JXJbnZSkkAFmE/z7ldnMlz\nZ2P6+Onp9tdlcu/bNyd5VpLzktyRSaj96Bjj8zOf4+okX0xyY5KzM7m14Pn7OBc4Ja85CdutNuEv\nflVdlEmdCXPbhJ9xWJUNfkWfo2OMD5xup9eWBKAdb3kD0MgGT2ILZXIDoB3hBkA7aknas3KSLlSO\ne2dyA6Ad4QZAO2pJgDWgclwskxsA7Qg3ANpRSwIsmcrx8JncAGhHuAHQjnADoB3PubFVvFoJi+A5\ntPVncgOgHeEGQDtqSYAZKsceTG4AtCPcAGhHLQlsDZXj9jC5AdCOcAOgHbUkW8sN3X2oG9nJ5AZA\nO8INgHbUksBaUzmyHyY3ANoRbgC0o5YEVkblyLKY3ABoR7gB0I5aElgKlSOrZHIDoB3hBkA7aklg\nT9SMbBKTGwDtCDcA2lFLQrz9jcqRbkxuALQj3ABoRy0Jzakc2UYmNwDaEW4AtKOWhA2mcoRTM7kB\n0I5wA6AdtSSsIXUjHIzJDYB2hBsA7agl4RCoGdkEs6+ruuk/syY3ANoRbgC0o5aEA9r0+obtsG1v\n5WRyA6Ad4QZAO2pJOAOVI5tg2yrHvTC5AdCOcAOgHbUk7KCKZNXUjAdncgOgHeEGQDtqSYAVUT8u\nj8kNgHaEGwDtqCUBFkTNuD5MbgC0I9wAaEctCbAHKsfNYnIDoB3hBkA7aklg66kc+zG5AdCOcAOg\nHbUk0Ja6cXuZ3ABoR7gB0I5aEthIKkfOxOQGQDvCDYB21JLAWlAzskgmNwDaEW4AtKOWBJZO5chh\nM7kB0I5wA6AdtSRwICpH1pHJDYB2hBsA7aglgdNSObKpTG4AtCPcAGhHLQlbRM3ItjC5AdCOcAOg\nHbUkNKR+ZNuZ3ABoR7gB0I5aEjaIupFlqqpVn8LCmNwAaEe4AdCOWhJWSM0Iy2FyA6Ad4QZAO2pJ\nWBKVI6yOyQ2AdoQbAO2oJWEfVI6w3kxuALQj3ABoR7gB0I7n3GAHz6fB5jO5AdCOcAOgHbUkW0Xl\nCNvB5AZAO8INgHaEGwDtCDcA2hFuALRjtSTtWSEJp1dVqz6FpTC5AdCOcAOgHeEGQDvCDYB2hBsA\n7VgtSUtWSMJ2M7kB0I5wA6Ad4QZAO8INgHaEGwDtCDcA2hFuALQj3ABox03ctOHGbeAkkxsA7Qg3\nANpRS7LRVJEwv67vvj3L5AZAO8INgHbmDrequqyq3lJVn6iqB6rqGac45uVVdUdV3VtV76iqC3bs\nP7uqrq+qu6rqnqq6saoec5ALAYCT9jO5nZPkD5NcleQhT3hU1YuSvCDJc5NckuSzSW6qqrNmDrsu\nydOTPDPJU5I8Nskb93EuAPAQdZAn5KvqgSTfM8Z4y8y2O5L81Bjj2unH5yY5keQ/jDHeMP34U0mu\nHGO8aXrMhUk+nORbxhi3nOLPuSjJ8X2fKG1ZUALza7Kg5OgY4wOn27nQ59yq6glJzk/yzpPbxhh3\nJ3l/kkunmy7OZJXm7DEfTXL7zDFwWmOMBx8Ap7LoBSXnZ1JVntix/cR0X5IcSXL/NPROdwwA7JvV\nkgC0s+hwuzNJZTKdzToy3XfymLOmz72d7hgA2LeFhtsY49ZMAuqKk9umIfbkJL873XQ8yRd2HHNh\nkscned8izweA7TT3y29V1TlJLshkQkuSJ1bVk5J8eozxfzNZ5v+SqvpYktuSvCLJx5O8OZksMKmq\n1ya5pqo+k+SeJK9KcvOpVkoCwLz289qSFyd5VyYLR0aSn55uf12S54wxXllVj0zymiTnJXlPku8c\nY9w/8zmuTvLFJDcmOTvJ25M8f19XwFawMhIOpsny/z070H1uh8V9bmzCzymss4bhdnj3uQHAOvCW\nN6wt0xqwXyY3ANoRbgC0o5ZkragigUUwuQHQjnADoB3hBkA7wg2AdoQbAO1YLcnKWSEJy9HwJbf2\nzOQGQDvCDYB2hBsA7Qg3ANoRbgC0I9wAaEe4AdCOcAOgHTdxAzSyzTduzzK5AdCOcAOgHeEGQDvC\nDYB2hBsA7Qg3ANoRbgC0I9wAaMdN3Kzc7E2n3pUb5ufG7YcyuQHQjnADoB3hBkA7wg2AdoQbAO0I\nNwDaEW4AtCPcAGjHTdwAG8iN22dmcgOgHeEGQDtqSdaK15kEFsHkBkA7wg2AdtSSABvCCsm9M7kB\n0I5wA6Ad4QZAO8INgHaEGwDtWC25IfZyQ7OVVNCPv9f7Y3IDoB3hBkA7ask1c5DXU5z9f1UZwDYz\nuQHQjnADoB3hBkA7nnM7RN6fbD7e241t5TnzgzO5AdCOcAOgHbXkgqjNANaHyQ2AdoQbAO2oJafU\nisAqWSG5WCY3ANoRbgC0s1Hhdvz48YwxlvIAoI+NCjcA2AvhBkA7Vks25b3dYP35u7k8JjcA2hFu\nALSjlmQjePsbYB4mNwDaEW4AtCPcAGhHuAHQjnADoB2rJQEOkRu3D4fJDYB2hBsA7aglt0C315l0\nQzebpsPfu01jcgOgHeEGQDvCDYB2hBsA7Qg3ANqxWhJgCayQXC2TGwDtCDcA2hFuALQj3ABoR7gB\n0I7VkgALYoXk+jC5AdCOcAOgHbXklvH2N7BYHf4edWRyA6Ad4QZAO8INgHaEGwDtCDcA2rFaEmBO\nVkiuP5MbAO0INwDaUUsC7IEqcrOY3ABoR7gB0I5akja8ziRwkskNgHaEGwDtqCW3WLe3v4FF8/di\nc5ncAGhHuAHQjloSYIYqsgeTGwDtCDcA2hFuALQj3ABoR7gB0I7VksDWs0KyH5MbAO0INwDaUUvS\nkre/YTeqyN5MbgC0I9wAaEctSRJvf8N28LO9PUxuALQj3ABoRy0JtKWG3F4mNwDaEW4AtCPcAGhH\nuAHQjnADoB2rJYFWrJAkMbkB0JBwA6AdtSQP0e11Jr39TX8dfk5ZLJMbAO0INwDaWXi4VdVLq+qB\nHY8P7Tjm5VV1R1XdW1XvqKoLFn0eQG9V9eADdlrW5PanSY4kOX/6+LaTO6rqRUlekOS5SS5J8tkk\nN1XVWUs6FwC2zLIWlHxhjPGp0+z74SSvGGP8RpJU1bOSnEjyPUnesKTzAWCLLGty+/qq+kRV/UVV\n/VJVfU2SVNUTMpnk3nnywDHG3Unen+TSJZ0L0IQqkr1aRrj9XpIfTPIdSZ6X5AlJ/k9VnZNJsI1M\nJrVZJ6b7AODAFl5LjjFumvnwT6vqliR/leT7knxk0X8eAOy09Ju4xxh/W1V/luSCJL+TpDJZbDI7\nvR1J8ge7fa6rr746j3rUo75k27Fjx3Ls2LGFnS+wXlSQ7MfSw62qviKTYHvdGOPWqrozyRVJ/ni6\n/9wkT05y/W6f69prr81FF120zNMFoIGFh1tV/VSSX8+kivwnSX4syeeT/Or0kOuSvKSqPpbktiSv\nSPLxJG9e9LkAsJ2WMbk9LsmvJHl0kk8leW+Sbxlj/E2SjDFeWVWPTPKaJOcleU+S7xxj3L+Ec4Ev\n4XUmN4MqkoNaxoKSXZ8AG2O8LMnLFv1nA0DitSUBaEi4AdCO93PjjLq9txvry88Xi2RyA6Ad4QZA\nO2pJYGVUkSyLyQ2AdoQbAO2oJYFDpYrkMJjcAGhHuAHQjloSWDpVJIfN5AZAO8INgHbUkuxZt9eZ\n9N5uy9XhZ4TNZXIDoB3hBkA7aklgYVSRrAuTGwDtCDcA2lFLAgeiimQdmdwAaEe4AdCOWhKYmyqS\ndWdyA6Ad4QZAO2pJ9qXb60yyO99nNonJDYB2hBsA7aglId7+Zpb6kQ5MbgC0I9wAaEctCagiacfk\nBkA7wg2AdtSSsKVUkXRmcgOgHeEGQDtqSdgiqki2hckNgHaEGwDtqCU5sG5vf9PhdSY7fB/gIExu\nALQj3ABoRy0JTagi4R+Y3ABoR7gB0I5aEjaI6hH2xuQGQDvCDYB21JKw5lSRMD+TGwDtCDcA2lFL\nslDdXmfyMPl6weKY3ABoR7gB0I5aEs5g2W9/o4qE5TC5AdCOcAOgHbUkHAL1IxwukxsA7Qg3ANpR\nS8ICqR9hPZjcAGhHuAHQjloS9uh0N3SrImH9mNwAaEe4AdCOWpKl6Vzddbse6MbkBkA7wg2AdoQb\nAO0INwDaEW4AtCPcAGhHuAHQjnADoB3hBkA7wg2AdoQbAO14bUkORefXmQTWj8kNgHaEGwDtCDcA\n2hFuALQj3ABoR7gB0I5wA6Ad4QZAO27iZmncrA2siskNgHaEGwDtqCXZF5UjsM5MbgC0I9wAaEct\nSRI1I9CLyQ2AdoQbAO0INwDa8ZxbU55DA7aZyQ2AdoQbAO2oJdeMOhHg4ExuALQj3ABoRy25T+pD\ngPVlcgOgHeEGQDsbVUsePXp01acAwAYwuQHQjnADoB3hBkA7wg2AdoQbAO0INwDaEW4AtCPcAGhH\nuAHQjnADoB3hBkA7wg2AdoQbAO0INwDaEW4AtCPcAGhHuAHQjnADoB3hBkA7wg2AdoQbAO0INwDa\nEW4AtCPcAGhHuAHQjnADoB3hBkA7wg2AdoQbAO0INwDaEW4AtLPScKuq51fVrVX191X1e1X1r1d5\nPgD0sLJwq6rvT/LTSV6a5F8l+aMkN1XVV63qnADoYZWT29VJXjPGeP0Y4yNJnpfk3iTPWeE5AdDA\nSsKtqh6e5GiSd57cNsYYSX4ryaWrOCcA+ljV5PZVSf5RkhM7tp9Icv7hnw4AnVgtCUA7qwq3u5J8\nMcmRHduPJLnzFMc/YulnBMAmOWMurCTcxhifT3I8yRUnt1VVTT/+3VP8L197OGcGwIb42jPtfNgh\nncSpXJPkF6rqeJJbMlk9+cgkv3CKY29K8gNJbkvyuUM6PwDWzyMyCbabznRQTRYprkZVXZXkhZnU\nkX+Y5IfGGL+/shMCoIWVhhsALIPVkgC0I9wAaEe4AdCOcAOgnY0Jt65vj1NVL66qW6rq7qo6UVVv\nqqpvOMVxL6+qO6rq3qp6R1VdsIrzXYaq+pGqeqCqrtmxvc01V9Vjq+oXq+qu6fX8UVVdtOOYTtf7\nZVX1iqr6y+n1fKyqXnKK4zb2mqvqsqp6S1V9Yvrz+4xTHHPG66uqs6vq+unPxT1VdWNVPebwrmLv\nznS9VfWwqvrJqvrjqvq76TGvq6p/vONzHNr1bkS4NX97nMuS/GySJyd5apKHJ/nNqvrykwdU1YuS\nvCDJc5NckuSzmVz/WYd/uos1/SXluZl8T2e3t7nmqjovyc1J7kvyHUn+eZL/muQzM8e0ud6pH0ny\nn5JcleSfZXLLzwur6gUnD2hwzedkcgvTVUkesux8j9d3XZKnJ3lmkqckeWySNy73tPftTNf7yCT/\nMsmPZfJv9PcmuTDJm3ccd3jXO8ZY+0eS30vyMzMfV5KPJ3nhqs9tCdf6VUkeSPJtM9vuSHL1zMfn\nJvn7JN+36vM94LV+RZKPJvm3Sd6V5JqO15zkJ5K8e5dj2lzv9Px/PcnP7dh2Y5LXd7zm6d/ZZ8zz\nPZ1+fF+S75055sLp57pk1dc07/We4piLM3mZxcet4nrXfnLbwrfHOS+T34o+nSRV9YRM3ilh9vrv\nTvL+bP71X5/k18cYvz27seE1f1eS36+qN0yr5w9U1X88ubPh9SaTl9G7oqq+Pkmq6klJvjXJW6cf\nd7zmB+3x+i7O5FWiZo/5aJLb0+BrkH/4t+z/TT8+mkO83lW+/NZenentcS48/NNZnunra16X5L1j\njA9NN5+fyQ9Iq7cHqqorM6kxLj7F7m7X/MQk/zmTav1/ZFJRvaqq7htj/GL6XW8ymVbPTfKRqvpi\nJk+B/Lcxxq9O93e85ll7ub4jSe6fht7pjtlIVXV2Jj8DvzLG+Lvp5vNziNe7CeG2TW5I8o2Z/Ibb\nVlU9LpMQf+qYvIh2d1+W5JYxxn+ffvxHVfUvMnn3+V9c3Wkt1fcn+fdJrkzyoUx+kfmZqrpjGug0\nVVUPS/JrmYT7Vas6j7WvJTP/2+NspKp6dZKnJbl8jPHJmV13ZvIcY6frP5rkq5N8oKo+X1WfT/Jv\nkvxwVd2fyW9yna75k0k+vGPbh5M8fvrfHb/Hr0zyE2OMXxtjfHCM8ctJrk3y4un+jtc8ay/Xd2eS\ns6rq3DMcs1Fmgu1rkvy7maktOeTrXftwG/O/Pc7GmQbbdyf59jHG7bP7xhi3ZvKNn73+czNZXbmp\n1/9bSb4pk9/mnzR9/H6SX0rypDHGX6bXNd+ch1boFyb5q6Tt9/iRmfxSOuuBTP/NaXrND9rj9R1P\n8oUdx1zPPvniAAABgUlEQVSYyS897zu0k12QmWB7YpIrxhif2XHI4V7vqlfd7HFlzvcluTfJszJZ\nVvyaJH+T5KtXfW4LuLYbMlkSflkmv8GcfDxi5pgXTq/3uzIJhf+d5M+TnLXq81/g12Hnask215zJ\n84r3ZTK1fF0mdd09Sa7seL3T6/n5TBYKPC3JP81kafhfJ/mfXa45k6XxT8rkl7QHkvyX6cdfs9fr\nm/79vzXJ5Zk0Gjcnec+qr23e683kKa43Z/IL2zft+Lfs4au43pV/web4wl6Vyfu5/X0mKX/xqs9p\nQdf1QCa/4e58PGvHcS/LZGnxvZm8j9EFqz73BX8dfns23Lpd8/Qf+T+eXssHkzznFMd0ut5zMnnP\nxlszub/rzzO5B+phXa45kyr9VH9//9dery/J2Znc53pXJr/w/FqSx6z62ua93kx+gdm57+THT1nF\n9XrLGwDaWfvn3ABgXsINgHaEGwDtCDcA2hFuALQj3ABoR7gB0I5wA6Ad4QZAO8INgHaEGwDt/H/c\noyl985TQzQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1a31a8090>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(20,10))\n",
    "plt.imshow(Lyot_stop_binned_crop_binary)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Write the binned, cropped Lyot stop to a dat (ASCII) file "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 474,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wrote binned, cropped LS array to /astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5flipPad19_N0128.dat\n"
     ]
    }
   ],
   "source": [
    "if flip:\n",
    "    LS_dat_fname_tail = \"LS_half_ann{0:02d}D{1:02d}_{2:s}flipPad{3:02d}_N{4:04d}.dat\".format(\n",
    "                         iD, oD, secobs_key, pad, N)\n",
    "else:\n",
    "    LS_dat_fname_tail = \"LS_half_ann{0:02d}D{1:02d}_{2:s}Pad{3:02d}_N{4:04d}.dat\".format(\n",
    "                         iD, oD, secobs_key, pad, N)\n",
    "\n",
    "LS_dat_fname = os.path.join(LS_dir, LS_dat_fname_tail)\n",
    "\n",
    "if not os.path.exists(LS_dat_fname) or overwrite==True:\n",
    "    #np.savetxt(LS_dat_fname, Lyot_stop_binned_crop, fmt='%.6f', delimiter=\" \")\n",
    "    np.savetxt(LS_dat_fname, Lyot_stop_binned_crop_binary, fmt='%d', delimiter=\" \")\n",
    "    print(\"Wrote binned, cropped LS array to {0:s}\".format(LS_dat_fname))\n",
    "else:\n",
    "    print(\"LS array {0:s} already exists, will not overwrite\".format(LS_dat_fname))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 475,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5flipPad20_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5flipPad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5Pad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5Pad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5Pad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5Pad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5Pad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5Pad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5Pad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5Pad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5Pad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5Pad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5Pad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5Pad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5Pad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5Pad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5Pad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5Pad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5flipPad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5flipPad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5flipPad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann34D90_wfirstCycle5flipPad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5flipPad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5flipPad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5flipPad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann36D90_wfirstCycle5flipPad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5flipPad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5flipPad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5flipPad13_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann38D90_wfirstCycle5flipPad19_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5flipPad01_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5flipPad07_N0128.dat',\n",
       " '/astro/opticslab1/SCDA/Apertures/InputMasks_WFIRST/LS/LS_half_ann40D90_wfirstCycle5flipPad19_N0128.dat']"
      ]
     },
     "execution_count": 475,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#os.listdir(LS_dir)\n",
    "glob.glob('{:s}/LS_half_*wfirstCycle5*_N0128.dat'.format(LS_dir))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 476,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#os.remove(LS_dat_fname)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [conda env:py27]",
   "language": "python",
   "name": "conda-env-py27-py"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}