{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "import matplotlib as mpl\n",
    "import datetime as dt\n",
    "import pandas as pd\n",
    "import f90nml\n",
    "import matplotlib.dates as mdates\n",
    "import netCDF4 as nc\n",
    "\n",
    "mpl.rc('xtick', labelsize=8)\n",
    "mpl.rc('ytick', labelsize=8)\n",
    "mpl.rc('legend', fontsize=8)\n",
    "mpl.rc('axes', titlesize=8)\n",
    "mpl.rc('axes', labelsize=8)\n",
    "mpl.rc('figure', titlesize=8)\n",
    "mpl.rc('font', size=8)\n",
    "mpl.rc('text', usetex=True)\n",
    "mpl.rc('text.latex', preamble = ','.join(r'''\n",
    " \\usepackage{txfonts}\n",
    " \\usepackage{lmodern}\n",
    " '''.split()))\n",
    "mpl.rc('font', family='sans-serif', weight='normal', style='normal')\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mesozooplankton\n",
    "      REAL(wp), dimension (1:3) :: zz_rate_mesozoo_sumpeakval !uM N magnitude of mesozooplankton summer concentration peaks\n",
    "      REAL(wp), dimension (1:3) :: zz_rate_mesozoo_sumpeakpos ! year-day times of mesozooplankton summer concentration peaks\n",
    "      REAL(wp), dimension (1:3) :: zz_rate_mesozoo_sumpeakwid ! year-days widths of mesozooplankton summer concentration peaks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "nml=f90nml.read('/data/eolson/results/MEOPAR/biomodelevalpaper/modParams/namelist_smelt_ref')\n",
    "nmlold=f90nml.read('/data/eolson/results/MEOPAR/biomodelevalpaper/modParams/namelist_smelt_ref_old') # for now just load the same one"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "#&nampismezo     !   parameters for microzooplankton\n",
    "winterconc       = nml['nampismezo']['zz_rate_mesozoo_winterconc']       #uM N mesozooplankton background concentration\n",
    "summerconc       = 1.0         # uM N mesozooplankton relative summer concentration\n",
    "sumpeakval = np.expand_dims(nml['nampismezo']['zz_rate_mesozoo_sumpeakval'],0)      #uM N magnitude of mesozooplankton summer concentration peaks\n",
    "sumpeakpos = np.expand_dims(nml['nampismezo']['zz_rate_mesozoo_sumpeakpos'],0)   # year-day times of mesozooplankton summer concentration peaks\n",
    "sumpeakwid = np.expand_dims(nml['nampismezo']['zz_rate_mesozoo_sumpeakwid'],0)   # year-days widths of mesozooplankton summer concentration peaks,0)\n",
    "\n",
    "#&nampismezo     !   parameters for microzooplankton\n",
    "winterconc_old       = nmlold['nampismezo']['zz_rate_mesozoo_winterconc']       #uM N mesozooplankton background concentration\n",
    "summerconc_old       = 1.0         # uM N mesozooplankton relative summer concentration\n",
    "sumpeakval_old = np.expand_dims(nmlold['nampismezo']['zz_rate_mesozoo_sumpeakval'],0)      #uM N magnitude of mesozooplankton summer concentration peaks\n",
    "sumpeakpos_old = np.expand_dims(nmlold['nampismezo']['zz_rate_mesozoo_sumpeakpos'],0)   # year-day times of mesozooplankton summer concentration peaks\n",
    "sumpeakwid_old = np.expand_dims(nmlold['nampismezo']['zz_rate_mesozoo_sumpeakwid'],0)   # year-days widths of mesozooplankton summer concentration peaks,0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rates ok\n"
     ]
    }
   ],
   "source": [
    "maxRate=nml['nampismezo']['zz_rate_mesozoo_r']\n",
    "# check that grazing rate has been unchanged since original run:\n",
    "if not maxRate==1.39e-5:\n",
    "    print('PROBLEM: GRAPH BELOW WILL BE INCORRECT BECAUSE GRAZING RATES ARE NOT CONSISTENT')\n",
    "else: \n",
    "    print('rates ok')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "zz_day=np.arange(1,367) # in model, nday_year starts at 1 on jan 1\n",
    "dts=[dt.datetime(2014,12,31)+dt.timedelta(days=float(ii)) for ii in zz_day]\n",
    "zz_day=np.expand_dims(zz_day,1)\n",
    "\n",
    "MesZoBar = winterconc + \\\n",
    "        summerconc*(np.sum ( sumpeakval * np.exp(-(zz_day-sumpeakpos)**2/sumpeakwid**2),1) \\\n",
    "                    + np.sum ( sumpeakval * np.exp( -(zz_day-sumpeakpos-365.25)**2/sumpeakwid**2),1) \\\n",
    "                    + np.sum ( sumpeakval * np.exp( -(zz_day-sumpeakpos+365.25)**2/sumpeakwid**2),1) )\n",
    "    \n",
    "MesZoBar_old = winterconc_old + \\\n",
    "        summerconc_old*(np.sum ( sumpeakval_old * np.exp(-(zz_day-sumpeakpos_old)**2/sumpeakwid_old**2),1) \\\n",
    "                    + np.sum ( sumpeakval_old * np.exp( -(zz_day-sumpeakpos_old-365.25)**2/sumpeakwid_old**2),1) \\\n",
    "                    + np.sum ( sumpeakval_old * np.exp( -(zz_day-sumpeakpos_old+365.25)**2/sumpeakwid_old**2),1) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQoAAACqCAYAAABVoAGTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXl4FVW2t98VZYoQZgQUCINCg4EMRGYFmYd2aLBRgVaZ\nbBFsbcX2aretXP2QwavQCKgt0oIIyiQQQCYZLgo2YQyCDIHLGCINiQRICMn6/qgKHiFzTp06Odnv\n8+QhqXOq6sc5Vb/ae+211xZVxWAwGHIjyG0BBoPB/zFGYTAY8sQYhcFgyBNjFAaDIU+MURgMhjxx\nxChEpKKIdBaR0Tm8PlpEQpw4t8Fg8D6OGIWqJgPxQJXrXxORikDL7F4zGAz+iRtdj5bAv104r8Fg\nKCQ+NQoRiQC2AeLL8xoMhqJxs8PHv94Q6mN1OaKB88A/b9hBxKSKGgwuoarZPsSdbFF0ASJEJFRE\n6ovIi6q6EKtFUTG3HVXVb37+/ve/u66hOGgyuoq/ptxwrEWhqh8BH3lsmmhvTwa6O3Veg8HgfUwe\nhcFgyBOnYxTFno4dO7ot4QYKq+n06dOcOXOGSpUqUa9ePUS8G1P2x88K/FOXP2rKDcmrb+JrRET9\nTVNxRlX59NNPmTBhAqdPn+b222/np59+QkR4/PHHee6556hRo4bbMg1+gIigLgQzDS6TlJTEb3/7\nW95//30mT57M2bNn2bVrF6dOnWLNmjX8/PPPNG3alA8++CDPYJahZGNaFAFKcnIy3bp1IyoqikmT\nJlGqVKls3/fDDz8wYMAAGjRowIwZM6hYMdcBKUMAY1oUJYyMjAz69+9PVFQU77//fo4mAdC0aVO2\nbNlCjRo16NChAydOnPChUkNxwRhFAPLWW2+RmprK5MmT8xWwLFOmDFOnTmXQoEG0bduWPXv2+ECl\noThhuh4Bxu7du+ncuTO7du2idu3aBd7/888/57nnnmPx4sW0adPGAYUGf8V0PUoImZmZDB8+nLFj\nxxbKJAAeffRRZs6cyf3338+aNWu8rNBQXDFGEUB88cUXZGRkMHjw4CIdp2fPnixcuJDHHnuMxYsX\ne0mdoThjEq4ChCtXrvDqq6/y8ccfExRUdP/v0KEDK1eupHfv3ly4cIFBgwZ5QaWhuGKMIkCYM2cO\nDRo08GrGX2RkJN988w3dunUjOTmZkSNHeu3YhuKFMYoAIDMzkwkTJjBp0iSvH7tJkyZs3LiRLl26\nkJyczCuvvOL11G+D/2NiFAHA8uXLKVu2LJ07d3bk+KGhoWzatIm5c+cyevRoMjIyHDmPwX8xRhEA\njB8/npdeesnRJ32tWrXYsGED27Zto1evXvz000+OncvgfxijKObs2bOH+Ph4+vbt6/i5qlSpwpo1\na4iKiiIyMpJ169Y5fk6Df2ASroo5zz77LJUqVWLMmDE+Pe/KlSsZPnw49913H2+//TY1a9b06fkN\n3sckXAUoly9fZs6cOQwZMsTn5+7Rowd79+6lWrVqNGvWjJdeeonjx4/7XIfBNxijKMbMnz+f6Oho\n6tWr58r5K1SowMSJE9m1axepqam0aNGCPn36MGPGDBISElzRZHAGR7oeHov8RKrqhGy2RwA7VHVt\nNvuarkc+6dixI6NGjfJJfCI/pKSksGjRIpYuXcqqVauoXbs2LVu2JCoqivDwcMLCwqhSxaz75K/k\n1vVwLEYhIvWB4ar6Xx7bhgGrVfWoiKxS1W7Z7GeMIh8cP36c8PBwTp06RZkyZdyWcwPp6ens3buX\nbdu2sW3bNnbv3k1cXBwVKlQgLCyMu+66i7CwMMLCwmjevDk332xSetzGb4zC47UI+7Wns3nNGEU+\nmDhxIj/++CMfffRR3m/2E1SVY8eOsWfPHuLi4tizZw87d+7k1KlTdOrUiT59+tC3b19TPMclcjMK\nJ9cIqA+MzeG1sUBIDq+pIW8iIiJ07dq1bsvwCgkJCTpr1ix96KGHNCQkRB955BHdunWr27JKHPa9\nl+397OuVwhCRvsDbQFXg5+x2ev3116/93rFjx2JXsdhp9u/fT0JCAvfee6/bUrzCrbfeysCBAxk4\ncCD/+c9/mD17Ng8//DCNGjXi9ddfp0OHDm5LDEjWr1/P+vXr8/XePLseInIf0BCIBA5jrVK+RlWz\nvck99hsG9AOewjKMfva+4+zjbNfsuyWal6aSzt///ncuXLjA//zP/7gtxTHS09OZM2cOr732Gnff\nfTcTJkwgNDTUbVkBTaFiFHYc4b+AVUAs1k1eBWgAPAxsU9Ub1g71glhjFHlw11138dFHH5WIClSX\nLl1i4sSJTJ48mVdffZU//elPXplGb7iRwhpFfVU9kstBI1R1h5c0eh7XGEUuHDhwgI4dO3LixIkS\ndcMcOnSIJ554gqCgID755BMaNmzotqSAo1CZmbmZhP26103CkDeLFi3igQceKFEmAdCoUSM2bNjA\ngw8+SKtWrZg9e7bbkkoU2bYoRCQEq5uRRRcnuhnZCjItilxp3bo1Y8aMoVu3G1JQSgy7d++mX79+\n3HvvvUyePJly5cq5LSkgKEyLIhp4GegCdAUaOaTNUABOnjx5retRkmnevDmxsbGkpKTQunVrDhw4\n4LakgCfb4VFVXSsi21Q12deCDDnz1Vdf0atXL0qXLu22FNepUKECc+bMYfr06bRr144pU6bQv39/\nt2UFLLnFKIxJ+BmLFi3ioYcecluG3yAiPP3003z99de8+uqrjBgxgtTUVLdlBSR5RsTsPAqDyyQl\nJbF161Z69OjhthS/IzIyktjYWBITE2nbti2HDx92W1LAkZ/Quamk6gesXr2a9u3bc8stt7gtxS+p\nWLEiX375JYMHD6ZNmzbMnz/fbUkBRX6MwgxB+AHLly+nZ8+ebsvwa0SEkSNHsnz5cv7yl7/wxBNP\ncO7cObdlBQSmRVEMyMzMZOXKlcYo8knLli3ZtWsXISEhhIWFmdXOvICpmVkM2L59O4888ogZBiwE\n//u//8uQIUOoX78+48aNo0WLFm5L8lsKXTNTRL4QkXn2zyr752tnZBpyYsWKFfTq1cttGcWS9u3b\ns2fPHvr06UP37t155JFH+O677zAPo4KRq1Go6u9VtT8Qr6rd1KpItd030gxZLF++3BhFEShdujQj\nR47kwIEDtGrVioEDBxIdHc2ECRP48ccf3ZZXLMhX10NEtgGdgGrAl6ra0jFBpuvxK86dO0doaCiJ\niYmULVvWbTkBQWZmJqtXr2bx4sUsWbKEoKAgIiIiaN68OXXq1KFWrVpUr16dSpUqUblyZSpXruyX\n5Qa9TZFL4dlTzp/CGgH50MkJYcYofs3cuXP57LPPWLp0qdtSAhJV5ciRI+zcuZPdu3dz8uRJTp8+\nzdmzZzl//vy1n1KlSl0zjayf0NBQOnToQNeuXalUqZLb/5Ui40rNzMJijOLXPP7447Rq1YoRI0a4\nLaXEoqpcvHjxV8Zx/vx5Dh48yPr169m8eTPdu3fn+eefp3Xr1m7LLTTGKIopmZmZ1KpViy1btlC/\nfn235Rhy4Ny5c8yZM4dx48YRGRnJpEmTimU1LrNSWDFl+/btVK1a1ZiEn1OlShVGjhzJwYMHad26\nNdHR0fzrX/9yW5ZXybVFISKh2W1X1aPOyDEtCk/GjBlDcnIy77zzjttSDAVg9+7d9O/fn969ezN+\n/PhiU2So0F0PEXnb408FKgPDVPWmPE6Y20phw7GK6x7JLihqjOIX2rRpw5tvvknnzp3dlmIoIOfP\nn+fBBx+kZs2afPbZZ8VigaNCdz1U9WVVfRn4AMskqpCPIjb2FPWsYryeDAc+UNWFWKMohhw4e/Ys\nP/zwA+3bt3dbiqEQVK5cmVWrVvHzzz8zePBgMjMz3ZZUJPLKzAwVkWnAS8DbdgJWrrU08yDao8x/\ngyIcJ+BZtWoVnTp1KhHj94FKmTJlWLBgAUePHuWFF15wW06RyKvzFO/xbz8ReVFEXizC+dSux2nI\nA5ONGRgEBwezZMkSYmJimDlzpttyCk1eHaeoIh7/+v7Ov7G6Iz8D53PaqaSvFJaRkcHXX3/N2LFj\n3ZZi8AKVKlVi8eLF3Hvvvdx11120bOlYYnOB8OpKYYUlm5XC+gIf2n+fx1pAaGc2+5X4YObWrVsZ\nOnQoe/bscVuKwYt8+eWXvPLKK+zYsYPy5cu7LecGTMJVMeP111/n0qVLjB8/3m0pBi/zxBNPUKZM\nGT744AO3pdyASbgqZpj4ROAyefJkVq9ezfLly92WUiBMi8LPSExM5M477+Snn36iVKlSbssxOMCq\nVat46qmn2Lt3L8HBwW7LuUZREq6uj6YJoNmtQu4tSrpRzJo1i8WLF7NgwQK3pRgc5NFHHyU0NNSv\nAta5GUVeox5ZSVZzsQrWmPqZDhMTE0Pv3r3dlmFwmHfffZewsDAGDhxIs2bN3JaTJ/mtR9EXa2nB\n81iZlUcdE1SCWxRXr16lRo0axMXFUbt2bbflGBxm6tSpzJs3j/Xr1yPi/jO4yMFMVV0A/AU4B3zp\nRW0GD7Zs2UJoaKgxiRLCU089xblz51iyZInbUvIkPyuF/c5O434ZqwxetPOySiYxMTFmtKMEcdNN\nNzFhwgReeukl0tPT3ZaTK3nN9TiHlSCVjBWf+GM2AU6Dl1i+fLmJT5QwunfvTr169fjwww/dlpIr\neY16ZFcxRU2MwvucOHGC8PBwzpw5w0035TqL3xBg7Nq1i27dunHgwAEqVqzomo6iTDM/Ys8WrW//\nRGDFKgxeZsWKFXTv3t2YRAmkRYsW9OjRg3fffddtKTmS38zMhlgTxO4GVjsnp+Ri4hMlm9dee40p\nU6bwn//8x20p2ZLf4dHQrF+Bvqo60TFBJbDrkZaWRo0aNTh8+DDVqlVzW47BJYYPH07VqlVdS8Ly\nxroeWSXxlF9qVFR0wjBKolGsWbOG1157jW+//dZtKQYXOXbsGBEREezbt48aNWr4/PyFMgo7yWq1\nR0Wq618fpqofeU/mteOWOKN4/vnnqVq1Kn/961/dlmJwmVGjRlG6dGlXCioXtbhuRaxiuElAJaya\nmZWBvzgx+lHSjEJVady4MXPnziUyMtJtOQaXOX36NM2aNXMlO9cbXY/6WDUuzzm5nKB9rhJlFPv2\n7aNbt24cO3bML9J43SI1NZWbb74522rV8fHxpKamEhwcTHBwMOXLl6dcuXIB+3m98MILXL16lUmT\nJvn0vKZwjR8zduxYTp06xT/+8Q+3pfiMhIQEjh49ypkzZ0hMTOTs2bOoKgMGDKBevXo3vH/jxo0k\nJCRw6dIlLl68yIULFxAR/vCHP1CrVi0X/gfOkpCQQNOmTX3eqjBG4ce0atWKt956iy5durgtxWds\n2LCBlJQUbr31VmrUqEH16tUpW7ZsvlsIqsqlS5coU6ZMti2QnTt3Urt2bapXr15sWx1utCqMUfgp\nJ0+eJCwsjDNnzgRUkRpV5dSpU2RkZFC3bl2fn3v58uUcOnQIVaVZs2aEhYVx6623FivTcKNVUWCj\nEJGvgX6qekFEIoDDOY1+5HDCHFcEs49X2T732mz2LTFGMW3aNDZv3szs2bPdluIVUlJS2L17Nzt2\n7CAzM5MOHToQHh7uihZVJTExkbi4OOLi4qhatSoDBw50RUth8XWrojBG0RfoglV7YqeIzFPV/gU4\n4Wh7359FZLqq/tHjuPGqukNEhqrqP7PZt8QYRY8ePRg6dCj9+vVzW0qRSE1NZdmyZRw+fJgmTZoQ\nHh5O3bp1/eYJrqpcuHCBkJDitaSMr1sVhTGKF1V1oj08+j2AvQxgfk/4har+3v59lap2s3+vCKzF\nWt/jL9m1UkqKUSQnJ1OnTh1OnjxJhQoV3JZTJFSV7du3c9dddxW7lc2Sk5MJCQnxG1O7Hl+2Kgoz\nKeyIiITb645WpeDrhOa0IlgDrOUJkwDH6m4WB5YtW0aHDh2KvUmAdYFFRUUVO5MAq9DttGnT2L17\nt1+uDzp69GhmzZrFqVOnXNWRbc1Mu6JV1u8fici2Ah43pxXB+tvmsy63uhYlYaWwefPm0b9/vntz\nrpOZmcnevXvJyMhwLe7gBP369SM+Pp5NmzaxYcMGOnXqRLNmzfymhVGzZk2efPJJxo0b5/VWRZFX\nChOR+7BiCUcLI8AjmHke2IZV+KYf1szThvbbzqvqumz2DfiuR1JSEnXr1uX48eOu1h/ID6rKoUOH\nWL16NWXKlKFz586Ehoa6LcsR4uPjWbt2LVWqVKFv375uy7mGr2IVhYlR/E5VF4rIUHtTfHY3tROU\nBKOYOXMmX331FYsWLXJbSq6cOXOG1atXk5SURNeuXbnzzjv95knrFKpKSkqK33UJ//znP5OZmcl7\n773n2DmKUly3pT0yUdn7skouxaXbsWnTJu68806efvppGjduHPAmAdbN4m8mAfDSSy/x6aefcvr0\naVfOn1OLIgLoDCxQ1SMi0jcrbuHR2gjXbBYZLrKgAG9RnD17loYNG3Ly5Em/XKjWkD1paWns2rWL\nqKgo16qQOd2qKHCLQlV3qOpE2yTqY41SZNHQjmE0cEBrwLNw4UJ69OhhTKKYkZaWxoEDB5g+fTrx\n8fF57+AAbrYqsjUKEbkvq6qVXTfTM4NyO1aQsorj6gKQzz//3K+6HWlpaXz99dckJye7LcWvCQkJ\nYcCAAXTp0oUlS5awcOFCUlJSfKqhZs2aPPHEE4wbN86n54WcYxRRWSMeIhIhIi96vBZvp2Sfz3ZP\nQ44cOXKEuLg4vynJv2/fPqZOnUpqaiqlS5d2W47fIyI0btyYESNGUL58eaZNm+Zzs3CrVZHr8ChW\notXbwO+zqlmJSIidmm1iFAXkjTfe4OzZs65PKU9OTmbFihWcPXuWPn36BOxwp9MkJye7MrztVKyi\nMMOjq4AvAFT1nyJyX9bwqN26EKCKOrCqeaAaRWZmJo0aNWL+/PmuVrJKS0tjypQpREVF0b59+2yn\naRv8m6y8ir1793q1HkdhjCJr4Z8uWGX666tqd/u1YcA8IDq72Z9eEBuQRrF+/XqeffZZdu3a5fow\n4+XLlylXrpyrGgKZc+fOUaWKsyE8J1oV3iiFV1FVk7N+x8q6jHUiCStQjWLQoEFERkby/PPPuy3F\n4CBXr15l+vTp1KpVi+7duzs2upWQkECzZs3YuXMnderU8coxC9OimIY9xdwrCgpAIBpFYmIijRs3\n5vDhw44/aTw5ffp0QJaK83fS09PZsGEDO3bs4L777iMyMtKRVuTf/vY3jh49yqxZs7xyvMIYRQTQ\nEqsCdzWsrkd/+7W+2Jma2dWT8ILYgDOKt99+m4MHD/Lxxx/75HwXLlxg5cqVnD59mmHDhpluhkuc\nOXOGpUuXEhQUxIMPPuj1h0RKSgqNGzdm0aJF3H333UU+XlHK9dcHzmEFLo/Y2+7DmujV0ImK3IFm\nFBkZGTRs2JAFCxYQFRXl6LlUldjYWL755hsiIyO55557AqrEXnEkMzOT2NhYmjRp4khq+IwZM5gx\nYwabNm0qcqulKHM9AH4PeE6l64YVozDzP/LBsmXLqFmzpuMmkZyczCeffMKuXbt4/PHH6dy5szEJ\nPyAoKIjo6GjH5o88/vjjXLx4kfnz5zty/CxyGx49jDUtfLvndPOsoVIRCS3sNPRcBQVYi6J9+/aM\nGjXK8WzMK1eusHfvXsLDw10fVTH4lnXr1jFkyBD27dtH2bJlC32cwsQoKmKlaGfN5+iSlTNhl8dT\n4N8FKY9XALEBYxSbNm3iySef5Mcff3RtIpHBf5k/fz533HEHzZs3L7K59+vXj6ZNmzJmzJhCH8Pb\nw6MRdnHca9u8SSAZRe/evXnggQcYPny4V4+rqqbVEACcOnWKpUuXUq5cOXr37k3VqlULfayTJ08S\nHh7Ohg0baNq0aaGOUZRgpiNp2rkRKEaxc+dOevXqRXx8fJGag55kZmby/fffs3fvXgYPHmzMIgDI\nzMxk69atbNq0idatW9OuXbtCtz6nTJnCvHnz2LBhA0FB+Qk//pqiGMUwoKJdkftFrPkfa7Krnu0t\nAsUoevXqRc+ePRk1apRXjnfq1CmWLVtGmTJl6N27N9WqVfPKcQ3+QVJSEsuXL6dJkyaFTvHPyMig\nXbt2DBkyhGHDhhV4/6IYReesNG3bNCrh0BwPj3MWe6PYsGEDTz75JPv37y/yrMyUlBTWrl3LoUOH\n6NKli1f6swb/JOu6L8r3GxcXR6dOndiyZQsNGzbMewcPimIUv7s+YOlpHrnsl9tKYcPs7Q2yS9gq\n7kahqrRt25ZnnnnGKytT7d+/n2PHjnHvvfcWy3L4Bt8zadIk5s6dy6ZNmwo06a8oeRRVswrYeFA/\nm/ddz3CsFPCFeKwJYmd1HlbVdU5kdfoDn332GWlpaTz22GNeOV6TJk3o1q2bMYkSzP79+9m/fz/5\nfYCOGjWKkJAQ3nzzTa9pyHPUQ0S+4JecCrCK2kzIax/NfqWw6VhZneeBpOxaJsW5RXH+/HmaNm3K\nV199VaiUWjOaYciOI0eOsGLFCsqXL0+3bt2oWbNmnvucPn2ayMhIPv30U7p27Zqv83hjeLQz0BU4\nlJ+WgIjMA4bZBW6uN4rpaq1nem37dfsWW6P44x//SFBQEFOnTi3QfmfPnmXdunXcfvvttG3b1iF1\nhuJMVir4xo0bqVevHp06dcpzOHXjxo08/PDDbN68mUaNGuV5jtyMIl8dGPvJX5DaEzmtFHbI87A5\n7VwcVwpbuXIlMTEx7N69O9/7JCUlsWnTJvbv30+bNm2Ijo52UKGhOJOVCt6iRQu2bt1KTEwMgwYN\nyrUFes899/DGG2/wwAMPsGXLlhvSyIu8UlhRyWGlsL7Ah1gxi8NYXY+AWCnszJkzREREMGfOnHyZ\nWkZGBjExMezfv5/IyEjatWtnZngaCkRBuqlPP/00Bw8eZNmyZbnm9BS56+FLiptRXL16lV69ehEd\nHc1bb72V7/127NhB48aNCQ4OdlCdoaSRmpp6gxlkZGQwYMAALl26xIIFC3KcLGiMwiFUlREjRnD0\n6FGWLl1q6k8aXEVVmTZtGlWqVKF169bUq1fvWqsjPT2dvn37UrZsWWbPnp1tfo8xCocYP348s2bN\nYvPmzYSEhPzqtStXrhAXF8fly5dp166dSwoNJY309HR27tzJ999/fy2u0bx5c0qXLk1qair9+/fn\n8uXLLFiw4IaYhTEKBxg7diwzZsxg3bp112oWqionT55kx44d/PDDD9StW5fWrVtTv35+Uk8MBu+h\nqhw5coTvv/+em2++mX79+gFWV/mZZ55hy5YtLFy48FfZm8YovMjVq1d5+eWXWbZsGevWrbu2DH1G\nRgb//Oc/SUtLIyIighYtWtzQyjAY3CAzM/NXk8RUlalTp/LGG28wderUayZijMJLJCYmMnDgQFSV\nuXPn3jCOnZiYSPXq1U3SlKFYMHXqVMaNG0dYWBgffvght912W5FK4fmc48ePk5GR4baMa6gqM2fO\npGnTplSqVImpU6dmm+xSo0YNYxKGYkPv3r157733SEtLo3Hjxrm+1y9bFLfddhtt27alZ8+eNGzY\nkDp16rjSzz927BgzZ85kxowZXLlyhaeeeoquXbvSrFkzV5aSMxicQFXZuHEjHTt2LF5dj5iYGN55\n5x1iY2Np06YN99xzD88++yy33HLLr9576dIlLly4QPny5SlXrlyBinWkp6eTkpLChQsXKF269K/y\n548fP868efN4//33CQoK4oUXXmDo0KFmIV9DQFNsYxQnTpzg888/JyYmhtjYWKKiomjZsiWRkZE0\nb94csPLZL168yOXLlyldujTBwcE0bdqULl263HDsH3/8kZiYGK5cucLVq1cpX748wcHBBAcHo6ps\n376d1atXk5CQwP3338+wYcNo3bq16U4YSgTF1ig8uXDhAt999x2xsbFs376dPXv28H//93+UL1+e\nevXqUbNmTSpUqEBwcDAVK1akQoUKiAhBQUEEBQWRlpZGUlISycnJXLx4kYSEBE6cOEFCQgL169cn\nPDyc8PBwOnXqRMuWLU0xXEOJIyCMIjtUlcTERI4ePUpiYiLJycnXzODq1atkZmaiqmRkZFCmTBkq\nVKhw7ad27drcfvvt3HbbbaZLYTAQwEZhMBi8R1FXCjMYDCUcYxQGgyFPjFEYDIY8MUZhMBjyxBiF\nwWDIE2MUBoMhT4xR5EF+i4/6En/UBEZXQfBHTbnhiFGISEURGS0ivxORiGxeHy0ixaJYgz9+of6o\nCYyuguCPmnLDqRZFtiuFwbUK3S2xyvkbDIZigFNGEe2x4nmD615ribXuh8FgKCY4ta5HTiuFRQDx\nWC2OL1X1aDb7mvxtg8ElirRSWCHIaaWw+vb2aHv7DcsT5iTUYDC4h09XClPVifZrX2C1KAJyRXOD\nIdDwu9mjBoPB/zB5FAaDIU+MUfgxIpLzirIuIn5aG1BEzJqODuGKUYhIiIi8KiJd3Tj/9YjIHSIy\nVkTquq0FQESCRGQoMF1EQl2Wcw0RuUlE/gD80W0t1yMibYHZItLLbS1Z+Nt1BYW/93xqFPYNEATc\nBZwBHhSRW/LYzUk9N4tIOaARcBZ40i0tth4RkWpYI0PrgXXAH9zUlIX9vTUCrgJNRaSly5IAEJFq\ndpbvKWA60FtEXF0i3t+uKyj6veczoxCR7sDHwH8BJ+0RjxO4cCPYT8bRwBvAi6q6ApgGVBOR9r7W\nY2u6FVgJDADeVtVDwDKghoi4tsqxiNQXkbnAc0ATVZ0DfAP0d0uTrauUiEyzdU3ACsyvxzIMV8zV\nH68rW1eR7z1ftiiigaHAASzBAIuAxiJSx4c6AGoBQar6KnCriAxR1UvA10BP23l9TR1giqpOAq6I\nyGBVPYf1GfURkdtd0ATQGetGnA2MsOMmK4DLLncdbwF2qepfgU3AaHv7HOBOEYkSEV+XUvfH6wq8\ncO/5RLjdR7sFUFX9EqguIl1UdT+QBnwiIg1zPYh3+Q1WExpgCvCwiFRQ1WVAJFZswNddokZAM/v3\nScAjIlIB2A/04peb1GeISCmgLnBZVROxcmJGq+plYAcwTkR+60tNHtwONBORIFWdjXVjdlXVI1jT\nBv4G3ObvP+1cAAAGpElEQVRjTX53XXnr3nNq9uivouKqegwI5pe+2higl4iUBxKAx1T1sBNaPPV4\n6NoMtBWRcPsDWwY8YH9gi4CRqnrRKT05sAS4W0QiVPVH++8woBLwhqq+oqqpTgoQkVoi8pCIlAZQ\n1XTgMPbEPvtJ+RsRqQJUBoao6lInNdm6bhhlUdU4+9cn7H//H3CHiFQGFmBdU8d8ocujpeD6deXU\nvefVhCvxqLVvCw7O+mDsD+sFrCbs3cAxe3apTxCR6sAlDz39sG7E5cBvgX+p6kEfablJVTPsvuN2\nVf3J3v4A0BhIB6KAF1T1jMd+jq5lICLRQAQQq6qxHtv/gdV8LgfcqapvOaUhD32/AVLtVgMiUh8r\nA/gA0A2YqqqbPN7vyOd13XVeHUhX1ST7b1euK6fvPa8YxfVfiN0c7YMVMHnT4z9wJ1YTTIEvnLro\n7eZopsff92M56svAIVXNsLc3A3oAG1R1mxNaPDR4fpHlVPWyWEON+1T13x7vqwjcD2xT1X1OaspG\nYzjWfJybgHV2jAQRuc3e3hb4h931cFqL5+dVCasrcQdwEespmPVaCNAKuElVV/pKk8e2AcBVVZ3n\nsc2V68r+25F7z+sp3CLSGXgaeBtoj/UFvnP9zesEYg0tZnpc4E2xnjZDgd+o6p/s7T5ZZchukoqH\nMZXButifBWYAdwJLVDXJfq9e96U7plOsmbw3q+q/RaS0ql4RkQlAVWA+VlziG6d15KBNVFXtJ2Np\nLOPspKojROSvWNfUG77W5aGvAVYuyRisa3yvqh7Pain6Wo+HLsfuvULHKEQkwm6qIiKVRGSciHSy\nX85y9n8Bj4lIDR+YRBSwGGt8+GYRGQ/8FeiLNTltq4i0AetudFKLracOMBer+YmIPIv1eVTGmhQX\nitVcztKUeb0uh3XeixUMbAY8KiKtsfr51e2f2+xAmK8+r9oi0jfrfGLlaXwIPIbViogVkbuBk0Ap\nEbnZaV1iDcG2y4rZiDWSMgXoYL+lNdZDqJut2ycm4ca9V5RgZkesC60WVpM0Hiu62gJIAsKxvuA+\ndsTc63gGbuw+9V77/HdjjWU/hxX5ro31f60ovhuiSsX6PNrbLZttWE+g5kB5YB4Qh/XE9Bken9lS\noAaWKfwM3K6q54FYoCfwLXDcF3pE5E1gGLDL46V2WEG4y8AGLIN4Gsv4p6nq1euP5UVNWddIVeBR\nrBEpsFo2X2ElLLVS1TVYM6NLiW+HYjvi43uvwDeNx4W2BKiJ9SGuAj7Fau40wroAN2IFeU57Q2g2\nOsoDx0Wkp1iRbrASgZYB3bFu1CFY/bUQYK2qrnSqZWO3Yt4XkVvtp91PWJ9JHNBDVb/F6iN2xnoC\nhWCZhaMjGdeT9RS2I93Hsb6vFOAWEVkE/IA1mhHvo2Z9EPATVqDtZ7ESzwCaAF9ifZ8tsG7IzcCn\nqnrSKTF2PCRJRF63N30FdLH7+A1UdbUdCzlqtzT+rKrTfdGacPPeK7BRXHehHcFyr3pYF39brCfS\nx6p6wMkLTVVTsLLeBvFLXc5oLLdPAAYDl4D3gP92yrA89FzF+vLGA/3sL7UF1pdWSqwRjoZY3aCP\nVTUZ62ZoIO5NZpqHldKbiJ1QpapfqG+HhgUr92AO8Ayw2A7yLsNKQvsMqK2q32E9OW+1A75OcREr\ncFoX+DNQCiuQewk4LyIvi8hm4CtVvWJ/jz7BzXuvSMFM+2Z4FyswFwQkqGqCl7Tl5/y1gUew+vlv\nYjUN47AivndgFcdJ96Ge32D1qZtiXWSPYz2hk4FOwFxV3e3x/hD9pbaozxGRDsB/A0PVShl3DY+R\noFrAX4BXVPWSiFTMuhntYPAVH8QmmgH3YF1HlbEeRN9iBQnT4NqDyjV8fe8V9UnWHmvcPUVV472g\np6AkYGUIVsV6KnbCGpL6N+4U8D1kn7ccVp8xAmti11bgW/1lXDvIDl66aRKlsPI1ujjZ3y8AqSIS\nBozE0lXKHtVI9vi80nyk5SDWiFQVVf2XiPwEvAVMcNsgPPDpvVeUUY+sC62zSyaBHW/YjmV4VYBB\nqrrRDS22nnSsoCVYXYxlQKKqpqvqxaw+ptMjQPnB1rTFT0wCrMBlZ6xY0khVTfZoavv081LVK1jf\nYz0Raa2qK1Q10qmgfEFx494LiFJ49lDkKTfHsD20CFayzTFV3eu2nuKEZ16EWzkSHlqCsIawj/qD\nsbtNQBiFIXBw2yAM2WOMwmAw5ImpmWkwGPLEGIXBYMgTYxQGgyFPjFEYDIY8MUZhMBjyxBiFwWDI\nk/8PwQRT7PcZYYYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fa24c435080>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig,ax=plt.subplots(1,1,figsize=(4,3))\n",
    "ax.plot(dts,MesZoBar_old*maxRate*24*3600,'--',color='grey')\n",
    "ax.plot(dts,MesZoBar*maxRate*24*3600,'-',color='k')\n",
    "ax.set_ylabel('$G_{max}^{MESZ}$ ($\\muup$M N d$^{-1}$)')\n",
    "#ax.set_xlabel('Date')\n",
    "fig.autofmt_xdate(bottom=0.3, rotation=30, ha='right')\n",
    "yearsFmt = mdates.DateFormatter('%b %d')\n",
    "ax.xaxis.set_major_formatter(yearsFmt)\n",
    "fig.savefig('/data/eolson/results/MEOPAR/biomodelevalpaper/figsMod/mesozoo.eps',dpi=400,transparent=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAN4AAACcCAYAAAAQ9xXQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXt0Vdd54H+fpKsXoMe9WA8k8xA2xmFwDAaHxHYH18ZJ\nOgmxjUNsd6VN61fT6ep0Wjt2ms6ym0nHTtLlTFYzK36kk5VOHkBquw+vaYxDhhQHLxPA79gBJEMA\nITB6IhDS1eWbP/Y5h6Or+zj3qSuxf2tp6eronr3PuXd/Z3/7e21RVSwWS3Epm+oLsFguRKzgWSxT\ngBU8i2UKsIJnsUwBVvAslimgohCNikg9sApYqapfd/6+F+gE3lPVVwvRr8UShPjx6Ry7BzM+O4Af\nU+DxWpAZT1UHgS4g7By6F3hSVZ8F7itEnxZLUOLHp4hsADpV9Weq+h2KMF6LpWquVtUh53VHkfq0\nWIKyDugQkVtF5AZgVaHHa0FUzUSISJ3vZiyWUmO3qr4mIluB/kKP10ILnji/d2Gm9SGgP+EbRWwI\njaVgqKokOOweOxB3/JekGa+5UkhV80ZghYgsBJ4GPi0idwOPJjtBVQv28/DDDxe0/WL1MZPupVif\nV8DxuU5EbgUeA54KMl5zoWAznqo+jbkhl68Xqq8LkWg0Sk9PDy0tLYRCoam+nGlHgPFZ0PFq/XjT\nlJ6eHrq6uujp6ZnqS7FkQdGMK1PN2rVrZ0Qfbj8tLS0A3u9C9FFoivV5FRIR+W1gMbAS4/frAn6q\naQwzkkIHLioioqVyLZaZhYigiY0rubS5AvgisBXYw3m/YAfwaYyV9DtJzy+VwW4Fz1IoCiR4v6Wq\n/+68nuR6EJEVmiLixa7xLJYMEZEHgA+JyLdFpI6JRhoAUgkdXEBrPIslj+xV1W3ghZs1ZtqAFTyL\nJXO6RORWVX1WVZ8RyVyLtWs8y4wnfo2XKDvBOf4A8CQmomX6ZSdYLKWMTs6e8QtjmAyyExx3QsZY\nwbNYDKswMZpCZtk0WVlLreBZLngcn9xuzguROtbKIGS1PrLGFcuFjCtoizAq5mpggMyyE7Ka8azg\nWWYc27dvZ/v27ene5mUnqOqzzhrvC5gZ7CngPhHpJ012gutWyBRr1bTMeAoRuZKknzqAdHGaYGc8\niyVn3Kx1jNqpwGfSnWMFz2LJnb2q+lAmJ1jBs1hyp0FEfgn0AajqR9OdkFbwss03StDOCkxMm2S7\nILVYSpRGVV2dyQlJ/XgiskJEtmAciLuBh4BngEHga049ikA4gaSo6s8wpluLZSYhInK/Ux7w1kAn\nJLMkisgiVX0vRU8p843i3lsPbMP4Rx5MNFtaq6alUBTaqulOLC6q+kzac4ox2H1q5joAVf1igvfo\nww8/7P29du3aGVEawDL1FMudkAkJBc/xR4R9h25MlcaethORx1yrj4g8mkzw7IxnKQSlKHjJjCur\ncepGYHwTi3PsZ5NvOn4xx7YslpxIsqnOKmAF8Cpm3Bc0LSjVGq/eSZ8oCnbGsxSKRDOeiCwC7lXV\nLzo7Bb2oqgcdZ/hW4ClVHRKRJ1T1j/J9TUmtmsUUOotlKlHVpx2hW4GZ5TLaZMexaC50Xt8gIltF\n5P5U56RNC8o20c9imYZsBB6E83GXAbkauM1xsd2rqjcBc1OdECRypaQWpRZLOgJmJ4BvbDs2iEeB\nCJlvWtIH/BT4I0wUS1qhTetOEJHfdhzfBcWu8WYOpbavQ5I13j3AbZjSDldhNivpAvZiBPA+jNDt\nVtXX0rS/CJNm1Occ+gywySkdkficAIJ3QzFCvKzgzRwOHz5MV1cXHR0dXHzxxVN9OYUqaLsBY5BJ\nGDopIvc4G6MkPr9UBrsVvJnDdJjx8tTuY0A9xiAzADQAl2CCRR5U1YNJz0012J1YTfcNbtFODRJ9\nnSlW8CyFogghY4sw1s++wGGUQQa7P9okWeRJrljBsxSKUoxcCVplbJ2IzHEke10hL8hiuRAImgh7\nD2aHTHVeWywZU2prv6nEGlcsRWOqrJ1FWOM9inFB/AzYpap/nO4cW/rBUjQKvYvtFNKI2aTyHud3\nWqzgWYpGKBQqFb9eouwELxsB40jPJDthD7BYVV8Vka4g15BS8NzAz3hS+ScsllJHVQcdAbnROeRu\nUjIkIk8AB+L+Tped0AE8KiK7MeFmaUk34/k7VMyUeg9QHqRxi2WasNq3XVcHpnjR3/r+TkcjpibR\nPc7vtKQUPF/W+CJM1HYY45m3WGYU8fuYJ9rXPAV7gA5H1Uxap8hPEFXzQefPx3JRMZ2g1E7nArMu\nI2GxpCOL7IRdnM9GyHTTEoAtwEYnB+9koI7ThIydw+yQ2YVvOyLfNBwIJ6C0P1WWg3UnWApFgOyE\nfsw6rx9T9qGLzLITNmM2Ovk0JqTy82mvKY3grUh0PNMaFM4CdTfmRgYSZTtYwbMUiiL48Z5Q1T9y\ngqYXqWpueyfkucjLblV9zalpkTDN6JFHHvFe2/J+lmnEV53fT2IqrqelWHU178eUfX9NRF5IlN1g\nZzxLoSjFbbqKJXh+B+VAorWeFTxLoSiCqjlhm64gqqaN1bTMeIogeF7B5qCkcyfEb0PrSnTe8/Es\nlmlM3rfpasT4MzZhisCUVDKhxVIiZLxNV9AM9A2YBNh+TAzbwawuL3UfVtW0FIQiqJpbME74LoBU\n1cVcAmWgq9l26EHMVPrjHK7RYpmJbMYRuqAE2RH2VsxsN4CZ7b6e5hSLpeTx71BM7puUDGAmpkZM\nQmz6/tNErvRh4tb8F1IQ44pVNS2FIl7VdJZOXU5Q8z2YsnxeGlCmm5SIyGbXhZDMTx1PuhnvqgTH\nrHRYpjs/BbY5lsiHgKfj0oIyRURkDmbpFmgtmS5k7D2nVXfjkgaM2pk2CNRiKWE6gC9gxvIXAc0w\nDSieB4G/x0xK9wY5IWjph8UYoYtgN5a0lDgB0oI+4zi8f+YENmeaBgQY+4djwVyBMbCAidU8mPbc\ngO6Ehe5LYEOmaUFBsGs8S6FIsMa7kvO7HLupQIHTgHztbFDVZ3y7HQOeFyD1uQEF7zG3Tc6bTevz\nKYBW8CyFophB0kHV1VRbMee0G0qmWMErLBdyMdkCblriT6Z9CjNjvhikwkLSNZ4zhT7mZBYk3A0l\nD9duKRI9PT10dRllpRRK7M0QOlT1oyKyH7jKcUdsAbIXPJhU7MjdDSXv6ztL4ZnBxWRLgW0+zbAz\nyAk2LWiakUhlvJDVyCAUUNXcjbHyr3N+C3BDkIBpW0l6mpFIZZxuauQMelC4G/hs8R3bnOiN8RRV\n8ETkAZzQnGL2O5NwVcVIJMLhw4dpaWmZdmrkdHtQJCOXmkQJsxNE5AUnBAYRWeHWksgFX736cK5t\nTXei0SiHDx8mGo1mfK67/0Bvby9dXV309PR4x6bL7NHS0kJHR8e0eVAUgmRpQU8BXxORKx2pzofb\nYBUB68rPdNwnfk9PT9ZtTOfBO90eFIUgmaq5SFU/77gTOgiotybDScHYTZrSZxdKeb98qIalsvPO\ndMVf2RyTY5pLWlDm/SeyJDrO806nHN89wG1BUh2SdmJy+gYxjsatiRyM1qppKRRJ0oK8yuZ+20OW\naUF3Azdh/NyadVqQP9ZMVZ92zKZZo6rPOmu8L+TSjsWSJ9YBu0WkATMhrMoxLWidqm7M5ISEguek\nAXW5tVXyMfWq6iCQ9axpsQQl4KYl/srm/TmmBYlTtDlwzZVkquatzix1t3OoK1ER2nxiVU1LoUig\navorm28FtgL/qKoH/dnkGbS/gYmb+uRc7GiVsx5rzORCLPklF/eDJSFPA+sc28NjODv9OBNNoJop\ncQxgNnENXBIlmVXzPeep8NX4f/hmwyuD5i1ZcmOmOJxLBWfZE1+0K5ciXveq6k1gfOBA2hkvmXHl\nVZwCR06A9IDv34udNWADYAWvCEy3yJQgzKCwMThfc0XIpeaK37ji1F3xby+7F5N3lI31x2IBZtws\nnnHNlWRrvKtci6YTMna/739dzowYuDaFZSKZrtnyEeky1cTf83SOvEnAClXd6BhlbgxyQrI13h6n\nzsp9mMXnKt//ep3fgfKOLJPJ9Gk/E1TN+HueKZE3TlmUG0RkNUbNXEkOibAPYVIdOlV1UET8Qnav\niAgm2Nmu8bIgU0GaCYN0Jjw8kuBaQZ+C8yUx05FM1bwPs12yOPuX+8s8DGK2nP1pdtdZuhTLbH8h\nBgnP1Ht2LKS/dASuT0S+HeS8oFXG6p0O/Lu77smnU70UHOiHDx+mq6uLjo6OaT/DWM5ThN2C7sGo\nmGHgniARMMkiV76NCRotmipZCoI3w0zcFodkgucGR2PWZhlnJ/giVgT4DCZkbGuQCSmZ4K3AGFTq\ngbmYNCF3U4YNOJEsQcqYBaUUBK/UcR8MkUiE3t5e+4AISCLBczS3pzDLqE+TRXaCiNzgvJwwcIMI\nXkoHuuM872Ni1ng/prDL4kTnWgqHaxkcGBigr6+PgYEBli5daoUvO9zEbAFWZ5OdoKrbwFM11+Gk\nBQHZCV4cGzEzn1vW7yaMPrs36AVa8oO/3gpAX18fPT09dj2aIQkSs3PdtOTGfKUFbcXouy9iKuMe\n9P17q6r+zLefgqVI+N0KS5cu9dajpcxUrJsDpAUtwmhxqzHhkFltWuIjb2lB9c6FuNPujepsRuk4\nDBVjQk3bQVAu1DXeTDfolIKlOMUabwum7MMWsti0xNeWf9MSzVrwEjTsdyesULOTpncsyPkYnXoF\n8KqrG8e9Z9oJXj6EJtXALCWhzPZaSuEeiuBOqMNYNesBglRbT5mP52xnhF/AXFNrUKFz2IiJgvlb\nZtCeC4WuFpao/SBO/kIEAmR7rzPVcR7H0xg18xICGh3TGVdWi8iNqvq3Ph32p5kuQt1dhZxF7YyJ\n8Sx0tbBE7QeJ8yxE5P8MDvnKB/2quk1E1mHWj2lJqWqKyA1xJtMGIOyu9zJFRB4FHk0kuCKiDz/8\nsPd3puX9pkKliUajHD16FIC2trai9BvkPlO9pxRUv2JTBFVzkaq+57jfVmqAjSnTzXj17gvfrHVD\n8renvLgNmEyHCMZ6NAl/Xc1MmYr8rqNHj3Lo0CEAysvL0/aby4Yjfud5OlLNoqWSBzcTHgBO6Qj3\n9QrnZSBDRTrBi4jIwjh3QqCp1I9P6Low/r+sZsxUTKUq1NDQEKjfoBuO+AclGAEfHBxkYGDAc577\n358JU/k5+e8r6AOgxAXUnUUztgqm2x/vaRHZ4qQFvegczrjwkTP1pp1+c2EqUmfa2tooLy8PPCgS\nDfpIJMLAwMCEmcw/KAFvVm1oaKC2tpY5c+ZkLThTmWLkagixWIy2tjYg/QOgVGboRLgqpeNic4Ok\nNwU5N23kiqpudNTLm4ADvtCaC55MB3Gi9584cYK+vj7mzJnDwoULgckCGovFvPcfOnSIBQsW5DwL\npFsHFnLtevz4cZqammZSEvCKuGJHad0Jgbbpcgwsk3xvFwrFUnf8/fgHpSuQ0WiU8vJyRkZGOHbs\nGGNjYyxenF3IbKqZpKenJ6O1a7J7iRfetrY2jh8/ztmzZ9m/fz8f/OAH07YT9OE2xSrpoK8Grbiv\nUyURpKurWbIUwleVrM1kPqyg15DqfW1tbXR0dNDW1pbWV+YOwuHhYQAGBzNxpU4kEonQ0NDA2NjY\npOtqaWmhvb2dhoaGQMaceKLRKO+++y6HDh3i0KFD3v2EQiHmzp0LwOzZs7O+9kRk4mcUkXoRuUFE\n7nd+14vIAyJyq89IkgmbMckEfZg0o37ShJ5N2x1h/esFd0bIlWSzQDJ1J936IxaLceDAAfbt20dd\nXR3PPfccO3fuZHx8nIqKCqqrq2lqaqK5uZnm5mba2tpoaGhg1qxZ3qaTqZ7ep06dYv/+/SxcuDDj\np3xvb69nrKmsrJxw/aFQiMrKSgYGBujt7aW2tjajtnt6eujr66OhoYH6+voJn1traytnzpyhtbU1\nozbTkaFKuhETg7zNV0naSwvCFKfNhAFMYEgjxl2WXV3N6UwuKkeyLy+ZupPo/adPn+a5557j+eef\nZ+vWrTQ2NrJmzRrWrl1LRUUF69evJxQKce7cOUKhEO+99x5Hjx5lz549dHd3c+TIET784Q9z9913\n8/jjj7Nv3z7a2tpob2+nvb2diy66iNmzZ1NWZpSV7u5uxsbGmDdvHjU1NVRVVWFK4qS/V3ftmGiw\nZru2ikajxGIxFixYkHB92Nvb6wllpgKdikzW2wkCOlb7wryyKVuZn4K20wG/RdFPupkwlWDmYiw5\ndOgQ3/rWt/jud7/LmjVruPnmm3n88ceZN29ewj7dGM3Vq1dPciOEw2EOHz7MZz/7Wbq7uzl69ChH\njhxhx44dvP/++1x33XWsWbOG06dPs3PnTjZt2sT777/P2bNniUajVFdXU1NT4/2uqamhtraWSCRC\nU1PTpJ/+/n6am5u56KKLUjrlgxhc3PVhR0dHwvck2ko632uygJuWgJn5HgSezsOmJRkVtA0UJF0M\nkgVJZzqDHTx40LP8JRK8INHyqfo8c+YMnZ2dLF68mJqaGnbs2MGTTz5JVVUVTU1N3HfffSxaNNHV\nmajPVMKYLpL/9ddfZ2BggKqqKlpaWiYIQiwW4+zZs4yMjDAyMuK9PnPmDL29vZw4cWLCz/Hjx73f\nJ0+epK6ujubmZm655RZuuukm3nzzTU6ePMnixYuZP38+YNaW4+Pjnkrq/6moqCAWizFnzhwqKysZ\nGRmhtbWVxsaJXqig9xrv08wmIidJdsIGjIssAmwgt01LFmG2O1DgoSCVxkp+xsvU0drU1ER5ebn3\nRI0vk5BMhRoaGuLdd99l6dKlDA4OJu2zs7OTvr4+jhw5wle+8hVWrlzJpz71KZqbm5MOokR9Jppd\n06l37j22trYyNDTE6Ogohw4dmmB5LC8vZ9asWcyaNSvp+ckeYufOnfOE8/jx457xRlV59913GRoa\nYnx8nNdee42BgQGi0ShjY2OegcZ97f6sW7eOO+64g8cee4xt27Zx5513ekad2tpaFi9eTENDQ8J7\ndYn3aaayxAb19/kCOu7FBHQ8CtwnIv1ksWmJI2gZJcLOmBkvfqZzn6jhcJi+vr60T9ZXXnmFs2fP\nUl1dzcqVKyf0OTQ0xDvvvENVVRVbt24lFAqxa9cu7rrrLkKhEO3t7VRWVia9xnyZut17qqioYHx8\nnPLyclpbW5k/f35G2kBDQwMf+MAHCIVCBTXDu203Nzdz8OBBjh075v1v165d/OAHP+DAgQNceeWV\nXHvttVx77bVcc801E2bHQs14+cTxc2dkXCl5wQuKX/Da2tq89UhTU1OgwkD79u3j2LFjVFVVccUV\nV3gLf1Vl+/btnjHj9ddf56abbmLp0qVp1VqXoDl3kHhgHT9+nEOHDtHb24uIUFNT41lGE/X9jW98\ngx/96EeMj4979wDw5S9/mbq6OgDvWtxr+/nPf84bb7xBVVUV1dXVhMNh1q9fz5o1ayZ9bsePH2d8\nfJxIJEJ1dfWk+020Hjxz5gxvvvmm93Bbvnw5tbW1DA8Ps2vXLnbs2MHu3bsJh8N0dnZy9dVXc/31\n13PdddelnRXTUQTB89RTEXlBA2zFPC0EL9nTLNGgjUQinjro5rkFCUyORqO8+uqrjIyMEA6HCYVC\nbNmyhR/+8IfcfvvtXH/99dTU1NDU1OQNJneAxWIxysvLkxod/OvC2tpafv3rX3s+p4qKCubPn8+e\nPXtYtmwZ1dXVnlC417l582Y2b95MJBJh3rx5fPzjH6epqQmA9vb2SU70/fv3c+LECaqrqz0Lp4jQ\n0tLC6OgowIR7eO211+jp6eHs2bOMjo4yNjbG3Llzqa+vT/iweOSRR3jqqafo7e0lHA5z88038/bb\nb3P//fezfv16Ojs7OXLkCGAEfHR0lP7+fs6cOUN1dTVnz55NuN51AwNEhJdeeont27fz8ssvs2TJ\nEtauXcv111/PkiVLWLRoUUazcxEEbwtwF8aw8o+uhTPlOaUseO4XEovFPEuZfxD4ZxJXwNz3hsNh\nry5J/HvGxsY4cuTIhNni5MmT7Nmzh8HBQb73ve/R1tZGfX09t9xyC6tWraK7u9sLVPZfh3sNYNYf\nr7/+Or29vfT29rJx40Y+97nPTZrxnn32WSKRCHv37uXEiRNcfvnlRKNRrr76as9p3dvb6w3EeOHa\nu3cvp06dAoxfrKamJrCRxv1M6+vrOXTokPcw8ONX01NVMRsbG+Ptt99maGiIUCjEokWLaG1t9a5P\nRLjyyiv55je/ybZt21i/fj2bNm1i48aNPPPMM3zpS1/iE5/4hNefK5Tl5eUAzJ07l7q6Og4cOMCO\nHTsoKyvjIx/5CNu2beP06dOsXr2a5cuXs2TJkpSCWIy0IM4bVx6MSypIfE4pCd7Ro0eZPXu2ZyE7\nevQoXV1dLFiwIGEwsn/mct0IdXV1NDY2TlBx9u/fT319Pd3d3USjUYaHh5k9ezZvvfUWL7zwAm++\n+SbXXXcdd911FwcPHqS6upqWlhZeeeUVtm/fTn9/P6tWreK2226bNBij0Sjbtm3jwIEDnlUwEokQ\niURYtmwZixYtmjTDJlpr+XH/X1VVxejoqLfudGeoXbt2eWpkQ0PDpIdBKtO/27Y7yMPhMMuXL5/Q\nf9B1XzIBfemllzwfoX/N3NXVRWNjI/39/dTV1bFw4UIaGxsnPQxefPFFnnnmGU6ePMnJkydRVW65\n5RbuvPNOKiqMPXDHjh28+uqrvPXWWxw8eJBIJEJra6sXdTN//nyuuuoqOjo6uOyyywoieE7h5y+o\n6infsfuDlH4oKcHz/11XV8cdd9zB7Nmz2bFjB6Ojo5SXl1NRUcHJkyc5fPgwc+bM4WMf+xg/+clP\n2LhxI7fffjsAzz77LC+99BKrV6+msrKSDRs2sHfvXlauXMmxY8cYGBhgfHycX/ziF9x2221cccUV\ntLS0cPDgQRYsWMAbb7zBkSNHiMViNDQ0EA6Hqa+vJxQKTYoSiVcj063Z3AeBqjI4OJhwVnJVtdmz\nZ3vhYfFGI4CKigpWrFiRcA2bbNZzBa+1tZXR0dGEM148QVR9//EdO3Zw7tw57++Ojg4ikQj79+8n\nFApx6tQpLr/8cm+9ma5fN+SuubmZaDQ66Zq/853vsG3bNnp6eryInOHhYe8eOzs7CyV4uzGznGdQ\nCVoQt6QET1UZGxtjeHiYw4cPMzAwwKxZs6iqqiIWizE+Pk4sFmNoaIjh4WHPOTw6OkpZWRmhUIiK\nigrq6uq8dsrLy4nFYtTX10+KbVywYAFNTU10dnZSW1vLkSNHPIuhq5q6azgwpnq/scZdF7qGDlcI\n3AEPTBj8/nWkK9CJZhVXONrb2zl16hSDg4OeuhmNRnnvvffo7+9POXhTCUummQeZVgrbuXOnF/9Z\nV1fHZZdd5n0uLkFmWrffZBpPUAqlaorIt1X1805a0ApMReqvqurn051bcn68yspKwuEwc+bMSanu\nRKNR3njjDYaHh1m4cCHt7e10dnayYMECBgcHqa+v58033/TM7gMDA1RXV3vGBTCOYHfddu7cOWpq\nahgZGfHWTP4ofZeBgQFPgN2ZE2B8fJxdu3ZxzTXXABNTemKxmDfgR0ZGKC8vp7a2Num9NTU1cerU\nKVpbWykvL2dwcNBb94RCIbq7uz0HfjKSReGEQiHKy8vp6upKmXmwfft2r/RGpuFjVVVVnuANDQ15\nD6pYLMbY2BgjIyMsXrx4Qh8w2Rfn79dVsxNFu0xhZoIAqOpDIrISU0G6PvUphqIInm+HocCbQqQL\n3/rNb37jqWH9/f1UVlbS19fHuXPnGBgYoK6uzhOKWCyGiHgzo6sGDQwM0N7eTllZmWe86OzspK2t\njV/96lfU1NTQ2trKqVOnGBkZ8YQoHA5767gzZ85w9uxZAG8m8d+DO8hPnTrlmd5jsRjd3d2e6goT\nB48/njFRaFxnZyciknUJ9yCC5BeKTEPp3O8F8IKk3Xt1AxXGx8e9Ptxjl1xyyYSKa/H9JnOSZ5os\nm814TMKT7gtV3QuscmoTpaVYaUH3YqK/n8UUDs0Z16pXVlbmDfyOjg4v3SRZoLB/7VFWVsb8+fNZ\nvnw5oVCIEydOMGfOHE/NPXbsGKOjowwPD09QN/v6+ryo/ebmZq+9mpqaSekpkUiEmpoa+vr6vAHp\nCsrY2Jh3rv88f8m/ROXxjh075gUGJEuDSZWKVOiSe66aHQ6HWbZs2YR+3nnnHUZGRnjnnXcmHdu/\nf3/K60qWyuT/vAKSl/GYSGDdAOx0FEvVzGpTiFTMmTOHwcFBmpubJ5jTo9EolZWVRCIRz5AiIuzb\nt4958+Z5zmQR8XxhMDH5s729HTA5Y62trd4s2tDQwKWXXuqpTmDWSbFYjOHhYV5++WU++clPAudn\nE9ctEA6Hqa2tZWhoiPLycqLRKCMjI17/8WpVqif3+Ph42hLuuaZNVVRUZB3EfPHFFye9/nA4THd3\nN+Hw+X1wGhsbOXbs2KR4zniSpTJlUc4i7+MxY1S14D+YRME65/XWJO9R+2N/CvWT6Xgs9E+xZry0\nm0IU0sFpscSR6yYlOVMUd4JvMZvVphAWSz4phfFYMn48i+VCYtoWO5pKRGRySP40RURKzpd7ITCj\nBE9ELhWRR0VkfoHaL3NKtz1RjI05JUjxlNza/wjwfRH5nQL3U9DvZTpSNMETkToR+ZKzo0q+264Q\nkRrMNkkngT/Ic/siInMxC/LtmAiF38tnH3H9lYvI75F5taug7c919nTrBp4A/pOI5K/y0Pl+Cvq9\n+Pop2NgqFAUXPGeWKAP+A3AcuFlEJtclyK7tchF5APhr4H5V/Tfg28BcEbk2T300Az8Bfhd4TFUP\nAM8DTSJyTT76iOuvDDNQx4EPiMiqPLYdciLq/wz4OmaNvx0jgHl7kBTje3H6KdjYKjQFFTwR+Sjw\n95hNSo6qqax7hPx9ya1Amap+CWgWkbtU9QzwAvBx50vJlYuBb6nqN4ExEflDVe0DngM+ISLteegD\nEVkkIpswQrFUVX8I/D/MTqP5Yhbwuqr+FbADeMA5/kNgiYhcJSLleein4N9LEcZWQSn0jLcauBvY\nx/kdgp6e3HB/AAAHP0lEQVQDLhORfOxAcTlmZgD4FvBpEZmjqs9jNpF4Ig9PwEuAZc7rbwK3O6Xc\n3gV+B/jjPBlbbsDMQt/3tflvwEgeVah2YJmIlKnq9zFCsU5NsZ4O4L8BbXnopxjfS6HHVkEpmOA5\nC+lZmKiBHwMXidld9l1gFPiuiGRU+N81NviMDr8APiIiVzrtPg98ymn3OeBPVPV0jrfyL8DVYvZ+\n/7Xz93LMJp1/rap/qapnc+lARELAfGBEVU8Au4EHVHUEeBX4qoh8MsM2JxlmVPUt5+XnnN//A7hU\nRBoxuzndqaq/yeL63e/FHU8F/V4KMbaKTd4EL/6Ldr7AWs4vqL8M/I6IzAZ6MF9yRtsyq6qKyEVO\nuzjqy4+AW0TkQ0AL8IqqdqrqU6o6lqK5RPdQ7vz+qNOP28f3gHUi8l+BNZj93N/W88mPGVkfRaRV\nRG4RkUqnjygmUv4+5+8vAZeLSBhTueouVf3XTPpQx0ErIpeLKU3g8jhG2P4AuB94U1X7VfX/qOqZ\nTO7Ffa/ve6lz/s7391LwsVVscnagi5wvluJ8QLXu08x56vwFRn26GviNBih95mu7TFXP+f5ej/mw\nH8JsGRZzji8DPgb8XFV353D9Nao64lgU31HVX/reVw+sx0Q6vJOkuaB9rsYkTu5R1T2+43+HWQfV\nAEtU9W9yuJcGjOp4KXAaMxjd/9UBHwLKVfUnWVz/pAI5IvK7wLiqbvYdy/p7SXA/eR1bU03Wghf/\n4Tuq0CcwC9yv+D6wJRi9XoEt8V9YkrbnAuccIwYi8gGMLn83cLmq/pdE15DBtZdh7t0V3CrMAP1T\n4H8DS4B/UdUB570ad6851SIUkSsxO+uWAz/z3Webc/wjwN85qmbQNsWZeQSoxDwkrlfVPxaRv8II\n2V/neu1xfXZgXB5fBq4F3lbVwyJS7n62WbZbsLFVKmStasZ9MDcAvw88DZwC/tw5Xqaq+1R1k6pu\nDih0VwH/hDENV4jI14C/wpTZ7gdeEZEPx19DUJyF9ybgk87ff4pRJRuBLcBCzCacbh/n4vvJpF8R\nWeHMcLiqJcY18UngDODfJK5bVV9S1a8FEToRmSemKrKr7q0CngLuxMxye0TkauAoEBKRiiw/s5CI\nXONevxjr57eA65y3rME8FG9yriVroXPvxdd33sZWKZGR4MUNogYR+aqIXO/821VZvgfcKSJNfjUx\nTbueDu+oXm9jBuTVGF/Qn2GsbfOca66X7E3SZzEL82udmXQ35ol9BTAbkzLyFmbWyAf/EWNJXAbc\nISJrMEaNi5yfNsdYEFigxfAV4B7gdd+/rsGsfUaAn2ME7vOYh9a3VXU8vq00/bifcQS4A2PhBTOb\n/jPGd/YhVf0pMIgR7qzcEYUaW6VKpoN3LWYQtWLUoS7MIP4gZo+wKzFP2k+osc6lxVkQHxaRj4ux\nroHxXz0PfBQjKHdhVI06YJuq/iQDoa4Qkf8lIs3OE/994B8wwvUxVd2JUVfc7abrMMKXq6XSfZj8\nK9CEEbIhoF1V+4E9wMeBncDhDJsvA97HrG+GxDj5AZYCP8Z8dh/ECMMvgH9Q1aOJGkpx/Q3AgIg8\n4hz6Z+BGR73rUNUXnfXhQWcm/HNVfSKH2W4teR5bpUwgwfMNon/BWKguwWzm9w8Y3f4SzOD6dyCq\nqscStZMIVR3GRDV8lvNp+KsxT9Me4A8xKtn/BP57Jm077Y871/w14DbnXj7oXGtIjCN2MUaN/XtV\nHcQM4A7JIYDYnb0c69phzGc0DMwSkeeAX2GslV1ZqEmC8ZP9EPjPwD85BqHnMQ7/HwDzVPVlzABu\ndoxDmXAaY5yZj1HvQpj15xmgX0QeEpFfAP+sqmPO55YxhRxbpUzGxhUR+Qzmg/i/mFCdv8CoAN9X\n1dFU56Zocx5wO2Zd9RWMKvMWZjF9KfBjNSb3rBCRyzHrng9gBtHvYwb+IHA9sElV3/C9P5e90hL1\nX4Yx438Xc19N+RhAct4K24rZNOMv1bgE6l1BcAxHY1mu7ZYBv4X5HhoxD8admJ12RsF7cOaFQoyt\nUiWbJ/oW4BsYNeAQsFFV0288nZoejKM4ApzACMPPHXP+L1OdGJADTjs1mOtegQl0fgXY6TNRlznG\nlLwJncM1OKqSo4rl66l9VkSWA38CRDEzuKjqoO9echmw+zEW3rCqfk9E3gf+Bvh6PgXORyHGVkmS\njYHiWszAHVbV1/LxwTjrtb2YB0EY+Kyq/nuu7fraj2KMKGBUyueBE6oaVdXTrrpTiAW7mKiUKHCj\nmgDrfFKLWZtuU9U/UdVBn4qb872ocXTvBhaIyBpV/TdVXVnANVbex1apkpGq6QyiqzBO5IwsZAHb\nvxhjUs/JHJ2kbcE4c3+jqm/nu/2pwu/zyqePztd+GcbFcrCQlsRCj61Sw5Z+mAEUQuAshcUKnsUy\nBcyo0g8Wy3TBCp7FMgVYwbNYpgAreBbLFGAFz2KZAqzgWSxTgBU8i2UK+P8SRarKiJfMzgAAAABJ\nRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fa218f03f98>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "##version with data\n",
    "df=pd.read_excel('/ocean/eolson/MEOPAR/obs/MackasZoop/1-s2.0-S007966111300061X-mmc2.xls',skiprows=1)\n",
    "df=df.dropna(0,how='any', subset=['Year','Month','Day']).dropna(1,how='any',thresh=200)\n",
    "yd=[(dt.datetime(int(r.Year),int(r.Month),int(r.Day))-dt.datetime(int(r.Year)-1,12,31)).days for i,r in df.iterrows()]\n",
    "df=df.assign(yd=yd).sort_values(by='yd')\n",
    "dtsdf=[dt.datetime(2014,12,31)+dt.timedelta(days=int(ii)) for ii in df['yd'].values]\n",
    "df['dts']=dtsdf\n",
    "\n",
    "fig,ax=plt.subplots(1,1,figsize=(3,2))\n",
    "\n",
    "p2,=ax.plot(dts,MesZoBar_old*maxRate*24*3600,'--',color='k',label='this study')\n",
    "p1,=ax.plot(dts,MesZoBar*maxRate*24*3600,'-',color='k',label='Moore-Maley et al. (2016)')\n",
    "ax.set_ylabel('$G_{max}^{MESZ}$ ($\\muup$M N d$^{-1}$)')\n",
    "fig.autofmt_xdate(bottom=0.3, rotation=30, ha='right')\n",
    "yearsFmt = mdates.DateFormatter('%b %d')\n",
    "ax.set_ylim(0,10)\n",
    "ax.xaxis.set_major_formatter(yearsFmt)\n",
    "ax2=ax.twinx()\n",
    "p3,=ax2.plot(df.loc[:,['dts']],df.loc[:,['Total Biomass']],'.',color='silver', MarkerEdgeWidth=0,MarkerSize=4,label='Mackas et al. (2013)')\n",
    "ax2.set_ylabel('Zooplankton Biomass (g m$^{-2}$)')\n",
    "ax2.set_ylim()\n",
    "ax.set_position(pos=(0.15,.16,.7,.8))\n",
    "ax2.set_position(pos=(0.15,.16,.7,.8))\n",
    "#fig.legend((p1,p2,p3),('this study','Moore-Maley et al. (2016)','Mackas et al. (2013)'),loc=2)\n",
    "#fig.savefig('/data/eolson/results/MEOPAR/biomodelevalpaper/figsMod/mesozooWithMackas.eps',dpi=400,transparent=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANsAAACoCAYAAABpNZsvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXmcFNXV978HRQWVRUSDS2QLedwSQVBfNfbgING4JAJq\nxDWIoMZPHqKCms2JxOCSvImKETdcoj4CKiYRFVChjRoV1DxqRFGQV+MCso0bRNTz/nFuTdf09FLV\n3dXTDPf3+dSnu6ur7r1VdU+d/VxRVTw8PJJHu9YegIfHpgJPbB4eVYInNg+PKsETm4dHleCJzcOj\nSti8tQfg4VENiEhnYCAwQFWvcvvOBJYAvYEZwBj3+y1VfbHSY0iUs4lILxFZKCLXi0jPJPvy8CgE\nVW0ElgLbAYjIcGCJqj6uqjdjhHaDqt4PjE1iDNXgbIeq6kdV6MfDIw4OAxaKSBegERgYcDyM01Uc\n1SC2E0REgeeTYM0eHmVgoar+U0TmAGtEpFOSjCFRYlPVt4CbANwFDU2yPw+PCBD3+WbW/gWYiPkR\nsCaJjhMlNqeATnNvixYX4Dieh0ciUFXJ2jUE6O/sBzcBY0SkN3A5sBAYKyJrgElJDSixDegJ1AOj\ngZ45/tck+3d9NLSlPubPnz9+/vz5On/+/PEb+7Uk3Icm3UfcLWkxchmwLMk+NkFMzfr02Ejg/Wwb\nGVKp1CrgqqIHetQcNoUIkvm+j5rrpxp91BzEybet07mIaksl1sOjbGTPrewIktDv3sBq4FE25ggS\nD49agWZFkGCEppgVshdViCDxxOaxSUJVHwPeAsao6u+AQZpxaCcSQeKJzWOThVrQxRUiMgNQEemU\nZH/eGumxqUEARORyYIqqLhORXsA9RIwgEZFDgT7AAEzHWwo8qkVCvbyBxKNNQETqgLrQrkuy55aL\naBqB6WRdgG4YwTyPiyDBCG2hqv4zRx/9gYuBOe6cQAfsDRznzrs57xg9sXm0RSQxt0TkEFV9wn1v\nEbQsIv0LWTG9zubhEQEiMh7Y3+VmdsIF2IdRzF3gdTaPiiGdTncDRgFTXaRLW8ILzoIZJJ52jduA\n52welcQo4Er32dawVESGAajqfcANcRvwnM2jkqjZIOk8ESRNESOYsSNvBIlzE7wV+n1f7DF4A4lH\nkmgt0TLX3HIm/jGqerHTwW5Q1Y9EZAqWTHpj8FtVzyrQ9qGq+njcMXkx0iNp1KpomR0xEieCpCQG\n4cVIj6RRy6JlM/N9jBokJYmDntg8EkUN5t8FXOk5MhEja4lXg8RzNo9NFzkiSHIhVw2SNcBvMQNJ\npBokgQsg9hirYSAJK6NZ+72BxCMRuGJS7YATgTmqujKhfjoBRBE/E+dsIZPrDIxFe3hUC7OBHbA4\nxooSW1BrEhMpFTih2DnVECMHYvKwh0fiEJEtgPPdzznAH1X1iwS6ekFVL4pzQtJ1I/tj0dQDkuzH\nwwNARHYHpgPvAL1c+sxgEVmgqp9UuLsuIrIAK6mAqn632AmRia3EHJ5emIVnEMZyW6QfiEhD6Od8\nVZ0fdUweHllYCUzEVJZOInITVoX7KODlCvfVVVUHxTmhqIGk7Bwe09mmAzOyj/MGEo8kICJHYLGL\nDwPjA4bg5nJXTM9aSBkFfkRkOuY+WArgapcUPicCsfVycWH5/i+Yw1OkbU9sHhWDiHQF/gCksArc\n4epaw4GlqvqiSyLtQihcq1B4Vp6+hod/R4mVLCpGFiI0979fmcajVrAL8DGwt/sM41HgMadnXQTc\npGUsEVXxQGTnQ9gutGtIIZExduees3kkhBx1IwMR8jBMjOwFnOk42xxVTXyFpWKcbRBOL3MD7JP0\ngDw8SkGECJITnKn+cVfsJ/ElorIRRWfr7ApcVr5zz9k8EkIOzrYPGWaxhggFfio+Jp/P5tEWkfTc\nEpELgHudL68euBALC/tdvnNi+dlKSZjz8Gij2A9ARNYCh6nqUCee5kWc5FHPgTw8MggW4xiIRZMU\nraYcJ1zLL8nr4ZHBFVjKzlz3+2asqnJexCE2z9k2QrTx8nJVhXNkz1XVj5z/OVw78j7nLM+LyMRW\nasKcR6sjqAECtZUxvdFBVe8TkctdCOISLMO7C9AX8+FdWOj8SMTm4sACMTIoTqlRIp09Wh01WwOk\nmnBc6WIg4O7HU0JsZJBW4yp19QZWF7JANhtDHNO/iExS1Yuzv5cKb/r3SAq5/Gyq+k9nyOiNRZKU\nHBtZCuKWsjtMRLZ1VH1YEgPy8EgCIaf1Ye574osfZiNu8uiZmNyv7rtHlZDL0OGNH/HgdK2uod9R\nS9dVBLGIzcm1ibNbj5zIZejwxg+HGNW1lrrv4VJ2sWMjRWQSVoXrceA5VT2n2Dm+lN3Gg1yGjo3K\n+JEkJ3YZ/vOD3yJySY7DumIGEWheyq5g6bo86IoZXM50n0XhYyM9qoZ0Oj0e48QTUqlUopy4CrGR\nZwJ9VPUiEbk8SvGfqKb/nrn2q+qyOAPc1FHum70N6GgbFScugt7AJBFZSMTqcVHFyLCephgLPRPY\nrNBJoZqRgT8idnZrG0O5OtZGraPVYCnyctAVy/g+030WRSRiy3LkXYgpln0jnDoQI86FQH2Uvto4\nyn2ztyXOsLHjeaC3q2lSsHRIgEg6mxMjg1CUK+KIj45AJ6jq2Tn+8zqbRyKogs7WGYtC6QwQJYok\nqlN7aehzhIhc4JLnisIFbF7pQr48PNoKbsTooS8Ry4VE1dn2LWU0Lpluiqq+5ThcrmMaQj99kVaP\nxOAsiEswG8IMyqgbCaxR1cdE5DCseFDx/pM0/bu6D92wKsrPZ2d6ezHSIynkiI0cjhHI4+53s2V+\nS6gb2SvERAZEMf55P5tHm0QOYpuCGepWA43Y2tonuP9KLmVXU0tGeXjUEBa6yP85wJpyYiNrdcko\nD4/EESE28s2s3+XWjYy9ZFRU03927JhgyaM+n82jJpFDjOxMxiCyljLrRjqxdF9iLBkVldimYG+B\ne4AXcPVIiq0DEKFdT2weiaAKfrZpgc4X+ZyYmdrDsaTRNZglZ1msEbZszxNbAbSBWMhWQxWILfaS\nUbEytZ1580KMdc4oYYwe8RDEQo5q7YF4tMA0MsEekRCnIvIwjKutxbhaWwkorWX4WMjaxVqM8XQl\nYj5cVJ1tNWa9CXvZvYHEo2ZRTZ1NRGZHMZCUE67lKyR7xEIb00FFRLbFVLFIRB01xeYt1/qhblcX\nTKRsEcnv4VEArZaP58KqZmAS2hWYka+c2MgLgVswpjMmyglxndp9MELrRqbGuUdCaGOcAFpfBz1U\nM4vZN4uNJGIhKxEZ5iyP/TEjCVjs77Ji58YltoDABBhe6ECP5iiRcDbqzGzIed2teR0niIhiiZ+D\ntLQ1tSXrM/t7XsQltoD6FVjqUhY6Ry2/vImjFMJpbU5QCdTECyO8EEY5sZFBdH/wGaeNosSWtXJH\ni1iwYit3eDQhNuHUACeoBKrywigWG+nm6TRHGGsoMTbStTMCC/W6ESPauap6c9FzI5r+L8fSv3Ou\n3FFqJIk3/XskhRyxkT0xm0MvbBHD1ZQQGxmscSEibwD7Op1vuqoeX+zcUgr+xFq5w8OjFuAYwrKs\n3eVIDY+FxMclBY908MmjHm0SSc0tVydyLub6mosZR+pVdVDRcz2xtUkT+yaPBImtf679Ufx0PnnU\nUBMWM4/aRwnO7yYUJDYRmQ2MUNWPHUUviWMqDVVE7g+8qLW7VHBbMLF71DgKipHO7D8E87T/M27C\nnDOTzlXVZbmKqtSKGOnR9lCLc6uYGNlLVc92i3b3JhOeEgmqGjgR+xPRYuORH163LB9BmBZm2Cgn\nNjI2ihHbW24t4otCzryiGak5cDyZ8uUeMREiso5Ag9vtdcuYCKk1M4DjKCE2MtTWaGAo5nPWslNs\nwoUnVfUmZ/aMBSeKXo4FL7fQ93xF5EgIDDgNbuuYTqe7ee4WGwOxyBGh9NjIAIdFcWSHUcxAciiw\nNIgQictqQ4S2FCsU1CLZVFUb4rS5iSJswAkI7zM8d4sMp8osxCL0AbTMNbXFrXcRuQZJMTGyizNu\njHa/l2pWCfFCcJxxU1+TrWyEYyTT6bS3nOZAhLqRvbBYyEFYyGG5dSOnETOBupg1cpiq3h/UQheR\n4VFqmkfuvAYtRh6VQWsbc3LNLaezTcd0tumUVzeynlANkkpwtrccq7wiR2cBIe4Td6AemwRqLlBA\nVRuBsCGjnHGNCVxZzh9dHrE5He1F12AvjP0G6ON0ui6AJzaPbLR1cTeoQSJETB4tJkY2M5Bk/VeP\npSnsGyWXJ0/7m6wY2dpiVltHFapr9cIkPiVimlkxMXJfzaxn1R+Lbg5Sa5aqrU9Vitl0o0aFCCWn\nmFVrRFhr46kh9A9M/86AWJThFCO2513S3VjMhD8w9F9w4zfFyJBK6CP5xKyS206IMGpO92ptuGTq\nehEZhImQA6gAsV2EWW2WqGqjiIQJa4yICGY+3dR0trL1kQIlD3K2HZGQkiCMtq57lYKgAvKNEH2B\nmWLENtZ9DnEhLb2AwM/WiPkaiibNtTUkWRukQNtRCKnihNFG6qCEQ7V6Y7aGRykxNtIxngVOjeos\nIteratEaqnFXsenszKfh9a5arJUdo71N1kASF153ioccNUjqMWPGGqAeE//KWVP7TEx83A44M0ok\nSjFr5PVuQImIieUSW2tNQNfvue7n5Fqf/JsioeZxavcCJrhMlukhA0fkNbVdCKJixHoCFq41JwrD\nKSZG3ggMEpEhwPZYyk2wmMBwzHtOqab/CqC1lPdzyUTfF41RLGeyJ2n5rDZam+id2HeFiMwAviox\nNjLwNSswJc6JkZza7o2wGmOZAdZgBU/6xOmwwmht5X1OxL5bTPYoE88dcztwJJBKp9OnlThJW/s+\nBUiM6CPUjbwcmOJifXthq+jGjo0Mqg04MfIwXIoNGVtGXsSpQXI8Vjsy8LMNxWTWF2K0UVG0ovI+\nGeNoRd/QIrJrv379Phg5cuQDt9xyyxbvvPPOWGDZ3Xffndppp52CLIh81zAKI7RF7nNUgWPzooaM\nHIkRvUvNmh/8FpFLsg65B4t6Og6YgFtTW0TWEHF9tSwMiZtiU0xnm4NZa+YCL4S95CJyqKo+LiI9\nfZHW5hCR7TDCOAMTv58F3gHWA9sAPTt37rzn0Ucfvc0jjzzy+MqVK+cCjz/88MMrOnToMIrMZAz0\nwruAY9kIda5SRcdyRc4qRJDEXua3GGc7DmO1vbG3wljNLIA4VEQOw1IVlpU66LYEEdkcI5CfAY9g\nk+VZVf0q+9h0Ot3t/fffP//pp5/+98qVKwcBE+66664uo0eP7vTkk08O7Nu375tf+9rXfoa9hUtJ\nAakVlCo61oSeWQBBiRAhYqpNMZ2tEfOnBeuzhTO1p6nqi84FsMnDha3dhYmXdar6apFTRvXo0ePi\nqVOnTkilUj8SEVm+fPm+TzzxxEXXXHPNFhs2bBj7/e9/f/VDDz1Ud8YZZ3z3iCOOqH/ppZf61tXV\nXQO8D6zNRcQ1iJlACnghnU4/CIxLpVJvRjivVvTMfJiLWSMjz/+otf4TSaOpZTEyjhjjrLV3Y7L/\n1VGIIJ1O9wX+SJ7JJyLtgH7AXjvttNOgYcOGDX3ooYfWH3LIIb1nzpy5VWNj49bYi3C121aFvoe3\nt4F/qeqqYteUhEsjnU6PxzjUImB3YFYqlTqq3HaLoQpi5DTMWn8cVoOkqFM7qoFkkIgMUdXfhVLB\nHy1mNg157QfoxrfgfSQxRkTGAJdi9TWfiNH+sZjRI+2yr5sRgSPY19x2L3DhiBEjLgEaTj/99Ia6\nurrLMEvYdnm2b2B1X3oBe4jIujFjxnw8cuTIvk8//XT/urq6icBiVf0y65ob3PeSyy6EiZoMZ5rn\n2h5XSps1iDWq+phTpXpFOSEqsS0NFVhtxMz9g8hRUyQMF9ayFKs9ubGhqBjjXjxnA99R1Tey/y/C\nSXLVFYHCE7xD8KmqXwAr3VYQLoZ1lyeffPKAnXfeeUxjY+O3unTpMmvt2rU7iMhLmEX5xQMOOGDB\npZdeOnGLLbb4khLFtyx3BalU6qrQy6RU10Wx/lrDdxckVN9Apq5JQUQltia5NFQLsj7W0CqEat3c\nfObyoP+TTz55W+BETD97J08zeYmomnVF1HSFd9Lp9GfAacCexxxzzIS6urobgX2wyVL3zDPP/HTE\niBG9R4wYsfrBBx/su2LFimcwK2pA1J8AnwMbQtsXrv1sQpsVup4kjR1RJZBm1bkx03/s2EgRGRb6\nHtT9L99AEkK3HCb+SKwzAVT8wcUk4FHAlfX19R/efvvt/VX13QLHRiKiGH6wrdxnF6cLxX3hBH67\nWcBUZwBLuw2ARx555GcdOnS4bLfddtv6kksu+TZwOOa+6I65LdpnbSIi64F1p5xyyuZnnHFG55df\nfvnDiRMnfrhixYrjgOduvPHGO/r16wfJvEyivqiOx6pzP+ZcWnMorW5koAfGXpEm6vpsN4nIdJdi\nE6yr3TVGP5VUVAPr1swKthmZgEeOHNlt6NChK9atWze4CKHFdiZHIPr93ecRwNejjDcLTRMzH5F2\n6NDhBmBDKpWaqqpFCdm5O7YCtvroo496LFmy5EfXXXfdghUrVnQD9gPOHT9+/M6nnnrqW7Nnz16x\nePHi/1HVz2OMuSCi3mNtWZ17UCgROnICtGaW972cTCDyPVHOjRxBoqrHO9FxKPBmDIPHEKB/Pud3\nCUVamwwLVC7DOdLb0UUlHHPbbbd9S1WXx+wjCooR/fbu8wvM/xaLU0SZmHFfEE53/MRtK4HzRo0a\n1eyYe++9d2L37t1/sX79+l8uXrz4chG5GbheVd+LM/4KIajOfVOZdSP7ZxX8Kbo4aKwlo5yRJNZK\nNO6NclOB/xvitEcCGc5RJpiI/BLzqwwuh9BiGE1y4SzMxTAmlUrVzIpAxV503bt3/yPw0UknnTTz\n4IMPHjthwoQey5cvf1lEbgcuV9UV5Y4hQt3IIHh+EmalLbduZGOonqoE3wsF5W9067MVIIyC4mU5\nnE9EfgGMJCah5emzWFByIaI/GNjBfVaF2CL63gq+6IJnlk6nx++2224XTJs2bUJdXd1emDV7kYhc\nB1yhqp+WOs5isZGSqc49BrO+TqK82MhwkdYbopyw0RFbAeQUL0MoyvmyicOZzH+NOS4Hq+oHMceU\nq89c3CsqV+6a9dkEEekHHAV0wlwEm7ntf1X1thzHfwMj2pWhbdXDDz8sQXymI5IovreoRoqm45w+\n+BMR+R1GBK+JyATgnsC6GQVRX6Kauzp3OQa2tVS4SGuroQROVOyBR5kQTZNeRP4A/Ano39DQcFpd\nXd3N6XS6WbRHhDG26DMPZ25xnIh0BQ7AfJp9gd7XXnttau+99wb4Vg5r5LbArmTEoi/dlk9E6wIc\ngumBTdtDDz20ZPjw4fu6Y65yY+r4ySefbHXfffc9kEqlWjSUfU357kuua1fVt4GRInIwcDXwYxH5\nyfz58/9f0IY7NN99bq0YysoWaW1lRLqJuUSwXA+7gN+sKWwK92BnzZo1HYva2Bo4tK6u7n9wTtp0\nOn0amUkQSXzKN3ZnydvRWTWzj+sHnAe8iVnP0n//+99v7dev39lbbrnlC9n9qurzwPMR7lEQobIA\n+FH2sffcc0934PTgXrjjfy0iJwMLbr311i+B1932GhZJlN1v030JObRnUiBzQVWfFJH93Jgeuv/+\n+98ZNmxYuJpbvvvcWjGUIjGLtNYysUW9ibkebKR1zNwk/CsWs0cqlTpKRP4KPAw8NXny5LP22muv\ns0JtjaM5gRUcY3iS19XVfYrplHsD33Kf3wSewHxZzaCqz2LJieF25mFO5XuxWMioEyzy2/+HP/zh\nh7mOUdU7ReQuTGf8ZmjrTUsin7py5cptn3322elHHnlkkNU+BLNk5x2DCx27WUSmz5gx49erVq3a\nY968edufccYZk+vr6yHH9bZirt6FwC2Y3ja2yLFADRNbjJsYNoyE1zGLYhofhRHaog0bNoxzcY6X\nAT9T1ZvS6fTVwE+Aa4LgWUfQHd1GtkFDRDqo6rpQ+2HCvBB4CdMrrwVejWgUCNr5EHMut0+lUkUX\n3wuhhfGoFIOR06WWuy1vHKjTd78EXl2/fv3nw4cP57nnnuusqne+/PLLf80limb18xHwUxG5Gpg0\nceLEf0ycOPFS4JM4BaqShFr5usqtz7aRIGwYKeq0zcJMIPW3v/3tjo8//vihrl27frJmzZrBqvqK\nm4wnuuOaOE8qlVrlwp6uXLZs2TZ1dXULMC4VbN1FpJuqbgiPxxkEDi3xGqdixHKk+/1azPNPcucu\nwAw+kJCuExDx/PnzJ9fV1f3moIMOmgKM7tGjx47jxo1bvWrVqp3+/Oc/v559/OzZs6dNmjTpnbBx\nxPllTxSRA7AX6K+cQeVmVf0kzriyg+JD1eFKWuZXEljFZmNAyTL7smXLxvbs2fPInXfeeciAAQO2\nPP744y+tr69/xf09CuiuqivWr19/Sq4+zzvvvOMwA8bLWC7by5jD/wsozJ2zIuOD/nK+JByBn4Yt\ndTQYWJd9TAlocd8qVJioSYRX1avS6fS/AXbdddee9913359SqdS8rFNHAVeuWrXqTGA7EXkZ4/6v\nAG9gybfPAIeLyL6Yu+BXInIvVoV4QRTrZY6g+DGUscwvbcxAEglZAb1B7lRHx32aTZp0Ot1t5cqV\n/33BBResXrZs2TGdO3fe4+yzz55zyCGHvAqMe/3113sPGTLkPKB3t27d9jzyyCPXPPDAA50HDhwo\nHTp0aLL+BX3ef//9LVJjwigyecOcBVrqnS0seel0+iWM2LYiBwr011QzJdd9CyHQrzqS4YA5kaOv\nfCL85NBpuQh7HjDnxBNPfHHt2rV3TJ8+fRdMpz0QOBU4BVcJwBliRojITsAvsZfPVyLyN+BB4KkY\nvrpyl/nd+AwkInIZlh/3BqaTrAVWquqGXJMn34Ry+8OGkStfffXVnnV1da9hetleJ5988rdHjx7d\n6fDDD39typQplzQ2Nj4wadKkz92571122WU7YjGHb65atWrWHXfc8Trw9q9+9aufktsR3ZRKQgTf\nXp48rzBHLmbh/FbWZ8H+gp0hp3K3EgOYm5CLg5FxEeyMOf9PyrIC5yLcYKyzsBIbQ88555xV06ZN\nuwoLEs7V59RUKrVKVd9zBPc55u4YjdV72VJEFgHPRLmWMsO1wgaSMVFOaHViw27UttjbWgG6dOny\n9pgxY77q2LHjV8cee2zf6dOnn1JXV/cvYNDIkSN3HTNmzBY33HDDZXV1dV8C7YDnxo4du8OJJ57Y\n77bbblsxc+bMr773ve99MWvWrNOxOMK1wL9nz549q1evXtssX778Z6oaiItNk/Hdd3PHFedJgWkW\nQZ/n2rLPa+IcqVTq1zQnpoAYdwYWY2/8bAScbXEeoim2TkA3bJLk41wtOGAOFDJCjcQMOPPT6fQ+\nqVTqzQLcNjhvJqZLthh3COPduLth60+gqt+Hply97YBdgJ5YyYhCJfEDLvQcJYRriRUunhA2kLi8\nxsrGRiYBVd0RmsoAbAt0+f3vfz+uT58+41555ZU7Fy5c+OTcuXPT2FvsH0888cR2ffr0GfDUU089\nBXzUs2fP9mefffagxx577NZFixYNWLp06d2HHnroR6effvrRY8eOvdGJX8EDnwkcW19f/z60fGPG\n4aREMMYU0NkOTKfT3fKIlv/tvl+WTqcfzWo/yCvcl4y5uQUHy9FfQCBzcvwXZbxh5LvuURihgfkn\nr8NW+Qz6bibaZ/WVU2QN3fv93K7+2cc4fW2ViOzt/m9xTAhNQfFYvO6YEsK1BgGPi0jYKNI3yomt\nTmwBXBmARqAxnU7/Bnhvr732mppKpVadf/75LY6/9NJLAUin05cAJ++///47ACPT6fQq99b/LfCF\n40pNhU5pLvZli125xLCColnMy5yMPax8NSBnAkcDO2EWx2wR9Dj3uzd5XBt5Xg6B6b8BW1AiqdV3\nziaT59jejSXoqyPxrZ/Bvb8C47h5SyoUi410x2QHxZdihV2gVr78chEZi6tBEuXEmiG2MIpN5CxO\ndaDbPRQ41709Z5LxhY3HJvc8zMK1gJbi1tQsnS9bj4KWolk48mQN0TnkOJq7KcI4CfiO+76a5gQ1\nCuMYADMKBCznejk0uUeKBDqXDHfdyzBiW4WJu0EAc1dMz7yCeIQe15VTDQiAql4kIgOwSsiRKmzV\nHLHFXIcshRHZHOBp918w0T6juejUE5sA84C70ul0QKzhNhuwCR5wTLBI91wT9I/YBO6NRYA3uP0t\nOKS7prvdWAN9rRDmkCnMGiBwpkNzC1828hlesvcVRIlugFuxl9/lWLpLB0zXCvBZlLZiZEG0Bpoi\n/FX1BWCgWCnyoqg5YqOIszWLA92F4xIhjvIZRlCXYW/SmzH9IQgVGkyGUH6EWSqDSdxAxiLY4PYd\nmE6nX8R8W5NdP33dOW+58zvQnAtlf44K9V8IYQNFMzN8AateM+SSCkq0Rpbi9P4dsCWmA22OSREN\nGGfbnTxiYAE3Qpy+q4Jczu8gC7wYapHYir2FA0JYBNyVFXwciC2XYZN7g7OIPe1+r3LnNbjjAmti\n8AZucMQUcJED3XkBoQQpJn/EiPZDt39d+A2cY8IHZvEjXH+thTiTuJRggQuA27B5tYjotSezx9Va\nwcWJouaILYLhYSb2UHcH/hiKwg/7fa7GxMYljguBcbvBWA7XYKzKVODzCog0PIZfOwIOojbmkXn4\n4zDxcXcc0WcPMkeEyFAskr8ByzvLhfCka2GGr0BlsciTuMT6KRMxV8w64NQYY2w2rigpO61Ywq5k\n1ByxRcCxWOT5YszgEYhbwdYBMxX3c9vh7nOZO/9JMg8o4IrBxJ6XdiWyyRg9nsOIbUMwAMctj8Gy\nA3bHDDPv0vzBhwmnozvuQzIvhFwTJtsgUHKERy4kHCE/CtjNfe+A3bOFkJdYsitCFxpXZCtxLaNd\nko2LrTc8XkSGhWrslYt5GDd5hIzyPYFMnX2wBw3w1osvvtjPfe/pPntkN+j8PldhE/lILO0mCP3a\nyvU3FHso3A9RAAAKIElEQVTAAY7F+Veeeuqpke7Y8P8zMUNHNzLFVbuHxgaZCTMqPI5cb2oRqUol\n4TL6mYoZD1ZiZbnD3PNc7DrvTqfT3Vwf12H3+roIbc/ExO+wQWsqMYoeJTQXYyFpzlZusGcuNGBc\n4h3czXZ6VkAcYcftskcffbRr//79u7jjd8UCh0eRFT7lJnhYPAwKrx6BccbsSJGpGCH1v+aaa7od\ndNBBV2X9fywZfe8KN65FmHO3WzYni3DdPwCGUyTCowLi1Q8wjhMLrq+zKPyMh2KENxj4X/c7SgZD\ni6yFErh0EnMxFpImtnKDPXOhwbX181QqFV5VJxz+MxqLJPhxY2PjfKwEwLsYsc3DTXhyRDdgfq7A\nLRBYLWeRVTrbfb8IYPny5fNzmKjDpnqwibUrRsifAVfFnTARj69F8WoyGWNTgNVZn0kjibkYC0kT\nm5YZ7JkLg7EJ26QTQIuJeFGwf9tttw1EuHUYJwwMKWEO0RTd4Igm0OWajChxuUTYVO8IO4gcKRRL\nWQmUVcS2Q4cOm5cbrJwNJ3mMJGMwGky0OMwAcY7NhyTmYswRqCa2YVbDnu77tBz/q9/8ltQWZy5W\nY4u0PlupCGXDrgEWagJrvHl4REEtzMVEic3DwyODRE3/Hh4eGXhiiwERyVmOIIF+El/62NWs9Kgi\nEic2EekkIj93y6Em1cc3RGSSiHw9ofbbuYUTprjEw0QgIpuJyKkk7AMSkQOBO0Xkewn2kegzcX0k\nPrcqicSIzU3QdsBeWJ3BH4jI1kVOi9vH5iLSAYvkWEmOCr9lti8isj2WPj8fy106tZJ9hPpqh13H\nF9ga2AOLnFJKH9uLSCfgPWAKcKSIdCxyWtw+En0mro/E51YSSITYROS7WDGUi4F31ZbR+TcVmqiO\nA4zH/FgXqOrDwPXA9q5mfCX62BELCTsJW9boTayC0w4iclAl+nD99BKRe7Dolf9S1bsxx/sJFeyj\nvaudMQ7zIYrLbH6PjeiZuH4SnVtJIinONgiL4lhMZpH7mcA3RWTXCrTfA2inqj8HdhSRM1T1M2A2\ncIR765WLXYHJqno18LmIjFLV1dh1HCUiu1SgD4B6jADuBM5xeuHDwLoKikdbY6vZ/AL4Oxb3CZbQ\n2k9E9hWRzcrsoxrPBJKfW4mh4sTmZPStMafiDKxC8BBVfQ34D3CriPQps5vdMXELLLrgOBHZVlUf\nxJZenVIBsaIvsKf7fjXwQ1cn8DXge2QIo2SISHusdN46tQUBFwLj1cqXvwhcISJHl9OHwy7AniLS\nTlXvxIjhMFchqjdWg3HnMvtI/JlUaW4lhwpEiUiOfVcDZ7jv3wb+L7b4+U+BHUrtI/TZEVtrax/3\n+1zgZCzIeAywRQWuqyNW5bZ/qI8DMQIcVmKbPbAA5S1C+04Drg79vhvTEUcFfZf7PNz+a4FR7vu+\nwDlYBvUpQMcynkm7pJ5JNeZWNbeSndoiIhrcETNVd1RXjda9Xc7HRKP9gLc1Qi30Iv11Bz4L9TEC\nq63/EFaR6nZVfaPEtjdT1S+dPvCCqn7o9n8fW6llAzZBz9fQyqPhexCxn0FYgPTzGlpmSUSuxcSt\nDkA/Vb2slOvI6mt3YL0a90JEemGTfjEWEPwnVf173GvJeu7dgQ2qutb9rsgzqfbcqhZiE1v2Q3Fi\nzlGYkvqb0E3qh4kPCkyPOSnbqZW2C34fg1m1LsJq6X/p9u+JJYemVXVhzsYiXIe4lWec2X2R2tpl\nwXGdgWOwEJ9FcfrI0ec+WPWpzYDH1XRARGRnt/9A4FrNrIJTyrV0wcTCbwCfAiND/3UC9gc2U9VH\nSu0jtO8k4AtVnRbaV5Fn4n5XfG61Ksphi5hyfy+2Osg47M0PTrQoob3tge1Cv/fAMhPOormolVNU\nitB+O2yiBb+3xMzHN2KrfJ4KdAkdK1nnR+4X42CD3Pct3GeQ8/Y9bNngsq4nfC5WYm1LrI7hn9y+\nXwCXlNtHVn+9sQyJbTCi2tXt36wS7Sc1t2phi2QgEZH+TgRCRLqIyBUiEmQcB2/I27HlWnfQEFeK\nCrEVSh7AfCabi8iVbrIMx4JHnxWR/wMunDt++7sC92DiDSLyEzfmrsB0LJN7KBD08VV2PzH7TWFG\niT3JLHv0Wyxbuzuws1P4Y1+PiOwktiA7qqrOJ3cjVv77U+B5sVU83wXai8jmJfTRXkQOEpEt3O99\nRWQymbqWB2BWwaFuHF/GaT/UT+Jzq1YQ1RpZh02cHpiosxSzCn0bq6O/D/aQj1KzqkWCk8cBUNNh\n/uXa3A/z14zDrGQ7ubF2LsOEvN6N+WAR2QOz/F2KFQ/dBqv9+AqwRYntA82u6W9YrZTuWD35XVR1\nDbZK5xFYnct3cjZSoG0R+Q1wJpbpHOAgzECxDivt9y5WnXg4cL26Jawi9hHc32B9uqBg0jHAXzAn\n8v6q+ihWwbp9mW6DOhKYW7WIghM3NHH+CnwNu/FzgDuwKlV9sQn1BKYovx+1YxHZBnhHRI4QW6wd\nzJn7IFYjfj22MslRQCfgMVV9JOqbzXHH60RkR/dm/9CN+xXgcFV9GpP767G3cyeM4NZHvYZcCDiI\nqi7BiKkv8AmwtYjMBF7FrGlLS+DQ7bCiQXcCH4k53gH+C6sC9iA2SRuBp4A71NbrjgSn760VkQa3\n6y/AEKcj9VbVuWq63jLH8c5T1SmlcLUk51atoiCxZU2ct7C3zG7Y5DwQe0vfoqqL404ctZUjr8dM\nz8EiEYOwN+cHmOn7M6w0wcS4N9u9zb+G6Rcj3MP9Nvbw2jvLYx9MRL1FVRuxSdtbKhekOw3TCVfg\nnNaqOl2jryGWDcF8WXcDPwYecEadBzEn/F3ATqr6D4xD7OgMPFHxKWZc+TpwHtAeM9x8BqwRkYtE\n5CngL6r6ubtnJSHJuVWriGyNdJP1D5iC3w74QFU/KKtzW2Prh5ie9BtMVHkFsz59A5ihtlxuqe3v\njukxe2CT5zSMszRiqfn3qOpLoeMrmjYvIt/BaimOVgv3qlS7gfW0B1Zh7Geq+pmIdA4IQES2BD4v\nQVfbEzgEewZdsRfh01hJ8f9A04uyYkhibtUi4rzBD8YsbJ+o6tIK9f8BFinRDXv7D8ZMxgvIrNlV\nDt507XTA9ID+WDDxs8DTmvHdtHMGkUoSWnvMPzckjs4UEevFlkg61/XR3pnNG0PX8p8S234Dqya2\nnareLiIfYhWmr6o0kYWQxNyqOUS1RgYTp76SN8PpXy9gRL8dcIqqPlHB9jeQKQq0BhO3VqjqBlX9\nNNAbkrBwuT6eSYDQwIwh9Zgee66qNobEsrKuRVU/x+7ZbiJygKo+rKoDkjJOJDW3ahE1URbBmeXf\nK9V8XKRtwfxBb6vqvyrdfmshy5EdK5IlQtvtMFfIso3Z1F5rqAli8ygNlSYyj2Thic3Do0rwNUg8\nPKoET2weHlWCJzYPjyrBE5uHR5Xgic3Do0rwxObhUSV4YvPwqBL+P81VkPFTh0xdAAAAAElFTkSu\nQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fa218dc7fd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "##version with data\n",
    "df=pd.read_excel('/ocean/eolson/MEOPAR/obs/MackasZoop/1-s2.0-S007966111300061X-mmc2.xls',skiprows=1)\n",
    "df=df.dropna(0,how='any', subset=['Year','Month','Day']).dropna(1,how='any',thresh=200)\n",
    "yd=[(dt.datetime(int(r.Year),int(r.Month),int(r.Day))-dt.datetime(int(r.Year)-1,12,31)).days for i,r in df.iterrows()]\n",
    "df=df.assign(yd=yd).sort_values(by='yd')\n",
    "dtsdf=[dt.datetime(2014,12,31)+dt.timedelta(days=int(ii)) for ii in df['yd'].values]\n",
    "df['dts']=dtsdf\n",
    "\n",
    "fig,(ax,ax1,ax2)=plt.subplots(3,1,figsize=(3,2.2),sharex=True)\n",
    "ax.yaxis.tick_left()\n",
    "ax1.yaxis.tick_right()\n",
    "ax1.yaxis.set_label_position('right')\n",
    "ax2.yaxis.tick_right()\n",
    "ax2.yaxis.set_label_position('right')\n",
    "for iax in (ax,ax1,ax2):\n",
    "    iax.patch.set_alpha(0)\n",
    "p2,=ax.plot(dts,MesZoBar_old*maxRate*24*3600,'--',color='k',label='this study')\n",
    "p1,=ax.plot(dts,MesZoBar*maxRate*24*3600,'-',color='k',label='Moore-Maley et al. (2016)')\n",
    "ax.set_ylabel('$G_{max}^{MESZ}$ ($\\muup$M N d$^{-1}$)')\n",
    "fig.autofmt_xdate(bottom=0.3, rotation=30, ha='right')\n",
    "yearsFmt = mdates.DateFormatter('%b %d')\n",
    "ax.set_ylim(0,5)\n",
    "ax.xaxis.set_major_formatter(yearsFmt)\n",
    "#ax2=ax.twinx()\n",
    "p3,=ax1.plot(df.loc[:,['dts']],df.loc[:,['Total Biomass']],'.',color='silver', MarkerEdgeWidth=0,MarkerSize=4,label='Mackas et al. (2013)')\n",
    "ax2.plot(df.loc[:,['dts']],df.loc[:,['Total Biomass']],'.',color='silver', MarkerEdgeWidth=0,MarkerSize=4,label='Mackas et al. (2013)')\n",
    "ax2.set_ylabel('Zooplankton Biomass (g m$^{-2}$)')\n",
    "#ax1.set_ylim()\n",
    "ax.set_position(pos=(0.15,.15,.7,.8))\n",
    "ax2.set_position(pos=(0.15,.15,.7,.63))\n",
    "ax1.set_position(pos=(0.15,.8,.7,.15))\n",
    "ax1.spines['bottom'].set_visible(False)\n",
    "ax1.spines['left'].set_visible(False)\n",
    "ax2.spines['top'].set_visible(False)\n",
    "ax1.xaxis.tick_top()\n",
    "ax1.tick_params(labeltop='off')  # don't put tick labels at the top\n",
    "ax2.xaxis.tick_bottom()\n",
    "ax.xaxis.tick_bottom()\n",
    "ax.spines['top'].set_visible(False)\n",
    "ax.spines['right'].set_visible(False)\n",
    "\n",
    "ax1.set_ylim(90,160)\n",
    "ax1.set_yticks((100,130,160))\n",
    "ax2.set_ylim(0,83)\n",
    "\n",
    "d = .015  # how big to make the diagonal lines in axes coordinates\n",
    "# arguments to pass to plot, just so we don't keep repeating them\n",
    "kwargs = dict(transform=ax1.transAxes, color='k', clip_on=False)\n",
    "#ax1.plot((-d, +d), (-d, +d), **kwargs)        # top-left diagonal\n",
    "a=.6/.15\n",
    "ax1.plot((1 - d, 1 + d), (-d*a, +d*a), **kwargs)  # top-right diagonal\n",
    "\n",
    "kwargs.update(transform=ax2.transAxes)  # switch to the bottom axes\n",
    "#ax2.plot((-d, +d), (1 - d, 1 + d), **kwargs)  # bottom-left diagonal\n",
    "ax2.plot((1 - d, 1 + d), (1 - d, 1 + d), **kwargs)  # bottom-right diagonal\n",
    "\n",
    "ax2.yaxis.set_label_coords(1.14,.63)\n",
    "#fig.legend((p1,p2,p3),('this study','Moore-Maley et al. (2016)','Mackas et al. (2013)'),loc=2)\n",
    "fig.savefig('/data/eolson/results/MEOPAR/biomodelevalpaper/figsMod/mesozooWithMackas.eps',dpi=400,transparent=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "with nc.Dataset('/ocean/eolson/MEOPAR/NEMO-forcing/grid/mesh_mask201702_noLPE.nc') as fm:\n",
    "    tmask=np.copy(fm.variables['tmask'][0,0,:,:])\n",
    "    navlon=np.copy(fm.variables['nav_lon'][:,:])\n",
    "    navlat=np.copy(fm.variables['nav_lat'][:,:])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(48.7, 50.1)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUUAAAEwCAYAAAAgm2MvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnV+sHsV9979jH9uI+thgbiorMT7HlvGfU/l1XoMJbyku\npGpruMI0ROplIUnTuwqH8KoRqFLlJuKivaj4W3FVFUyC1MgQFRs4mNZ+jQ2k6rGxXRwbGh+llYjD\n80CFATPvxTzL2bNndndmd2Z3Zvf7kVbPOfvszv6Zme/zm5nf/EZIKUEIIUSxqO0bIISQkKAoEkJI\nCooiIYSkoCgSQkgKiiIhhKSgKBJCSIqxpi4khKDvDyHEC1JK4SqtRi1FKWXuZnu87ZZNFwAefPDB\n3GtXeQYX90biIrZ89VF+295cI5rKPCGELLqWEPlC7/IehRBfpFd0zSJc3E/Va5OwqFqm2yJUsa7D\nqE47e9mNNZ/LqCtUeSTpJenXLRS689NCa3o/JG5iE8MEm7LaVzo70CKEKC2cr7zyinUB0aVJQSSk\nOzTafE7+tv2VrXqPpmnZipULISVxEquFmKVL1qLr5nMrlmJMhacOJtYq6Qax5XNs99skwfQp+ibv\nlzGvcGSPN+mLSaeVHhljAewGtmWIxEln+xSB8iZCUWHOWnm2aeUJJOkeseYtxVxPFKJYpxlaVZCq\n+i+mYaHrBvxR6xetiWKogpGuAC4GVGysTRIWfRDDUOthm7RqKYaWIekKUOXeTPqcul7JuoCNGHIw\nrXtE0XxuAt9ixYoTB30Uw648hyuCE8U2LKmmRhVZ+MLFxDpMhJD52G1aF0WXztouqXoPdUa8STuU\n/Sj2QQi7/nw2tC6KQJjCWKeQtH3vxBz+iJEshc7bQogJAM8COAbgBwAuAvgmgLMAzkkp3xodt1K3\n3zWuJ7NLKY2dt0n82ORpH8WQwSIUJjNabpVSDgBACLEHwGNSyoEQ4lEA3x4d882c/cHDQtAPGLSD\nmGLSfL5bCHGPEGIbgOsTgQQwmTomb39v6YOPG+ke/FEosRSllOcAPAEAQogXAVwUQqxICeAXh+bs\nN0ZnuifNW4oL8UVSvigGJKGsT/FeAM+MxO4iVN/iKgDJ/wl5+63IE8amcelszQrXDjYhvpg/8wnd\nEJmensb09LS39AvjKQoh1gJYB2ACwEEAvwLwLSjhOw7gAwC7ATye3i+l/JkmrdK3HEpGuB4NZ6Vr\nnq7EPWyTUOpjGa7jKbYSZLaIEDLCx3oxrIjNwwAd9QihLprQiSCzoeOjMMRSwAhJ6OuPCEWRdBL+\nCJGqBCeKoUyp8uFSQzed5gihDHWBPr7H4EQxNChihPSLzoji7Czw+OPq0zWJhedKIGkxNkMfrRwf\n9O09dmLhqtlZYN064OOPgaVLgX//d2DDhrbvihASI8FaijZ9i/v3K0EEgE8+AW67DRgOPd6cI2gx\nkljok7UYrCjacMcdykJM+K//Ak6caO9+SDhwKQhiSydEcfVq1WT+0peAJUuAzZuBLVv8X1dnzaaj\nM5tGamZl9QuF0Q19sRaDm9Giw/Qeh0NlIW7ZAoyPV72aPXnTAqvMIe1LwWuauouSEUWIPyqdn+an\nI8SM8A0rrj2mQYP5busRWn3kNL+eEFrBiwWT98YBLlIERTFgWHHbg8KZT9ct7Sj8FEOP7+YTNvvc\nkLy79PvURdHpazmzpcvvipYi6QRFYlclDdJfOi+KwyFw5EgcztxlsEmnxyZuoqm1HUpgkpDp6vuJ\nRhSrZMBwCNx8M/A7v6M+uyCMpD4UPFJENKII2AvjzIzyW/zsM+DkSc5yIfMpi7BOq7ycLv64RCWK\nwPwZI2VMTSlH7jqzXEJsfrPCuqOoPHWxwpNyonDezsPk3stmuQyHyqKcmlr4/XAI3HQT8PbbwKZN\nwOHDzc6UMaHvFVdXBqq+E5dp9Y02f6TpvJ3CxGIcHwduvDFfEG++WW1f+crCWIxHjyrBvHxZfb7+\nusObdwQtRkLcErUo1mVmZk703nkHuOWWsJrJCSYBJdikJm3SJYu616I4NQVMTMz9/+678wdjduxQ\nx4yNqc8bbmj+HpPC1qVC5wL+CBBfRN2nmFDnGWZnlYX47rtqMOa11+Y3tduKvJNQJexVHwTU5F3Y\nvgf2KdanjR+qXkbJKaPuM7QtfHlkKyRFcSFF74Si2DwURZsLeRRFoHsDDnkxGuuc32Xy3o3Ne6Ao\nuqHpusjR5xyaLLwh+i5m76lrPxJVqFsmKIj9pDOiCNg5dlelzI3HB2WDCpzOqMf0hyHv/fKHpRqx\n/5h0ShRdUbSGdJ4bj2vr0dTNZjgEnnpK3VOfpzPGXhFJQKQrn88NgGxjs+XCBSmXLZMSUJ8XLsz/\nfjCQcv169T0g5ZIlUh48KOXUlJSLF0u5erWUp09bX7YSg4GUmzbN3Qug7mMwmDumrffeVH6W5XNZ\nOnnntv2cXdiaYnQtZ1pFSzHDj34EXLqk/r50CXjuufnfj48Dr74KrF8/N6f6o4/mrMfZWeXTWKdZ\nbWp1Hj2qpiAmCAH8zd/MH0GXHW8CJgW56Puy84kfYrXeOy+Kthmzdu38/6+9duExq1cDb74JHDqk\n/BqvvHL+959+Crzwgt19JtTpH5yYaMfBvElirWgkIlyanUUbIjHlBwMpN29WTeHNm+c3RYvO2bBh\nrgmra3abcviwlGNjc03zI0eKrzs1pY6fnNRfs+333mZeFr2DKudwaya/quSvdKhVRn6KQog9AB4D\ncA2ArwF4H8A5KeVbqWO2AbgayvfxJU0a5RdqAJPnreLMPRwC09PA+fPA7t3KmqxCYimePKmfYVPl\nXrtmXZnkYUIVB/iuva+2scmvKjTuvC2EWAngcQDfA3A/gL+WUp4XQuyTUn59dMxuAD+XUr4lhLhH\nSvmkJp0gRBHwn0l1cTnDposV3Db/bKdKdvGdtY3POudaFE1W89sO4BiUKfwOgKtG+1OhFHAQwEtC\niGNQwklqkIQ7c4GU+QvE94WkQpa9h76/J6IoHGgZNYmPAxAAIKV8GMDkaP+vUodOAvgugF8DeMDP\nrZKqJH0lXcFlEFnSDDH94JRZihMAVgG4HsBFIcRLUFahAPBm6ri7pZTfA/CyEGKvlzt1SEhr1hZF\n/s5SNUBE9viYCmgaX3kW6/voK9PT05ienvaWvmmf4j4Az0JZjeugBlT2QQ287IYSynWjUy5KKV/W\npBOGCqVoWxiTQZWk/7BsUKWuKJqkGSJN5FMM76EL+MhLRslxSNuieOSI8kf87DPlCH7oUHFfYl7F\nrfMcMYhB3XxKnjEvnRjeQVeIQRQ777xdRNvr/7pYbdAG3UyZNvobddcr8huzQRcUpCgdCmKzxPC+\ne20pJrh8BzZ9hMnxNu43VeMsJisTnjoFbNy4cGVC08Jq2i+ZjHr77AdMp122hrPufNIOrssELUUP\nuKogVaboFa02aIrJ/ScrE372mX5lwjLLLPtd+v+ic3xSJ1Rc10bkYyL0HyQTP0ViyMyMsvqSEF6v\nv67mRZtajXVwaZXZpJPnB+nbSjSBokeqwOZzhjrvIz1F77rr1L5Tp8xGlk0xabbm3dtNN6moOps2\nLWw++7ovXx3rRZj0H5o2vYk/XJUNjj43RNX3kvQRfvgh8Id/ODey/NOfurEa61hJTS3QVWUFQps0\nbEXRpL+Rwtg8FMWeiGKCL6vRxQp1ttgOHlUVRVOxKxO57PU5CBMurnxtKYoN4EpMdFZjmT+ijjoV\n1UWXgKmDeZq8fkUb0TLFRXOYVmM7uPBDpSg2iCtxtAkHVkRT4phYhx99VF/QsxT1P1YdSTa9hkk6\nFMVmoShGiovmdJUYjdlma90Ka+rTmFiHGzeqfadPzwn6cAjs3w/ccUf1uJGuBz0oaHFTd1YWRbEl\nmnTxyGu2+qr06WfLTj/86U+B3/gNdR/DIbBuHfDxx8AVVwBnzxYLY15/pEnTusr9UxTjJCRRpPN2\ngAyHwD/+Y7PLlqadoLPTD2+4Yc7BfP9+JYiA+ixai6bImZ0DHCRNUHlfNOfU5YYA1opwtflkMJBy\n61a17soVV6jPrVsXrhXj/zmXS2DHguteuKDuC1CfRWvR2Kw3U+eZmn0v3HxtVRmd60yraCkGRnpW\nzOXLwH33KWss2w8pvTflPwRwFCtWzP8FX71aNZmfeKK86Wwa8KJOYI7kPbQd3IPUJ5T8Y59iBXy+\ns/RI9eLFarnUqaniEeumClOV556dBZ5/Hrj99oUC6uq+pSO3HtI+VcoY+xQDwKdVMj6uBPDv/k4J\n4uXL5X2KidnvG9vnHg6BXbuA73xHfdqsYU36SQg/bhTFGvgSx/Fx4BvfUBaiTazFpsTRlGyAjKef\ndhvLMTk/hIpE7AitrKZh87kGvt9d3bnKTbjvFGHaFUA3nH5S9INmU7fYfO4RVWItpqNr+xJtUxGq\n0hVgSqhWBsknO8ob6o8ZRbECoZr+Or/Atu+1rCugasUItUKRhVQdPGkLNp8rEKIgAsULYfkoZDbv\nIa8roE7TmcIYPkVlpGpcTF06LpvPFMUKhCqKJoEnkul1riPTVKVuf2KdNIg/TMqGSb6ZpsM+RU/o\nvNt1xzTB7Czw+OPq05SkD+/QoXy/xuT+6z/Hcgjx1VI3G90Kgtl+T9t7oSCGi0l+Bu9o73J6TNGG\nAKYRFW0hYTOVrg7V39dyCbwlgUvaKYgJ6SmLyXG6fbb3Vf/+ubVZh1ynOzqG0/xcE9Ivl03QhXaY\nArAZwFL8279dyh1Nzvopnjih35clW0jzKPueNIPvfGi6blIUU4QijHfcocJyAepz165272chMwBO\nArgE4G189at6nyHd3GfT+dBpiipcKHlGygm+2TyCAy0aQrA+ZmeVhbhrl30g1zprqpizHMAWACeg\ngkfo31sy6rxqFTA9rQR/fLy6U3o6DmMMFawPNPWjlXcdjj43SAjiaEvdNVXqontns7PA5CRw6RKw\nbBnw859Xj9idhqIYBkVi1dR1XIoim88FxFjpTPrsmuZHP1KCCKjPH//YTbox/mj1Cdf501R9pCh2\njCp9di7R9RutXYvC/+tAYSSuYfPZgNgqXtVAEr76f4ZDYMcO4MwZYMMG4OjR6qsZmhCjhR87TTWh\ndddin2JLxCaMVfDZKV434k8RyX2nr0lhbB5dHaEoFl2IohgFTYwWuiYtihTDdmjSUsxer5WBFiHE\nHiHECiHEhBDiXiHEnUKIbZlj7hVC3CqEuMfVzRE/6KbeJbh0xHVdIYru28f1iDldeveloiiEWAlg\nO4BrANwP4ICU8jkAD6SO2Q3grJTyZSnlk75utk26kum68GI6sXEljq7eW9FyqaRd2phZ5LM+mliK\n2wEcg5qH+A6Aq0b7J1PH/B6ASSHEbiHEbW5vkbgk67Lz+uvtiE2Z1ZelyNUoqZCc9hcWsRoShaI4\naiIfByAAQEr5MJT4bQPwfubw41LKH0NZk52kLJNjmMaUddmRstiv0YXQZN9LFauvzNWIYhgWTdQD\nIQQeeugh9wkXRYsAcCeA2wDsA3APgAkAKwCsBLA3ddx9AP7X6O9/zkmr9SgerraiaB1F34fCYKAW\npk9HrVmyJD9qjZTuItJIKeXhwypKDqCue+SI/X0X0Xb56NvWZh4k15IOo+SUjj6P+hT3AXgWympc\nB+Dq0b5rAOwG8DiAbwE4C+DXUsqXNekUXygydO9N5xoSA7buMnWtgMFAlgbDrUvoFnuXyCvvTeaB\npEtOu5gUgtiE0ZY6BV5K6dVvEaAoNkVROacoll2oQ6IIFFuKRcd0kSqFv4l3Q2F0Q9W8ilUUOfe5\nAqaFhJWyPfjuSVUoihXQVbi8ShjDiHRdkg5qG+q8FxN3nr5Y6T6pkq9Ac2XeVx5TFEnQZAXQxp2H\nwthtfAkvRbEhfP56Vln5zweurQqdAIYYL7ILZK3CWH5QfNznmPMUSaP4impdFVkxKIPOnUkngIkT\nd+LO03S8yK6RFZVYxDDBh6FBS7Eivqy+xHIyTd9XVOs6uKpYulksJmtb+7iXLsJ3o4eWYg3SiyhV\ntZCy6eWln4fPqNZ1yN53lXeTCGDWn3F8HLjxRrt76fpglw11XWxC8U30Bf0Ua5J+fyYFom6Byp7f\ndFTrqthWFh/lsgsVtg5132mZKLb5fum8HRAmTtwuzik63/fsEFfUecY27qFLuHifFEXHdFUUE2yb\ni66FMRZs1m2mtVgf24kGJi2ZumXXB5zR0gjLAdw4+rTHR4W2baqXceYMcP/96rMp2nb7qOqQHBuh\nO16HDC1FLcsBvAZgM4CTAG4G8KHRmaYWY51mt4s51mfOABs3AlICQgCnTqk+yTxMBn2qYvOOfF+z\nK7juu06fR0uxl0xBCeJSAJsAVHeGK+p/qdpczsZ/q8IjjyhBVOkBjz5afk7aXcj2/ouObdJy64Ml\n5Ot99sXKpihqOQ/gM6g4lpcBvNvq3QDuK/OOHcX/m2IqlGWVqSnrI12x0599qOxA+xZdDFAUtayF\ncuEUABYDuLZWaqaVrsnKefvtqrkshPrctcvPdWzE3LUwmoh0X8TQlj5Y1HmwT1FL0qe4CcDbsOlT\nBOachevEXGwiNqONK4+rCqIbLCobQHLhbNwXx2PA7TNWiX7UBnTJaYTlUH2JJ2AjiFnqDrzY5k8S\nNGFqKswQ/2WiqLuWK6fjsrS6IoqAW19Cm/ffBVFk8zmXDwEchUtBLMKFhVQUVst2SVEdbbnQuKTP\nzcIq9LF7gaLYMEX9hnUrbF5YLZcLyTdltbnsX+3TQEoain81KIotUdelRkfe2siuYxCqe67u3J4W\n/5Cszy6Jp67rxfcsq65AUfSIqeXnwu8Q0IfVGg6B995TgggAn34KrFlT+RIAVJpbtw4BvAo1IFVt\n1k/TTVnT6YVdpMoAS1d+IGzhQItn2ixYwyFw003KUkyzfHm9JvSRI6oproT2EoBboPpf69HEu7Kx\nUIua+iZ9wG2Lb5VR/ayFaWvRc6CFlNJmxZiZAd5+e+H+Dz+sN9853UzfunUZBoP/ByCOTnkbC6jI\nes9a93nHhILNMwMLrfi2Bb5JKIoNUHVqXF2mpoAvf1n/3WOPVU9X10zPzg7JCkWfm2Oh0HYwjlig\nKLZEnXnEpoyPA3v36r97/vl6Tegk+nXaFzKviZY3bzsmXLg0NU2frDuXUBRbIK+wui7EwyHwV3+l\n/+4//gOYnnZ6uXkzR0yPj4E6Lk0x/gCkifVHrA4UxcBwKYwzMyokWMKiRcDSpervzz8H7rqr+rKo\neZaTbeWJobL1aVnVuq2WGPKzDIpigLgSxvSAyNQU8OKLwF/8xdz3n3wCPPecfbouncEBd9bU7Czw\nwx+qzeUa2Hn+n6FTNEBUdk7VPOlCk52r+QWKi6CuutXwPvpo/jHXlgQA0s2l1llONqvr+SC9/jUA\nfP/7wLlzbtbAzltV0IaYHKPrlLsu+ILSUuw42QGR3/1dZe0sXqw+d+7MPzfPIgzRctq/f04QAWUF\nv/CCu/R1A0t16WpfXcyCCNB5Owpc55FpyLC0k/aSJcoFJ7EIQ1tBMGspLl3qzlL0jS8RaVNwmxZG\nhg7rKU0X8sRSPHlSWYSJT2KozM4C//AP6u8//uM4BDHBV5TxNmjDUmxcFIUQewA8BuAaAF8D8D6A\nc1LKt3THSSkHmjQoijVpo5CnLULAX6xGMp/YR4BjthRL+xSFECsBbIcSxPsBHJBSPgfggZzjVrm6\nOdI+SV8a4HbEOTRCc86u6h/oQhB1EwtMXXW6EK/SZKBlO4BjUKs4vQPgqtH+iZzjiCfankfdVV89\n1y5GrnERwzJ2oWqSQlEUQmwDcBxqBSdIKR8GMDna/6u844g/2irgTY84D4fAwYPAT36iPn0KVQyC\n72JmSVm5KVvoK4RmeROU+SlOQDWHrwdwUQjxEoCDUOL3Zt5xAJ50f6skTdNBWl346pmiC3k2NQUc\nPlzvunnr1ySCnwwoheBi5BLbMGlV6Yo1WmgpjvoOjwNYOdq1EsDvAfgjAHuFEBNCiPs0x5EO4sNX\nLyHdp5edngio/+tYcEVNZF3UH5v7bRqd1VbHkuuKmDlDFxfOxwbVJ8nN0xYzg4GUW7dKOTamPi9c\nkHJqSkpgbpuaUsfZpnv4sPo8cEDKxYtVWkuWSHnkiJv7Xb9e3W8buMh3l+Vq4fnLJXDj6NN7+Xem\nVZzm1xFcTAv0wXCowpQdPQr86Z8CGzYsPCbbp/fee6qp/PrrwP/8D3DllcANN8xZcCbLuCaW4YkT\nwMaNyqn78mX13XXX2TeR09dM3+877wC33AK8+Wbzbkqh5Lfe0kzWTt8M4CRs105vFZcKW7QhAGuq\nb1vbDAZSbtgwZ+0JIeXp0/rjtm5VFtzWrcUW4YULUq5cqdJbuTLfSjt8WFlygJSLFs23Op9+2v45\nspbs+vVz6Y2N1bM828BFOUq+0593owQujd7RxxLY4busO9Mqzn3uMG27YszMKEsqQUrgqacWHmfa\npzccKgvvgw/U/x98oP7X9eulR8t/8zfnf/ff/23/HFlL9tVXgfXrgbExdZ2uDc4kyBJrNL98zUBZ\niJcAvA0gwCH9PFwqbNGGACwnbuW//i4xtRRNOXx4vsWXbHlW2mCgvjt9Wsply9Sxy5bZ9wHmWbJJ\n+rZ9nSHgopyUn79cKgsxrj5Fzn3uKWX57qqPcjhU0WqOHgW+/W19n2IRZ86o+Ijvvw98/DHwL/+i\nFt5KuPpq4N13y/vzZmfVfezaVW1OdGgBMOriMkJ6CKPXkgEhSNM0VU7SnDmjBkXKGAy6IVRNYypm\nfRNF9inmshzAraOt2mLvXaKN/sm//3uz4/L8F0ObzxwSJnmZNCf7BkVRy3IA/wrgwGj7V1AYFU0K\n45/8SfkxV16pH+QIfT5zQhvCbSqIfYWiqGUKwCao17MIwEYAHR1erIBPqzEtEhs2AG+8AVx11dz3\nX/4ysGaN+nvHDuCXv9Q3nWOYz9yWcJdZgFUEsWzwIibYp6glsRSnRv/PAPg/iMb5tGFclaG0w/WW\nLXPuOcOhcuQG5jtxm6SVDpALhBUP8uBB4A/+QDmVZyObdw3fLQwOtDTCcqj4FoCKiEZBzCNdhuoE\nqiha/qAK2QC5OsFNjmtaLLNBL1wEvIgFT1HG3SXq0r+nxHxu3TePW/ObDTYzW2xJz3BJz31Oz1b5\n0pfq+VFWvZ+xMSkPHmzmuqHhsJw50yr2KZJgqBKtxpS8eJDpvsdf/AL4rd9yu2a0yf1s2aK6BfqI\nDLC/kc1n0gghFH6dA/ZwqETyF7+YO+6JJ4B77mnnfvqIi+a0ZJ8iiQ0pZbCRfM6cURbiJ58AV1wB\nnD1bPOuljT7ILuKyb5GiSKImKXPZStGmYJpOA8wbISf2UBQpiiRDYj1m9zVBVWvP9Qh5U4Rm3boe\ngXYpihxoIa3RxOwY3YyROk7TRQt4hTqtMMTZPSF2o3yBy6Hsog0BuIdwi2ury2CgljFYvHj+cgZ5\n7jk26WZDhmUD0YYUTqzu8/rEYVmhSw7pPnUtyaNHVZPx8mX1mcyKqbtcq24Br5CnFTa9PK0NMkCL\nkWu0kKDxMWLtY7lW02VS2+jba2p5Wtu8CiHkmA4OtJDgqVpGk6l0p06pxavyptFVFarsebOzapGu\n22/Xj2B3eeQ6T+CK8i7U0Wf2KXILYivqZ6pD2ZIBVfsCdYtZlaUTct9eXQIoP+xTJN1BZqyJbCGt\ng67/L41NX2B6dDl73vPPl6cTct8eSeFSYYs2BGCNcAt3a4uiIBSDgbLuBoO5ZU2zlmFyXvp/YG57\n442F16uz2FX6nkIigPLjTqtcJlZ4oQAqHrc4NtMK54oyF5upKSknJ+eELmn6Zs8bDOYLYrJlhTF7\nbVORC9ntJ4AyQ1Hkxi3ZfJDu/1u8WG0Yidz69fmClByT3r7xDf2xtiIXQ59ki2XAmVaxT5FET7I8\ngsvRzHT/36ZNahsbA9avB159dX6UnbJZLDt26Pfb+jayT7IhXCps0YYALApu/dpsyTZl081j01ks\nFy4stBQnJ/VWYJWgunX7JJugpbx2plX0UyS9oaisV/Eh1AWHkFKlc/ny3HFFgSO6HFOxSedsyYAQ\nhNhTVEmrTNPTNWenptS2ZAmwdKlqchc1dctchkwINRBFrNBSJL1DV+Z1q/+Zrhqoi+Z94oRaivW9\n9/xagaazdorOz5vNc+YM8Mgjqk/09turP0MTFqNLS5F9itx6uemIob8uyz/90/z+y2QBrKR/9MKF\nfJefotHv06elFGIu3euus38vuvfsMT/daZXRQcAeACsATAC4F8CdALalvl8J4DYA9wG4LSeN1isC\nN27J1gYXLkj5l38p5Xe+42bVwMFAyrVr54viT34yP2TasmXqUzeQU+Ti893vzk930aJqLkB579pD\nfjYnilCC9wyUID4KYO1o/77UMfem9r+Yk07rFYEbN8CfIBY5Yl+4MH+2ixDzhbHKTJXDh5VYJWli\nZCkeODB/HzSil1wzb/TbhaWYJZsHJvlkkafORNFkoGU71GrwEsA7AK4a7Z9IDpBSPiGlPC+E2Abg\nrEGahDROUuh9UBbdev9+4NNP0/cCPPWU2bl5TE2p/s+EzZv1S6UuWqQf7ClaUnbDBtVP+ed/Djz9\nNHDsmH2fYtp/VNevmM4LnTgV7fdKkWIC2AZlKe7BnCV452j/P2uO3wtgRU5arVsI3Lj5omy2SZGl\nWGemymCgrMODB+f7VybN502b5n/XJLZ5kbcvu73yyivywQcfzKbVTPMZSgBvA7APwD1Q1uEKKKHc\nmzl292j/RE5arVcIbtwAfbPMtoJnMXHETvoU/+zPFjadbZ24ywhl0MgkL0yOKztXOhTFUpccIcRK\nKFF8FsBxAOsAXD3adw2UGJ4D8AOopvObUsoHNOkUX4iQhpBy/iqCBnVAm0aWOo7YXXXi9u2Ok+Sl\n5BKnhLjBpPybiiIpJm85WxfC6VIUuUYL6S11hM3G0iSKWN4Tp/kRUpNQF2AKhbz34zqykStoKZLe\nYrL6nGml9bHqYJcIUfzyoKVIeouJiNkIXUwVn+TDgRZCUG0E2kW6faGBUWiGDiPEFa4txjRlszr6\ngM/n9vEjwtQSAAAMdElEQVSjwz5FQgxx6ULSJ7J+oXXT8g1FkRBL0hWTAtkOPt87RZH0nvTIcZ6D\ncR6m1mP2OPY12tHkDxH7FAlBvs9cXn9gdl8yb7YM0+Nc0JVlCpq2zCmKhBiSFsisuGX3p7+rKoJ1\nRG04VGu//PZvq882hTG2LgaKIiEV0DWzs+JXxyqsGmMx4ZVX1Hozn3+uPv/2b9sTxmwUGtvzEpoS\nV4oiIRWxraQ2bjlVVhdMc/78/P+//321wFUITekiYawinq6hKBJSEduKa+PSo1s+1Ya77lLnpjl1\nyl5cfZEVvzIhbLIJzhkthNRAV3/yRphtQ5DVjbF45gxwyy3AL3+p/p+asl8CtU1shJDxFAkJhLz6\nYyqMPutfsib0yZPA2rVqHZbVq71drjbZ0fy2RJHNZ0I8UNXv0SUzM6rJ/PnnwH/+J/Dee41degFF\nUx3T+79YEqDFEWuKIiGeyBPApoSxbr+kS3TuS1khzB7f1mBLo83ntn8BCPFBlTrUVOTu0NZ+Me1X\ntdWJaPsUKYqkq4QsjFUZDlUTfGrKjaDm1X0XvohR9ymGmPmENEGb/Yu21HUe1+ErZqVr2KdISEOE\n4oeXJm8qYV3n8TzyfBNDEUSgJVEM+ReSkCrYVOoqgSR8UGQN+h6kSXelhSSIQAt9ipl9jVybkKYw\nXUe6qC40VSePHFGC+NlnSvwOHVLBIxJ8DtK4rvtRD7Rk9jVybUKaos6yBXXTsCWxFE+eVNbga681\nP0JdRwPSTe/OiOJofyPXJ6QJXIhinXQSTEeOQ3DZqbsoGEWRkAiwrVe286KLSCzAROxeeAF49113\nrjW+MNUCXYuzU6I4+q6ReyCkKdpsRmf7Cq+9VoUS27KlnSayLUV6kGdYRe2nSAjxS3rk+NprgXPn\n5rvWdGWZAl8EIYp00SFE4WJ2x/i4sggPHQJefVWJZOJas2aNe6fsrhFE83n0fSP3QYgvXNYll6PR\n6cGUmZliN5wQYPN5BK1FEis+HLBdpjc+roRvfHyhU/aaNWxKZzGyFIUQewA8BuAaAF8D8D6Ac1LK\nt0bfrwTwTQBn0/szaRRaiqNjbO+fkGDw8cOuc/SuS2I5rlkD/P7vq5iLGze6icrtIohE25bigpW2\nNCtvrQTwDIAJAI8CWDvavy91zB4AK0Z/P5qTjjQBADduUW8m5TsEDhyQEpjbDh6sl95gIOXWrVKO\njanPwaBaOrbvePR/qZaZbibN5+0Ajo1u5h0AV432T6aOuV5KOdDsJ4RoKIpEHSuugkhIWdwdUfSd\nCwpFUQixDcBxAGJ0Mw8DmBztfz91qBRCrHBxQ74fmBAXlFVcH7gW0R07VDN38WL1ecMN9dJrKtK3\n7x+SsZLvJwCsAnA9gItCiJcAHIQSyTdTxx0bHTcAcDEvsYceeuiLv3fu3ImdO3dWuWdCokVqAi37\n6Dc0YXxc9SO6muaXuAI1MW0wrSWuKR1oGQ2i7APwLJTVuA7A1aN91wDYDeBxAN+CEsTjUsqfadKR\nNhnfpWYF6S5JmbZxoak6pa/JoBGuo25XwXTApZPT/HKO93g3hLghW6ZNrL6qwR+aXNclPXe6ramB\nbYliMH6KWdi3SGKkQSPDW9q+om7HQrCiSEgMVBGnrAuI6Uh0dnDHlzCGtDRqGwTbfE6d5+FuCHFL\nnXpk2yz2vQDWcAgcPQoIoUak2+5TNOmOcNl8Lht9JoQY4GoEOZuOLt3s4I7L0Wtdf2JbtNWFFnzz\nmX2LJBZctGryBLAp+t6fCEQgioTERNU+xirnZ89zIcp9708EIuhTTJ3v8G4I8Usd97Omzs0jb92W\nEHwXdfTGT1FzvsO7IcQ/VcXNxYCkj8GXEHwXdfTGTzEL+xZJl6lbvnVO5HWa1NklC/rU1xiNKBLS\ndYqCTOjELSt6ad/HouPKSKzC9JIFfeprjKb5nErHwd0Q0ix1y76pz172+Cr3kF0NMFmyIIQ1onX0\ntk8xlY6DuyGkWVz2p7uIIVCURmIpnjyprMKQ+g919F4UR2k5SYeQJnE90FhXHMuEMUSrUAdFERRF\nEic+vC98CmMs9Hb0OU0XMpL0D9vBDpNV9mydvesOwvSBKEURoDCSODERIN3ob1maFEd3RCuKhMRK\nmfi06RPYxQW1bIlaFGktki7SlE9gG4tvxUCUAy2ZdJ2nSYhvyupCdvS3aqCIItIhx1wM5LQFB1oy\nxJBphGQpE7nxcSWIMzPlfYq2adsc20ejg0FmCWmJouCw2QAMg4G0shjbWja1C0RvKQK0Fkm8pAc2\n0oLnYrDFds2XsnvsC50QRUK6Rt5gCwdH/BP9QEvmGl7TJ8Q36TpSNNXOtqxXWRArOS8bjKKp9adN\n4TS/4mt4TZ+QpjBZwc5lejZ0XRQ71XwOIYMIcQF/4NujU6JICMmfxucyfR/phkLnRJHWIukKulFf\nm5HgbCTuOgKWvq7unrpE50SRkL6jc/EpW+rAxH0nOTYvja7QSVGktUi6RAh+gqai2QU6KYqEdI1E\ndHyIT9oKrDoHOgThdkVnRbFLv1yEVBUcU7EyGTwpE+au1DnOfSYkcEIRm/R8ap2IdmW5AyNLUQix\nRwixYvT3NiHErUKI2zLHaPe3SYwZQkgal+sa1ZkLnT5Xt7pg9tzk/xib1KWiKIRYCWA7gFVCiN0A\nIKV8GcBE6hjtfkKIe+ougFVn0KTIKkzOzYpjbJg0n7cDODb6+yCAl4QQxwDcnzomb3/rxPprRUgR\nLnwP617PNDxZbGHMCi1FIcQ2AMcBJG9iEsB3AfwawAOpQ/P2E0JqkhderMkf+zwL02YQJ5YR6jJL\ncQLAKgDXQwnepJTyAQAvCyH2po67W0r5Pc3+eTz00ENf/L1z507s3Lmz4m3bQWuRxEpRua1TrrOR\nb8ooskzL0spaimUBJcrSm56exvT0tNF9V6E0Ss6oT3EfgGehrMZ1o68uAjgH4C4AB9L7R32L2XS8\nR8kpgqJIukKRkLhMU0eZSNscb3t+iegydFjFe2j1+oS4wEd/YgiiaHMf2Wu4FEX6KRJCghkM0bn7\n2Db169LZGS06Qsh0QkKl7sp+ec7bVepd+ryml2DolSgSEjNNiINt87iL9E4U+5KxpLvYuraYimld\n0S06tyztkOol+xQJCYysq03yf9V+v+Q8FwM0eYtXpYW6qviFIoy9sxSBcF4+IXno5hPXIW++ctEx\nNmT9EGNx1NZBS5GQgDBt5gLVRczF3GTdLJW8qYB1rtMGvbQUgbgyifQH00ALVdN2PU3Qxr8xFsuR\nliIhkVJ3+qorwyBrIeZZsrEYIr0WRc6JJiFiIyYhleEYBlFM6G3zOSGmzCL9JBTRs8X1YFFT9NpS\nJCQWigYsTAdeTGMfZo81CcrQJXpvKQL9yGjSDVyMOOela5J2rFarDb0WRZ8x2Qjxhc0Ite1sFpPj\nQxFGXyPaFMURtBZJV3jwwQed9uWFVDeacAzvtSgSQswIQRibslA50JJCN+dURyjNB9Jf6s4UqTKP\nui+R8xuNvN3IhQghvSPK5QgIISQG2KdICCEpKIqEEJKiF6IohFgphLhNCLEn8/99o88JIcRxIcQj\nQoi1mfP2CCHuFEJsa+v+bSl73tRxe4QQK1L/a99DDNR45ijzGDAq18n/9wohdqfOizKfazyvXR5n\nnTa7ugGYALB39Pe9ANaO/n4RwFoAKzTn7En2A3i07Wdw9byjz5UAnkn2p85Z8B5i2So+c7R5XPbM\nAG4dbdsA3NeFfK74vFZ53AtLMYuU8gkp5fnRr8bZ0e67hRD3ZH5JrpdSDkZ/TzZ7l+7Ied7tAI5p\nDte9h+iweOZO5DGw8JmllC8DOAfgm1LKhzOHR5/PFs9rlce9FMUUXwdwv5Ty/OgFPwngB6nvZbqp\n1QG+DuD+USE6DmCeG4OU8lzOe4iZwmdG9/IYGD0zoPIUwA+FEPuSLzuYz4XPC8s87pTz9qgfIe1j\n9OvRr8cXh2SO/WsA1wgh7gbwzOjX5GLq+GMAVgHI7g+Cqs8L1QRZBeB6qOd6cnTMvdC/h2Bw/cwI\nPI+BSs+8F6pc/1+o5uI5IcRE6pig89n188Iyj3vjpzgqCHcB+BaA/w1VWX4O4E0AjwFYB1VxDkK9\n9N0AHh8dfxHAcSnlz5q/82oUPa+U8gEhxEoA+wA8C+AlqOf9EVLvQUp5voVbr0zFZ442j4HScv0M\n1A/CVwC8AdW0jDqfKz6vVR73RhQJIcSEvvcpEkLIPCiKhBCSgqJICCEpKIqEEJKCokgIISkoioQQ\nkoKiSAghKf4/jEZznJ9GIq0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fa218f738d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig,ax=plt.subplots(1,1,figsize=(5,5))\n",
    "ax.pcolormesh(navlon,navlat,tmask,cmap=plt.get_cmap('Greys_r'))\n",
    "plt.plot(-1*df['Longitude (deg W)'],df['Latitude (deg N)'],'b.')\n",
    "ax.set_xlim(-125.3,-123)\n",
    "ax.set_ylim(48.7,50.1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}