RDA_effective_rainfall.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "## Although we have the classification scheme of banded and local rainfall, we can further decompose banded rainfall\n",
    "#into changes due to frequency and to intensity, and then see the net change in rainfall induced by each of those.\n",
    "\n",
    "#therefore, we create a new variable - f*i, effective rainfall, and we can calculate changes and compare to the banded rainfall\n",
    "#classification scheme\n",
    "\n",
    "#can also approximate changes as f*delta i + i*delta f\n",
    "\n",
    "#significance testing strategy will be to choose random sets of years and combine total combined front frequency and mean intensity\n",
    "#and thereby effective rainfall"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import time\n",
    "import os\n",
    "import netCDF4 as nc\n",
    "import datetime\n",
    "import matplotlib.pyplot as plt\n",
    "#import masked arrays - important because intensity shows up as a NaN at some points\n",
    "%matplotlib notebook\n",
    "import numpy.ma as ma\n",
    "import scipy.stats as st"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "51\n",
      "07\n",
      "51\n",
      "79\n",
      "80\n",
      "07\n",
      "80\n",
      "93\n",
      "94\n",
      "07\n"
     ]
    }
   ],
   "source": [
    "## LOAD MEAN FREQUENCY DATA FOR TIME PERIODS OF INTEREST\n",
    "loadfile = \"RDA_hov_periods\"\n",
    "ferretpath = \"/Users/Siwen/Desktop/Ferret/bin/\"\n",
    "loadname = ferretpath + loadfile + \".nc\"\n",
    "\n",
    "#import from netCDF file\n",
    "freqload = nc.Dataset(loadname, \"r\")\n",
    "\n",
    "#list of time periods of interest\n",
    "yr_beg = [1951,1951,1980,1980,1994]\n",
    "yr_end = [2007,1979,2007,1993,2007]\n",
    "\n",
    "#in this case, we're just going to load the smoothed plots of hovmoller frequency\n",
    "lats = []\n",
    "\n",
    "for yb, ye in zip(yr_beg, yr_end):\n",
    "\n",
    "    yb_last2= str(yb-1900)\n",
    "    print(yb_last2)\n",
    "    ye_last2= str(\"{0:0=2d}\".format((ye-1900)%100))\n",
    "    print(ye_last2)\n",
    "    \n",
    "    lats.append(freqload['lat_hov_smth_' + yb_last2 + ye_last2][:])\n",
    "    \n",
    "freqload.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "51\n",
      "07\n",
      "51\n",
      "79\n",
      "80\n",
      "07\n",
      "80\n",
      "93\n",
      "94\n",
      "07\n"
     ]
    }
   ],
   "source": [
    "## LOAD MEAN INTENSITY DATA FOR TIME PERIODS OF INTEREST\n",
    "loadfile = \"RDA_int_climo\"\n",
    "ferretpath = \"/Users/Siwen/Desktop/Ferret/bin/\"\n",
    "loadname = ferretpath + loadfile + \".nc\"\n",
    "\n",
    "#import from netCDF file\n",
    "intload = nc.Dataset(loadname, \"r\")\n",
    "\n",
    "#list of time periods of interest\n",
    "yr_beg = [1951,1951,1980,1980,1994]\n",
    "yr_end = [2007,1979,2007,1993,2007]\n",
    "\n",
    "#in this case, we're just going to load the smoothed plots of hovmoller frequency\n",
    "ints = []\n",
    "\n",
    "for yb, ye in zip(yr_beg, yr_end):\n",
    "\n",
    "    yb_last2= str(yb-1900)\n",
    "    print(yb_last2)\n",
    "    ye_last2= str(\"{0:0=2d}\".format((ye-1900)%100))\n",
    "    print(ye_last2)\n",
    "    \n",
    "    ints.append(intload['int_mean_' + yb_last2 + ye_last2][:])\n",
    "    \n",
    "intload.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "rainfall_eff = [f*i for f, i in zip(lats,ints) ] "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "5107\n",
      "5179\n",
      "8007\n",
      "8093\n",
      "9407\n"
     ]
    }
   ],
   "source": [
    "## SAVE AS NETCDF\n",
    "## SAVE FREQUENCY DIFFERENCES AND P-VALUES AS NETCDF\n",
    "outfilename = \"RDA_effective_rainfall.nc\"\n",
    "ferretpath = \"/Users/Siwen/Desktop/Ferret/bin/\"\n",
    "\n",
    "#netCDF output file formatting formatting\n",
    "fullname = ferretpath + outfilename\n",
    "\n",
    "##the following is only to make sure that we don't exceed the bounds of the number of rows of actual data.\n",
    "if os.path.isfile(fullname):\n",
    "    os.remove(fullname)\n",
    "    \n",
    "outfile = nc.Dataset(fullname, 'w')\n",
    "\n",
    "\n",
    "#create dimensions in NetCDF file\n",
    "days_dim = outfile.createDimension('time', None)\n",
    "lat_dim = outfile.createDimension('lat', 22)\n",
    "dim_tuple = (\"time\",\"lat\")\n",
    "\n",
    "#run through loop and save all differences of interest, as well as p-values\n",
    "for yb, ye, rain in zip(yr_beg, yr_end, rainfall_eff):\n",
    "    \n",
    "    #convoluted line of code below guarantees that 2007 will become '07' and not '7'\n",
    "    yb_last2= str(yb-1900)\n",
    "    ye_last2= str(\"{0:0=2d}\".format((ye-1900)%100))\n",
    "    years = yb_last2 + ye_last2\n",
    "    print(years)\n",
    "\n",
    "    #create variables inside netCDF file\n",
    "    rain_out = outfile.createVariable(\"rain_effective_\" + years, 'f8', dim_tuple)\n",
    "    \n",
    "    #fill with actual values\n",
    "    rain_out[:] = rain\n",
    "    \n",
    "outfile.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0.00272904  0.02105263  0.02362573  0.0297076   0.03204678  0.03298246\n",
      "  0.02175439  0.02152047  0.0194152   0.02315789  0.03321637  0.04046784\n",
      "  0.04421053  0.04888889  0.04959064  0.04538012  0.03789474  0.03345029\n",
      "  0.02783626  0.0205848   0.01578947  0.01481481]\n",
      "[37.90512512286032 33.09937650437558 32.62319777050344 32.27439730156701\n",
      " 31.4248218687748 30.81390648987485 27.981854063979544 25.891627473446828\n",
      " 24.665686081931327 26.049530006474722 26.44275616017386 26.39948539585312\n",
      " 26.229321411604467 25.583714726345526 24.61021146486886 23.485019439717135\n",
      " 23.071599384347625 22.633778298903174 22.359608015292316\n",
      " 22.194744497165775 22.830153643391167 22.72927693982279]\n",
      "[0.10344478591034004 0.6968289790394858 0.7707468946949351\n",
      " 0.9587949607716983 1.0070644669057653 1.0163183193151704\n",
      " 0.6087280533216604 0.5571999362706685 0.4788893438129357\n",
      " 0.6032522738341515 0.878332485320395 1.0683300522766292 1.1596121045130396\n",
      " 1.2507593866213367 1.220436217673029 1.065752928960263 0.8742921871963311\n",
      " 0.7571065021621413 0.6224078020631079 0.4568742726901961\n",
      " 0.3604761101588078 0.3367300287381154]\n"
     ]
    }
   ],
   "source": [
    "print(lats[0][200,:])\n",
    "print(ints[0][200,:])\n",
    "print(rainfall_eff[0][200,:])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "application/javascript": [
       "/* Put everything inside the global mpl namespace */\n",
       "window.mpl = {};\n",
       "\n",
       "mpl.get_websocket_type = function() {\n",
       "    if (typeof(WebSocket) !== 'undefined') {\n",
       "        return WebSocket;\n",
       "    } else if (typeof(MozWebSocket) !== 'undefined') {\n",
       "        return MozWebSocket;\n",
       "    } else {\n",
       "        alert('Your browser does not have WebSocket support.' +\n",
       "              'Please try Chrome, Safari or Firefox ≥ 6. ' +\n",
       "              'Firefox 4 and 5 are also supported but you ' +\n",
       "              'have to enable WebSockets in about:config.');\n",
       "    };\n",
       "}\n",
       "\n",
       "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n",
       "    this.id = figure_id;\n",
       "\n",
       "    this.ws = websocket;\n",
       "\n",
       "    this.supports_binary = (this.ws.binaryType != undefined);\n",
       "\n",
       "    if (!this.supports_binary) {\n",
       "        var warnings = document.getElementById(\"mpl-warnings\");\n",
       "        if (warnings) {\n",
       "            warnings.style.display = 'block';\n",
       "            warnings.textContent = (\n",
       "                \"This browser does not support binary websocket messages. \" +\n",
       "                    \"Performance may be slow.\");\n",
       "        }\n",
       "    }\n",
       "\n",
       "    this.imageObj = new Image();\n",
       "\n",
       "    this.context = undefined;\n",
       "    this.message = undefined;\n",
       "    this.canvas = undefined;\n",
       "    this.rubberband_canvas = undefined;\n",
       "    this.rubberband_context = undefined;\n",
       "    this.format_dropdown = undefined;\n",
       "\n",
       "    this.image_mode = 'full';\n",
       "\n",
       "    this.root = $('<div/>');\n",
       "    this._root_extra_style(this.root)\n",
       "    this.root.attr('style', 'display: inline-block');\n",
       "\n",
       "    $(parent_element).append(this.root);\n",
       "\n",
       "    this._init_header(this);\n",
       "    this._init_canvas(this);\n",
       "    this._init_toolbar(this);\n",
       "\n",
       "    var fig = this;\n",
       "\n",
       "    this.waiting = false;\n",
       "\n",
       "    this.ws.onopen =  function () {\n",
       "            fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n",
       "            fig.send_message(\"send_image_mode\", {});\n",
       "            fig.send_message(\"refresh\", {});\n",
       "        }\n",
       "\n",
       "    this.imageObj.onload = function() {\n",
       "            if (fig.image_mode == 'full') {\n",
       "                // Full images could contain transparency (where diff images\n",
       "                // almost always do), so we need to clear the canvas so that\n",
       "                // there is no ghosting.\n",
       "                fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n",
       "            }\n",
       "            fig.context.drawImage(fig.imageObj, 0, 0);\n",
       "        };\n",
       "\n",
       "    this.imageObj.onunload = function() {\n",
       "        this.ws.close();\n",
       "    }\n",
       "\n",
       "    this.ws.onmessage = this._make_on_message_function(this);\n",
       "\n",
       "    this.ondownload = ondownload;\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._init_header = function() {\n",
       "    var titlebar = $(\n",
       "        '<div class=\"ui-dialog-titlebar ui-widget-header ui-corner-all ' +\n",
       "        'ui-helper-clearfix\"/>');\n",
       "    var titletext = $(\n",
       "        '<div class=\"ui-dialog-title\" style=\"width: 100%; ' +\n",
       "        'text-align: center; padding: 3px;\"/>');\n",
       "    titlebar.append(titletext)\n",
       "    this.root.append(titlebar);\n",
       "    this.header = titletext[0];\n",
       "}\n",
       "\n",
       "\n",
       "\n",
       "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n",
       "\n",
       "}\n",
       "\n",
       "\n",
       "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n",
       "\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._init_canvas = function() {\n",
       "    var fig = this;\n",
       "\n",
       "    var canvas_div = $('<div/>');\n",
       "\n",
       "    canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n",
       "\n",
       "    function canvas_keyboard_event(event) {\n",
       "        return fig.key_event(event, event['data']);\n",
       "    }\n",
       "\n",
       "    canvas_div.keydown('key_press', canvas_keyboard_event);\n",
       "    canvas_div.keyup('key_release', canvas_keyboard_event);\n",
       "    this.canvas_div = canvas_div\n",
       "    this._canvas_extra_style(canvas_div)\n",
       "    this.root.append(canvas_div);\n",
       "\n",
       "    var canvas = $('<canvas/>');\n",
       "    canvas.addClass('mpl-canvas');\n",
       "    canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n",
       "\n",
       "    this.canvas = canvas[0];\n",
       "    this.context = canvas[0].getContext(\"2d\");\n",
       "\n",
       "    var rubberband = $('<canvas/>');\n",
       "    rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n",
       "\n",
       "    var pass_mouse_events = true;\n",
       "\n",
       "    canvas_div.resizable({\n",
       "        start: function(event, ui) {\n",
       "            pass_mouse_events = false;\n",
       "        },\n",
       "        resize: function(event, ui) {\n",
       "            fig.request_resize(ui.size.width, ui.size.height);\n",
       "        },\n",
       "        stop: function(event, ui) {\n",
       "            pass_mouse_events = true;\n",
       "            fig.request_resize(ui.size.width, ui.size.height);\n",
       "        },\n",
       "    });\n",
       "\n",
       "    function mouse_event_fn(event) {\n",
       "        if (pass_mouse_events)\n",
       "            return fig.mouse_event(event, event['data']);\n",
       "    }\n",
       "\n",
       "    rubberband.mousedown('button_press', mouse_event_fn);\n",
       "    rubberband.mouseup('button_release', mouse_event_fn);\n",
       "    // Throttle sequential mouse events to 1 every 20ms.\n",
       "    rubberband.mousemove('motion_notify', mouse_event_fn);\n",
       "\n",
       "    rubberband.mouseenter('figure_enter', mouse_event_fn);\n",
       "    rubberband.mouseleave('figure_leave', mouse_event_fn);\n",
       "\n",
       "    canvas_div.on(\"wheel\", function (event) {\n",
       "        event = event.originalEvent;\n",
       "        event['data'] = 'scroll'\n",
       "        if (event.deltaY < 0) {\n",
       "            event.step = 1;\n",
       "        } else {\n",
       "            event.step = -1;\n",
       "        }\n",
       "        mouse_event_fn(event);\n",
       "    });\n",
       "\n",
       "    canvas_div.append(canvas);\n",
       "    canvas_div.append(rubberband);\n",
       "\n",
       "    this.rubberband = rubberband;\n",
       "    this.rubberband_canvas = rubberband[0];\n",
       "    this.rubberband_context = rubberband[0].getContext(\"2d\");\n",
       "    this.rubberband_context.strokeStyle = \"#000000\";\n",
       "\n",
       "    this._resize_canvas = function(width, height) {\n",
       "        // Keep the size of the canvas, canvas container, and rubber band\n",
       "        // canvas in synch.\n",
       "        canvas_div.css('width', width)\n",
       "        canvas_div.css('height', height)\n",
       "\n",
       "        canvas.attr('width', width);\n",
       "        canvas.attr('height', height);\n",
       "\n",
       "        rubberband.attr('width', width);\n",
       "        rubberband.attr('height', height);\n",
       "    }\n",
       "\n",
       "    // Set the figure to an initial 600x600px, this will subsequently be updated\n",
       "    // upon first draw.\n",
       "    this._resize_canvas(600, 600);\n",
       "\n",
       "    // Disable right mouse context menu.\n",
       "    $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n",
       "        return false;\n",
       "    });\n",
       "\n",
       "    function set_focus () {\n",
       "        canvas.focus();\n",
       "        canvas_div.focus();\n",
       "    }\n",
       "\n",
       "    window.setTimeout(set_focus, 100);\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._init_toolbar = function() {\n",
       "    var fig = this;\n",
       "\n",
       "    var nav_element = $('<div/>')\n",
       "    nav_element.attr('style', 'width: 100%');\n",
       "    this.root.append(nav_element);\n",
       "\n",
       "    // Define a callback function for later on.\n",
       "    function toolbar_event(event) {\n",
       "        return fig.toolbar_button_onclick(event['data']);\n",
       "    }\n",
       "    function toolbar_mouse_event(event) {\n",
       "        return fig.toolbar_button_onmouseover(event['data']);\n",
       "    }\n",
       "\n",
       "    for(var toolbar_ind in mpl.toolbar_items) {\n",
       "        var name = mpl.toolbar_items[toolbar_ind][0];\n",
       "        var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
       "        var image = mpl.toolbar_items[toolbar_ind][2];\n",
       "        var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
       "\n",
       "        if (!name) {\n",
       "            // put a spacer in here.\n",
       "            continue;\n",
       "        }\n",
       "        var button = $('<button/>');\n",
       "        button.addClass('ui-button ui-widget ui-state-default ui-corner-all ' +\n",
       "                        'ui-button-icon-only');\n",
       "        button.attr('role', 'button');\n",
       "        button.attr('aria-disabled', 'false');\n",
       "        button.click(method_name, toolbar_event);\n",
       "        button.mouseover(tooltip, toolbar_mouse_event);\n",
       "\n",
       "        var icon_img = $('<span/>');\n",
       "        icon_img.addClass('ui-button-icon-primary ui-icon');\n",
       "        icon_img.addClass(image);\n",
       "        icon_img.addClass('ui-corner-all');\n",
       "\n",
       "        var tooltip_span = $('<span/>');\n",
       "        tooltip_span.addClass('ui-button-text');\n",
       "        tooltip_span.html(tooltip);\n",
       "\n",
       "        button.append(icon_img);\n",
       "        button.append(tooltip_span);\n",
       "\n",
       "        nav_element.append(button);\n",
       "    }\n",
       "\n",
       "    var fmt_picker_span = $('<span/>');\n",
       "\n",
       "    var fmt_picker = $('<select/>');\n",
       "    fmt_picker.addClass('mpl-toolbar-option ui-widget ui-widget-content');\n",
       "    fmt_picker_span.append(fmt_picker);\n",
       "    nav_element.append(fmt_picker_span);\n",
       "    this.format_dropdown = fmt_picker[0];\n",
       "\n",
       "    for (var ind in mpl.extensions) {\n",
       "        var fmt = mpl.extensions[ind];\n",
       "        var option = $(\n",
       "            '<option/>', {selected: fmt === mpl.default_extension}).html(fmt);\n",
       "        fmt_picker.append(option)\n",
       "    }\n",
       "\n",
       "    // Add hover states to the ui-buttons\n",
       "    $( \".ui-button\" ).hover(\n",
       "        function() { $(this).addClass(\"ui-state-hover\");},\n",
       "        function() { $(this).removeClass(\"ui-state-hover\");}\n",
       "    );\n",
       "\n",
       "    var status_bar = $('<span class=\"mpl-message\"/>');\n",
       "    nav_element.append(status_bar);\n",
       "    this.message = status_bar[0];\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.request_resize = function(x_pixels, y_pixels) {\n",
       "    // Request matplotlib to resize the figure. Matplotlib will then trigger a resize in the client,\n",
       "    // which will in turn request a refresh of the image.\n",
       "    this.send_message('resize', {'width': x_pixels, 'height': y_pixels});\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.send_message = function(type, properties) {\n",
       "    properties['type'] = type;\n",
       "    properties['figure_id'] = this.id;\n",
       "    this.ws.send(JSON.stringify(properties));\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.send_draw_message = function() {\n",
       "    if (!this.waiting) {\n",
       "        this.waiting = true;\n",
       "        this.ws.send(JSON.stringify({type: \"draw\", figure_id: this.id}));\n",
       "    }\n",
       "}\n",
       "\n",
       "\n",
       "mpl.figure.prototype.handle_save = function(fig, msg) {\n",
       "    var format_dropdown = fig.format_dropdown;\n",
       "    var format = format_dropdown.options[format_dropdown.selectedIndex].value;\n",
       "    fig.ondownload(fig, format);\n",
       "}\n",
       "\n",
       "\n",
       "mpl.figure.prototype.handle_resize = function(fig, msg) {\n",
       "    var size = msg['size'];\n",
       "    if (size[0] != fig.canvas.width || size[1] != fig.canvas.height) {\n",
       "        fig._resize_canvas(size[0], size[1]);\n",
       "        fig.send_message(\"refresh\", {});\n",
       "    };\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_rubberband = function(fig, msg) {\n",
       "    var x0 = msg['x0'];\n",
       "    var y0 = fig.canvas.height - msg['y0'];\n",
       "    var x1 = msg['x1'];\n",
       "    var y1 = fig.canvas.height - msg['y1'];\n",
       "    x0 = Math.floor(x0) + 0.5;\n",
       "    y0 = Math.floor(y0) + 0.5;\n",
       "    x1 = Math.floor(x1) + 0.5;\n",
       "    y1 = Math.floor(y1) + 0.5;\n",
       "    var min_x = Math.min(x0, x1);\n",
       "    var min_y = Math.min(y0, y1);\n",
       "    var width = Math.abs(x1 - x0);\n",
       "    var height = Math.abs(y1 - y0);\n",
       "\n",
       "    fig.rubberband_context.clearRect(\n",
       "        0, 0, fig.canvas.width, fig.canvas.height);\n",
       "\n",
       "    fig.rubberband_context.strokeRect(min_x, min_y, width, height);\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_figure_label = function(fig, msg) {\n",
       "    // Updates the figure title.\n",
       "    fig.header.textContent = msg['label'];\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_cursor = function(fig, msg) {\n",
       "    var cursor = msg['cursor'];\n",
       "    switch(cursor)\n",
       "    {\n",
       "    case 0:\n",
       "        cursor = 'pointer';\n",
       "        break;\n",
       "    case 1:\n",
       "        cursor = 'default';\n",
       "        break;\n",
       "    case 2:\n",
       "        cursor = 'crosshair';\n",
       "        break;\n",
       "    case 3:\n",
       "        cursor = 'move';\n",
       "        break;\n",
       "    }\n",
       "    fig.rubberband_canvas.style.cursor = cursor;\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_message = function(fig, msg) {\n",
       "    fig.message.textContent = msg['message'];\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_draw = function(fig, msg) {\n",
       "    // Request the server to send over a new figure.\n",
       "    fig.send_draw_message();\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_image_mode = function(fig, msg) {\n",
       "    fig.image_mode = msg['mode'];\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.updated_canvas_event = function() {\n",
       "    // Called whenever the canvas gets updated.\n",
       "    this.send_message(\"ack\", {});\n",
       "}\n",
       "\n",
       "// A function to construct a web socket function for onmessage handling.\n",
       "// Called in the figure constructor.\n",
       "mpl.figure.prototype._make_on_message_function = function(fig) {\n",
       "    return function socket_on_message(evt) {\n",
       "        if (evt.data instanceof Blob) {\n",
       "            /* FIXME: We get \"Resource interpreted as Image but\n",
       "             * transferred with MIME type text/plain:\" errors on\n",
       "             * Chrome.  But how to set the MIME type?  It doesn't seem\n",
       "             * to be part of the websocket stream */\n",
       "            evt.data.type = \"image/png\";\n",
       "\n",
       "            /* Free the memory for the previous frames */\n",
       "            if (fig.imageObj.src) {\n",
       "                (window.URL || window.webkitURL).revokeObjectURL(\n",
       "                    fig.imageObj.src);\n",
       "            }\n",
       "\n",
       "            fig.imageObj.src = (window.URL || window.webkitURL).createObjectURL(\n",
       "                evt.data);\n",
       "            fig.updated_canvas_event();\n",
       "            fig.waiting = false;\n",
       "            return;\n",
       "        }\n",
       "        else if (typeof evt.data === 'string' && evt.data.slice(0, 21) == \"data:image/png;base64\") {\n",
       "            fig.imageObj.src = evt.data;\n",
       "            fig.updated_canvas_event();\n",
       "            fig.waiting = false;\n",
       "            return;\n",
       "        }\n",
       "\n",
       "        var msg = JSON.parse(evt.data);\n",
       "        var msg_type = msg['type'];\n",
       "\n",
       "        // Call the  \"handle_{type}\" callback, which takes\n",
       "        // the figure and JSON message as its only arguments.\n",
       "        try {\n",
       "            var callback = fig[\"handle_\" + msg_type];\n",
       "        } catch (e) {\n",
       "            console.log(\"No handler for the '\" + msg_type + \"' message type: \", msg);\n",
       "            return;\n",
       "        }\n",
       "\n",
       "        if (callback) {\n",
       "            try {\n",
       "                // console.log(\"Handling '\" + msg_type + \"' message: \", msg);\n",
       "                callback(fig, msg);\n",
       "            } catch (e) {\n",
       "                console.log(\"Exception inside the 'handler_\" + msg_type + \"' callback:\", e, e.stack, msg);\n",
       "            }\n",
       "        }\n",
       "    };\n",
       "}\n",
       "\n",
       "// from http://stackoverflow.com/questions/1114465/getting-mouse-location-in-canvas\n",
       "mpl.findpos = function(e) {\n",
       "    //this section is from http://www.quirksmode.org/js/events_properties.html\n",
       "    var targ;\n",
       "    if (!e)\n",
       "        e = window.event;\n",
       "    if (e.target)\n",
       "        targ = e.target;\n",
       "    else if (e.srcElement)\n",
       "        targ = e.srcElement;\n",
       "    if (targ.nodeType == 3) // defeat Safari bug\n",
       "        targ = targ.parentNode;\n",
       "\n",
       "    // jQuery normalizes the pageX and pageY\n",
       "    // pageX,Y are the mouse positions relative to the document\n",
       "    // offset() returns the position of the element relative to the document\n",
       "    var x = e.pageX - $(targ).offset().left;\n",
       "    var y = e.pageY - $(targ).offset().top;\n",
       "\n",
       "    return {\"x\": x, \"y\": y};\n",
       "};\n",
       "\n",
       "/*\n",
       " * return a copy of an object with only non-object keys\n",
       " * we need this to avoid circular references\n",
       " * http://stackoverflow.com/a/24161582/3208463\n",
       " */\n",
       "function simpleKeys (original) {\n",
       "  return Object.keys(original).reduce(function (obj, key) {\n",
       "    if (typeof original[key] !== 'object')\n",
       "        obj[key] = original[key]\n",
       "    return obj;\n",
       "  }, {});\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.mouse_event = function(event, name) {\n",
       "    var canvas_pos = mpl.findpos(event)\n",
       "\n",
       "    if (name === 'button_press')\n",
       "    {\n",
       "        this.canvas.focus();\n",
       "        this.canvas_div.focus();\n",
       "    }\n",
       "\n",
       "    var x = canvas_pos.x;\n",
       "    var y = canvas_pos.y;\n",
       "\n",
       "    this.send_message(name, {x: x, y: y, button: event.button,\n",
       "                             step: event.step,\n",
       "                             guiEvent: simpleKeys(event)});\n",
       "\n",
       "    /* This prevents the web browser from automatically changing to\n",
       "     * the text insertion cursor when the button is pressed.  We want\n",
       "     * to control all of the cursor setting manually through the\n",
       "     * 'cursor' event from matplotlib */\n",
       "    event.preventDefault();\n",
       "    return false;\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._key_event_extra = function(event, name) {\n",
       "    // Handle any extra behaviour associated with a key event\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.key_event = function(event, name) {\n",
       "\n",
       "    // Prevent repeat events\n",
       "    if (name == 'key_press')\n",
       "    {\n",
       "        if (event.which === this._key)\n",
       "            return;\n",
       "        else\n",
       "            this._key = event.which;\n",
       "    }\n",
       "    if (name == 'key_release')\n",
       "        this._key = null;\n",
       "\n",
       "    var value = '';\n",
       "    if (event.ctrlKey && event.which != 17)\n",
       "        value += \"ctrl+\";\n",
       "    if (event.altKey && event.which != 18)\n",
       "        value += \"alt+\";\n",
       "    if (event.shiftKey && event.which != 16)\n",
       "        value += \"shift+\";\n",
       "\n",
       "    value += 'k';\n",
       "    value += event.which.toString();\n",
       "\n",
       "    this._key_event_extra(event, name);\n",
       "\n",
       "    this.send_message(name, {key: value,\n",
       "                             guiEvent: simpleKeys(event)});\n",
       "    return false;\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.toolbar_button_onclick = function(name) {\n",
       "    if (name == 'download') {\n",
       "        this.handle_save(this, null);\n",
       "    } else {\n",
       "        this.send_message(\"toolbar_button\", {name: name});\n",
       "    }\n",
       "};\n",
       "\n",
       "mpl.figure.prototype.toolbar_button_onmouseover = function(tooltip) {\n",
       "    this.message.textContent = tooltip;\n",
       "};\n",
       "mpl.toolbar_items = [[\"Home\", \"Reset original view\", \"fa fa-home icon-home\", \"home\"], [\"Back\", \"Back to  previous view\", \"fa fa-arrow-left icon-arrow-left\", \"back\"], [\"Forward\", \"Forward to next view\", \"fa fa-arrow-right icon-arrow-right\", \"forward\"], [\"\", \"\", \"\", \"\"], [\"Pan\", \"Pan axes with left mouse, zoom with right\", \"fa fa-arrows icon-move\", \"pan\"], [\"Zoom\", \"Zoom to rectangle\", \"fa fa-square-o icon-check-empty\", \"zoom\"], [\"\", \"\", \"\", \"\"], [\"Download\", \"Download plot\", \"fa fa-floppy-o icon-save\", \"download\"]];\n",
       "\n",
       "mpl.extensions = [\"eps\", \"jpeg\", \"pdf\", \"png\", \"ps\", \"raw\", \"svg\", \"tif\"];\n",
       "\n",
       "mpl.default_extension = \"png\";var comm_websocket_adapter = function(comm) {\n",
       "    // Create a \"websocket\"-like object which calls the given IPython comm\n",
       "    // object with the appropriate methods. Currently this is a non binary\n",
       "    // socket, so there is still some room for performance tuning.\n",
       "    var ws = {};\n",
       "\n",
       "    ws.close = function() {\n",
       "        comm.close()\n",
       "    };\n",
       "    ws.send = function(m) {\n",
       "        //console.log('sending', m);\n",
       "        comm.send(m);\n",
       "    };\n",
       "    // Register the callback with on_msg.\n",
       "    comm.on_msg(function(msg) {\n",
       "        //console.log('receiving', msg['content']['data'], msg);\n",
       "        // Pass the mpl event to the overriden (by mpl) onmessage function.\n",
       "        ws.onmessage(msg['content']['data'])\n",
       "    });\n",
       "    return ws;\n",
       "}\n",
       "\n",
       "mpl.mpl_figure_comm = function(comm, msg) {\n",
       "    // This is the function which gets called when the mpl process\n",
       "    // starts-up an IPython Comm through the \"matplotlib\" channel.\n",
       "\n",
       "    var id = msg.content.data.id;\n",
       "    // Get hold of the div created by the display call when the Comm\n",
       "    // socket was opened in Python.\n",
       "    var element = $(\"#\" + id);\n",
       "    var ws_proxy = comm_websocket_adapter(comm)\n",
       "\n",
       "    function ondownload(figure, format) {\n",
       "        window.open(figure.imageObj.src);\n",
       "    }\n",
       "\n",
       "    var fig = new mpl.figure(id, ws_proxy,\n",
       "                           ondownload,\n",
       "                           element.get(0));\n",
       "\n",
       "    // Call onopen now - mpl needs it, as it is assuming we've passed it a real\n",
       "    // web socket which is closed, not our websocket->open comm proxy.\n",
       "    ws_proxy.onopen();\n",
       "\n",
       "    fig.parent_element = element.get(0);\n",
       "    fig.cell_info = mpl.find_output_cell(\"<div id='\" + id + \"'></div>\");\n",
       "    if (!fig.cell_info) {\n",
       "        console.error(\"Failed to find cell for figure\", id, fig);\n",
       "        return;\n",
       "    }\n",
       "\n",
       "    var output_index = fig.cell_info[2]\n",
       "    var cell = fig.cell_info[0];\n",
       "\n",
       "};\n",
       "\n",
       "mpl.figure.prototype.handle_close = function(fig, msg) {\n",
       "    fig.root.unbind('remove')\n",
       "\n",
       "    // Update the output cell to use the data from the current canvas.\n",
       "    fig.push_to_output();\n",
       "    var dataURL = fig.canvas.toDataURL();\n",
       "    // Re-enable the keyboard manager in IPython - without this line, in FF,\n",
       "    // the notebook keyboard shortcuts fail.\n",
       "    IPython.keyboard_manager.enable()\n",
       "    $(fig.parent_element).html('<img src=\"' + dataURL + '\">');\n",
       "    fig.close_ws(fig, msg);\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.close_ws = function(fig, msg){\n",
       "    fig.send_message('closing', msg);\n",
       "    // fig.ws.close()\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.push_to_output = function(remove_interactive) {\n",
       "    // Turn the data on the canvas into data in the output cell.\n",
       "    var dataURL = this.canvas.toDataURL();\n",
       "    this.cell_info[1]['text/html'] = '<img src=\"' + dataURL + '\">';\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.updated_canvas_event = function() {\n",
       "    // Tell IPython that the notebook contents must change.\n",
       "    IPython.notebook.set_dirty(true);\n",
       "    this.send_message(\"ack\", {});\n",
       "    var fig = this;\n",
       "    // Wait a second, then push the new image to the DOM so\n",
       "    // that it is saved nicely (might be nice to debounce this).\n",
       "    setTimeout(function () { fig.push_to_output() }, 1000);\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._init_toolbar = function() {\n",
       "    var fig = this;\n",
       "\n",
       "    var nav_element = $('<div/>')\n",
       "    nav_element.attr('style', 'width: 100%');\n",
       "    this.root.append(nav_element);\n",
       "\n",
       "    // Define a callback function for later on.\n",
       "    function toolbar_event(event) {\n",
       "        return fig.toolbar_button_onclick(event['data']);\n",
       "    }\n",
       "    function toolbar_mouse_event(event) {\n",
       "        return fig.toolbar_button_onmouseover(event['data']);\n",
       "    }\n",
       "\n",
       "    for(var toolbar_ind in mpl.toolbar_items){\n",
       "        var name = mpl.toolbar_items[toolbar_ind][0];\n",
       "        var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
       "        var image = mpl.toolbar_items[toolbar_ind][2];\n",
       "        var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
       "\n",
       "        if (!name) { continue; };\n",
       "\n",
       "        var button = $('<button class=\"btn btn-default\" href=\"#\" title=\"' + name + '\"><i class=\"fa ' + image + ' fa-lg\"></i></button>');\n",
       "        button.click(method_name, toolbar_event);\n",
       "        button.mouseover(tooltip, toolbar_mouse_event);\n",
       "        nav_element.append(button);\n",
       "    }\n",
       "\n",
       "    // Add the status bar.\n",
       "    var status_bar = $('<span class=\"mpl-message\" style=\"text-align:right; float: right;\"/>');\n",
       "    nav_element.append(status_bar);\n",
       "    this.message = status_bar[0];\n",
       "\n",
       "    // Add the close button to the window.\n",
       "    var buttongrp = $('<div class=\"btn-group inline pull-right\"></div>');\n",
       "    var button = $('<button class=\"btn btn-mini btn-primary\" href=\"#\" title=\"Stop Interaction\"><i class=\"fa fa-power-off icon-remove icon-large\"></i></button>');\n",
       "    button.click(function (evt) { fig.handle_close(fig, {}); } );\n",
       "    button.mouseover('Stop Interaction', toolbar_mouse_event);\n",
       "    buttongrp.append(button);\n",
       "    var titlebar = this.root.find($('.ui-dialog-titlebar'));\n",
       "    titlebar.prepend(buttongrp);\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._root_extra_style = function(el){\n",
       "    var fig = this\n",
       "    el.on(\"remove\", function(){\n",
       "\tfig.close_ws(fig, {});\n",
       "    });\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._canvas_extra_style = function(el){\n",
       "    // this is important to make the div 'focusable\n",
       "    el.attr('tabindex', 0)\n",
       "    // reach out to IPython and tell the keyboard manager to turn it's self\n",
       "    // off when our div gets focus\n",
       "\n",
       "    // location in version 3\n",
       "    if (IPython.notebook.keyboard_manager) {\n",
       "        IPython.notebook.keyboard_manager.register_events(el);\n",
       "    }\n",
       "    else {\n",
       "        // location in version 2\n",
       "        IPython.keyboard_manager.register_events(el);\n",
       "    }\n",
       "\n",
       "}\n",
       "\n",
       "mpl.figure.prototype._key_event_extra = function(event, name) {\n",
       "    var manager = IPython.notebook.keyboard_manager;\n",
       "    if (!manager)\n",
       "        manager = IPython.keyboard_manager;\n",
       "\n",
       "    // Check for shift+enter\n",
       "    if (event.shiftKey && event.which == 13) {\n",
       "        this.canvas_div.blur();\n",
       "        event.shiftKey = false;\n",
       "        // Send a \"J\" for go to next cell\n",
       "        event.which = 74;\n",
       "        event.keyCode = 74;\n",
       "        manager.command_mode();\n",
       "        manager.handle_keydown(event);\n",
       "    }\n",
       "}\n",
       "\n",
       "mpl.figure.prototype.handle_save = function(fig, msg) {\n",
       "    fig.ondownload(fig, null);\n",
       "}\n",
       "\n",
       "\n",
       "mpl.find_output_cell = function(html_output) {\n",
       "    // Return the cell and output element which can be found *uniquely* in the notebook.\n",
       "    // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n",
       "    // IPython event is triggered only after the cells have been serialised, which for\n",
       "    // our purposes (turning an active figure into a static one), is too late.\n",
       "    var cells = IPython.notebook.get_cells();\n",
       "    var ncells = cells.length;\n",
       "    for (var i=0; i<ncells; i++) {\n",
       "        var cell = cells[i];\n",
       "        if (cell.cell_type === 'code'){\n",
       "            for (var j=0; j<cell.output_area.outputs.length; j++) {\n",
       "                var data = cell.output_area.outputs[j];\n",
       "                if (data.data) {\n",
       "                    // IPython >= 3 moved mimebundle to data attribute of output\n",
       "                    data = data.data;\n",
       "                }\n",
       "                if (data['text/html'] == html_output) {\n",
       "                    return [cell, data, j];\n",
       "                }\n",
       "            }\n",
       "        }\n",
       "    }\n",
       "}\n",
       "\n",
       "// Register the function which deals with the matplotlib target/channel.\n",
       "// The kernel may be null if the page has been refreshed.\n",
       "if (IPython.notebook.kernel != null) {\n",
       "    IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n",
       "}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAoAAAAHgCAYAAAA10dzkAAAgAElEQVR4Xuy9CZheaVnm/0uquiqppPJVOiHpbgk0KK0giqgMcIEyyi6bC4i4gH+UxQFHbUdUQJBNBAZGRxpEcRlQgREFUXQyso/N4KA9jmzSNnRj2nSnTbprSSqp6qrU/7rP9z7pt7+u5VS9Jyf11Xef6/rqbM+73c97zvur52zb8GQFrIAVsAJWwApYASswUApsG6jWurFWwApYAStgBayAFbACGADdCayAFbACVsAKWAErMGAKGAAHzOFurhWwAlbAClgBK2AFDIDuA1bAClgBK2AFrIAVGDAFDIAD5nA31wpYAStgBayAFbACBkD3AStgBayAFbACVsAKDJgCBsABc7ibawWsgBWwAlbAClgBA6D7gBWwAlbAClgBK2AFBkwBA+CAOdzNtQJWwApYAStgBayAAdB9wApYAStgBayAFbACA6aAAXDAHO7mWgErYAWsgBWwAlbAAOg+YAWsgBWwAlbACliBAVPAADhgDndzrYAVsAJWwApYAStgAHQfsAJWwApYAStgBazAgClgABwwh7u5VsAKWAErYAWsgBUwALoPWAErYAWsgBWwAlZgwBQwAA6Yw91cK2AFrIAVsAJWwAoYAN0HrIAVsAJWwApYASswYAoYAAfM4W6uFbACVsAKWAErYAUMgO4DVsAKWAErYAWsgBUYMAUMgAPmcDfXClgBK2AFrIAVsAIGQPcBK2AFrIAVsAJWwAoMmAIGwAFzuJtrBayAFbACVsAKWAEDoPuAFbACVsAKWAErYAUGTAED4IA53M21AlbAClgBK2AFrIAB0H3AClgBK2AFrIAVsAIDpoABcMAc7uZaAStgBayAFbACVsAA6D5gBayAFbACVsAKWIEBU8AAOGAOd3OtgBWwAlbAClgBK2AAdB+wAlbAClgBK2AFrMCAKWAAHDCHu7lWwApYAStgBayAFTAAug9YAStgBayAFbACVmDAFDAADpjD3VwrYAWsgBWwAlbAChgA3QesgBWwAlbAClgBKzBgChgAB8zhbq4VsAJWwApYAStgBQyA7gNWwApYAStgBayAFRgwBQyAA+ZwN9cKWAErYAWsgBWwAgZA9wErYAWsgBWwAlbACgyYAgbAAXO4m2sFrIAVsAJWwApYAQOg+4AVsAJWwApYAStgBQZMAQPggDnczbUCVsAKWAErYAWsgAHQfcAKWAErYAWsgBWwAgOmwCAC4GuBJwD3BE4CHwdeBNyY+f4scAZYAKTREvBQ4HMD1j/cXCtgBayAFbACVmALKjCIAPga4L3AZ4Ax4K3A/YAH9gDgI4GPbkGfu0lWwApYAStgBazAgCswiADY6/IHANcAFwNTaacigI8CPjLg/cPNtwJWwApYAStgBbagAgbA7uXf5wP37okA3gxcBHwF+E3g7VvQ/26SFbACVsAKWAErMIAKDDoAKsr3PuB7gb/O/P8dwCeBReDRwB8Cvwi8raePSL/LgJkB7DtushWwAlbACqxfgXHgaLq3fP2p66XYAYzUM13Taj7dE7+moQ36S4FBBsAnAu8EngV8YA23vQx4DPDwHruv6nl4pL+879paAStgBazAhVDg7sC/nqeCd+yG03rCsaFJV8PuZQhsSM1NlM2gAuAPAW8GngZ8qIY/BICPBR7WY7tH9w0eOXKEPXu02L/Ti1/8Yn7lV36lfxsAbIU2yAFboR1boQ1bxRdbpR1boU9NT09z6NAhuaQDTJ+nE241Lv0MMFpYwBzwX7p5nM/6FtbSyTeqwCAC4AuBVwJPAq5eRjg9DSxd9JSwHgbR08DvAgSBVy0HgFNTU30PgFdeeSVvetObNtqPNkW6rdAGCbkV2rEV2rBVfLFV2rEV+pQAsNMRS51XoKoA8KWArgOXTHoX2qsNgCUSbuq0gwiAgrrbAf1zoyne8/f4BIS6NPx6QCF6vQdQD4G8BfjtZTxZHWgGwM3Rx7fCAOHBenP0paiF+9Tm8cdW8EWbAPiKhgDw5QbAzXMQNFyTQQTAJiXcMgB4+PBhHvtYXeXu32krtEHqb4V2bIU2bBVfbJV2bIU+ZQDs3/FlK9bcAFjm1S0DgGUyOLUVsAJWwAqspUCbAKgvHjRxCfgljgCu5da+3W8ALHOdAbBMP6e2AlbACgyMAm0C4OuAnYXKngZ+3gBYqOLmTW4ALPONAbBMP6e2AlbACgyMAgbAgXF1XzTUAFjmJgNgmX5ObQWsgBUYGAXaBMA3NhQB/FlHALds/zQAlrnWAFimn1NbAStgBQZGgTYB8L82BID/0QC4ZfunAbDMtQbAMv2c2gpYASswMAoYAAfG1X3RUANgmZsMgGX6ObUVsAJWYGAUaBMA9fLaJh4C+Q+OAG7Z/mkALHOtAbBMP6e2AlbACgyMAm0C4G8BY4XKzgLPNQAWqrh5kxsAy3xjACzTz6mtgBWwAgOjgAFwYFzdFw01AJa5yQBYpp9TWwErYAUGRoE2AfB3G4oAPtsRwC3bPw2AZa41AJbp59RWwApYgYFRoE0AfGdDAPgjBsAt2z8NgGWuNQCW6efUVsAKWIGBUcAAODCu7ouGGgDL3GQALNPPqa2AFbACA6NAmwD47oYigD/gCOCW7Z8GwDLXGgDL9HNqK2AFrMDAKNAmAL63IQB8qgFwy/ZPA2CZaw2AZfo5tRWwAlZgYBQwAA6Mq/uioQbAMjcZAMv0c2orYAWswMAo0CYAvh/YVajsKeC7HQEsVHHzJjcAlvnGAFimn1NbAStgBQZGgTYB8C8aAsAnGgC3bP80AJa51gBYpp9TWwErYAUGRgED4MC4ui8aagAsc5MBsEw/p7YCVsAKDIwCbQLg4YYigI91BHDL9k8DYJlrDYBl+jm1FbACVmBgFGgTAD/cEAA+0gC4ZfunAbDMtQbAMv2c2gpYASswMAoYAAfG1X3RUANgmZsMgGX6ObUVsAJWYGAUaBMAPw7sLlT2JPAIRwALVdy8yQ2AZb4xAJbp59RWwApYgYFRoE0AvLohAHyYAXDL9k8DYJlrDYBl+jm1FbACVmBgFDAADoyr+6KhBsAyNxkAy/RzaitgBazAwCjQJgD+n4YigP9u+Qjga4EnAPcEdKVYV5xfBNy4gjPvBrwB+HZgP3Ac0OeKfxmYT2l0tfmjKT9tEp/cBtxjYDpIyw01AJYJbgAs08+prYAVsAIDo0CbAHhNQwD4zcsD4GsAfW74M+mTw28F7gc8cAVn3gt4OvAe4Hrg3sD7AD2sfGUGgB8BhoGlgekUF7ChBsAy8Q2AZfo5tRWwAlZgYBTYQgDY67MHAGLOi4Gpmg79KeBHM2hUBFAAOAIs1szDZgUKGAALxAMMgGX6ObUVsAJWYGAUaBMA/xEYL1R2BvjG5SOAvTnr8u/zU2SvbqkfBI4Bz+6JAOoy8miKLr4K+ETdDG23PgUMgOvTq9faAFimn1NbAStgBQZGgTYB8AsNAeB91wbAR6XLud8L/HVNZ/5SAsYHAUdTmgPAQeBzwM60/9XAgwHxrKeGFTAAlglqACzTz6mtgBWwAgOjwBYEwCcC7wSeBXygpiMV1XsmoI+MXLdGGl0S/iTw0pp522wdChgA1yHWMqYGwDL9nNoKWAErMDAKtAmA/7zBCKAew/1Y8ogez/29lSOAPwS8GXga8KGaTrwKeHSCvyM10ughkU8BL6lha5N1KmAAXKdgPeYGwDL9nNoKWAErMDAKtAmA1w/BnsIRfnoJ7tV9HKMDTGeOeiHwSuBJgN45vdY0BLwj3VKoS8a69693egxwLfAVYAfwPECvm9G7qPWAiaeGFSjsHg3Xpv+yMwD2n89cYytgBazABVFgCwHgWeB2YC4JKZbQq1sen4DwEPB54HFpXe//U3BR9gs9aTSOatJl3uekJ4lPp4dABJl6x6Cn86CAAbBMVANgmX5ObQWsgBUYGAXaBMAjo81EAA91Ea83AjgwPtvKDTUAlnnXAFimn1NbAStgBQZGgTYB8OaxZgDwklkD4FbtoAbAMs8aAMv0c2orYAWswMAoYAAcGFf3RUMNgGVuMgCW6efUVsAKWIGBUaBNADy+p5kI4P7uox++BLwFe6kBsMypBsAy/ZzaClgBKzAwCrQJgFN7Yc/2Mmmnz0LnNgNgmYqbN7UBsMw3BsAy/ZzaClgBKzAwChgAB8bVfdFQA2CZmwyAZfo5tRWwAlZgYBRoFQD3NxQBPO4I4FbtoAbAMs8aAMv0c2orYAWswMAo0CoAXtIQAN5sANyqHdQAWOZZA2CZfndJvW3bKxrOUdkN18wz3k8a5pHuomXy0LaVJr3DNJ/0XfM8zyhH26MMLY+ldeWtd6xqUl6xrPWwUR55u6JMzbVdeegXZedt6W2furH2yzbXQOXm9YoytS3af2tKG+9yTXnnxSprfQdAXxQ4k/Zr2ygwkb4JENXbnZosM71/LG9ipI1qzQDx+tm7uEJGqpsmVUY/aaO5EmqKjNQuLesXy6F9bqv0sV93xode4Z/QJPLR/tx3eZ/J+4hsevtMb1/sLWO1Pp23Ybk+2pu3bFbKbznbu+a5tPTyFY8G77hDAQOge8NmUsAAWOYNA2CZfndJbQAMeDMAGgDzw8MA2PCp5oJk1yoA3r2hCOCNlVR+CviC9JjzW6gBsExfA2CZfgZARwC7XwV1BLDnWHAEsOFTy6bIrlUAvLwhALzBALgpOs95qIQBsExUA2CZfgZAA6ABsLpM7EvADZ9KNmV2BsBN6ZaBrZQBsMz1BsAy/QyABkADoAGw4bPI5s2uVQC8N+zRvbcF0/QidL7sCGCBhJs6qQGwzD0GwDL9DIAGQAOgAbDhs8jmza5VALyiIQC81gC4eXtUWc0MgGX6GQDL9DMAGgANgAbAhs8imzc7A+Dm9c0g1swAWOZ1A2CZfgZAA6AB0ADY8Flk82bXKgDer6EI4OcdAdy8PaqsZgbAMv0MgGX6GQANgAZAA2DDZ5HNm12rAPgNDQHgZwyAm7dHldXMAFimnwGwTD8DoAHQAGgAbPgssnmzMwBuXt8MYs0MgGVeNwCW6WcANAAaAA2ADZ9FNm92rQLgAxuKAP5fRwA3b48qq5kBsEw/A2CZfgZAA6AB0ADY8Flk82bXKgB+K+yp+xXMFSSbXoDO3xkAN2+PKquZAbBMPwNgmX4GQAOgAdAA2PBZZPNmZwDcvL4ZxJoZAMu8bgAs088AaAA0ABoAGz6LbN7sWgXABzcUAfxbRwA3b48qq5kBsEw/A2CZfgZAA6AB0ADY8Flk82bXKgA+rCEAvNoAuHl7VFnNDIBl+hkAa+q3bdsralqeD7O6N8Is9BQe6S6Cc6AWJtq20tT7XdedmaHyjHK0PcrQ8lhaV976Pqwm5RXLWg8b5ZG3K8rUXNuVh35Rdt6WqE5sUzfWsmxzDVRuXq8oU9ui/bemtMojm/JilbU+SbUInEk22jYKTGAA3AIAKK8uLb38fBy8WypPA+CWcmffN8YAWOZCA2BN/doDwICfqJggphfWltsm+xy0tJ6D32rA1wtUAVEBewF6vdCWA5dsBXc5LKo+vTAZ9VquzNkegMtBTcsBdDnIKZ+LU2a5fbS9B+xkqbPGUo/j40wS3x4V2GmK5gj2AgBDlv3AQeBEBohKI0hUVecSWyrPYOJwg2xmUhlhF80K+U9mdQwA3ZHyj/qHLKrbuTbl8D2dNJUfAsY1F/jmfay38OhjeZ/K+5Ps836o5dimfHv/Gentn/k/CcsdhJG+t0/neff+Y7RaH18unyi3W5YBcO2TYasA+IiGIoAfr9rVAXQweNpCChgAy5xpAKypnwEwqCcfdA2ABsCALgNgzVNJX5u1CoCPbAgAP2wA7OtOt0rlDYBlnjUA1tTPAGgAdATQEcCap4sta2YA3LKu7cuGGQDL3GYArKmfAdAAaAA0ANY8XWxZs1YB8NGwp86dK6uoPX07dP7aEcCt2iENgGWeNQDW1M8AaAA0ABoAa54utqxZqwD4+IYA8K8MgFu1QxoAyzxrAKypnwHQAGgANADWPF1sWTMD4JZ1bV82zABY5jYDYE39DIAGQAOgAbDm6WLLmrUKgE9qKAL458tGAF8LPAG4J6Dn7fWs8IuAG1dxnt4NcFVKdxb4IPBCYCpL81TgVcA9gBuAlwLv27Id4gI3zABY5gADYE39DIAGQAOgAbDm6WLLmrUKgN/dEAC+f1kAfA3wXuAz6f1VbwXuBzxwFecJ+HQQPD29TOo9wCngu1OaBwMfA54BCDufAvwB8HDgmi3bKS5gwwyAZeIbAGvqZwA0ABoADYA1Txdb1mwLAWCvjx6QIE0vFc0jemEXEb1vBD6bNmr5H1K0T5HD303vG/y+LPM/TW8Kfc6W7RQXsGGDCIB1QtfqmL8BfAswCfw2sNynLAyANTuvAdAAaAA0ANY8XWxZs1YB8KkNRQAV51v7RdC6/Pt84N4rOO/JwLtTtDA30Svdddn3LxJAKir4uszgFwEB4bdu2U5xARs2iAC4Vuh6N3Bt+m/klcAVgJ6D+s/Ar/f4ygBYs/MaAA2ABkADYM3TxZY1axUAnw57RsqknJ6HjpBsdQB8VLpP73uB7ktj7jr9MPAG4NKeXTcDVwJ/BFyXbN6W2QgqtV/jsKeGFRhEAOyVsDd0/az0H8hlgG5U1fQfgZ8E7mMA3FgPNAAaAA2ABsCNnT22TqotCIBPBN4JaNz8wCqeWi0CqAif7g/UfX6OALbY3Q2A3SeX8tD1m4D7Ao/P/PBQ4G/Sf0H5F0YdAazZWQ2ABkADoAGw5uliy5q1CoA/uLEI4OF/Bf00zS/CVV+sFpf7FvAPAW8GngZ8aA2n6R7A6wEFXOIewAi+6EniuAdQY6ouCcfkewDP49Ew6AC4XOj67cCu9CRSSP91wOeAQ8DRzB8GwJqd0wBoADQAGgBrni62rFmrAPjMjQFgLn51CfgdywKgXt+iW6SeBFxd02F6slcHgcBR7PGu9AqZ70np9RTwR9PYq4igngJW6d/mp4BrKrxOs0EGwJVC15siAlgKTEtLL6/dFUrLql1QK4YCrfz7R7f3rKsSy22L7Xkl8wG7zjeVVLamhTQP6NsJGAANgFsbANd7eK/nHLXevDer/RYCQN0epRPpXNJaLLGUrpwJCBUs+TzwuAwQ9R5ARQw19spWQCiQnM78pcvBr07vF9R7AF8MdF9E46lxBQYVAFcLXT8TeD2Q3wP4U6mjLnsP4Ate8AJGRrp32z72sY+tfqVTKZSt5+RaWlZpW++cPiAq3xpAlW/L7bQ/1vNBdrWa9UJiDm6RTuC2HCzmZUc+Yad5Xh/l1Tvw53WNsgJcVWZeF+WnX+SzGpSeTnZ5XpFf3pbYrwC2ljUHglO1rM1DwHjKckcmhZa1T5OekY+0arb2adKjVLlMsl8E9MxfTPuTXdxUofQxpKgpso+8tC/Sajn2x3LIk2Vf5aU8VHa4KThcZWpZQ1Ken8rQduVfdbvwpzbMph1avjUvKRmHn/L+mvuut8/2VjbKi+29fbT3OIi69eYTQmr7cmlie4gRDl/un5yof09zq9W87y93jC6X5o5t6zlHrZ7T5t57+PBh9NM0Pz/PVVfpXcjLXlJtqiHdK1PPbigCqJeznN/6NtVu57NOBQYRANcKXWvo0l0P6vZ6YljQp3D0G9t8CrgUytZzci0ta519bg1zA2BXoBhQDYAVTBoAV4C5laBwORjNATP+Scmj1QbAZs9ld82t1QjgcxsCwN8yAJ7vfnGh8h9EAFwrdC1f3B94S3oPoF5qqbec6/M0vdN5uwewFMoMgGsdUo4AdiN/jgA6AqhjpffWiTh+HAFc60yynv0GwPWoZdvzrcAgAmCTmhoAm1SzyssRQEcAs0u8EsMRwOxSdH7A+RJw46ef85xhqwD4fNgzWtag6Tno/KYjgGUqbt7UBsAy3xgAy/RbJrUB0ABoAOz2Ad8D2Pjp5QJn2CoAvrAhANRjG74H8AL3nPNTvAGwTFcDYJl+BsBzCvghkOpSrB8CyY4JA2Djp5cLnKEB8AI7wMXfSQEDYFmHMACW6WcANACCnwLOnvTO77kzADZ+ernAGbYKgD/dUATw1xwBvMDd5rwVbwAsk9YAWKafAdAAaACs+kCAnwGw8VPKJsqwVQD82YYAUO+/8CXgTdSLmquKAbBMSwNgmX4GQAOgAdAA2PhZZLNmaADcrJ4ZzHoZAMv8bgAs088AaAA0ABoAGz+LbNYMWwXAn28oAvg6RwA3a38qrZcBsExBA2CZfgZAA6AB0ADY+Flks2bYKgC+GPbEl3k2KMj0Gej8igFwg/Jt+mQGwDIXGQDL9DMAGgANgAbAxs8imzVDA+Bm9cxg1ssAWOZ3A2CZfgZAA6AB0ADY+Flks2bYKgC+rKEI4CsdAdys/am0XgbAMgUNgGX6GQANgAZAA2DjZ5HNmmGrAPiKhgDw5QbAzdqfSutlACxT0ABYpp8B0ABoADQANn4W2awZGgA3q2cGs14GwDK/9w0Abtv2irKWtpban4LrSq3PYmjSe+Hi3XB6MbB+oVH+omDZnk5p8q+K7EzbNFc+kV5zdV/Zag7IJLLWfAgYT8l1M3kk17L2aZrMsleV46Zzfb83d6XsF4EzWUfan77z6xdBD+R7AEtPKUtL3dBUP02tAuBrGooAvqRSuANM95PWruvaChgA19ZoNYstCIB1AUyyLGe7XkEDdErzEtDkgKN6BCCtlHfsz1++q/rkgJXDV97mvLy8zb32eZ3iW2c5xEXagLsAtqhDbxrtj/SxHGVoLghMQFe1Y+zOkKcjPuDskrRrAtBH4wV7ArhdaVl2AXSaB7zly7EtlzgkVJ5zPf1hJoGgNisflSEAVD6Rl5pxKkunvNWs3F2Cyeg6IVEuZSwrr5BUeWhZc+WvfNV2Tb11UF20TbYqK6/30nJ9PLSPfVG5WJ/NEmlfwHr4Ll/P90cy7c/7qbZrPf9HIS97pX2yyftUb3/N25b3wbWO67wPrmXbq81a9svvNwCuqFt3XHpdQwD481U5BsCNddNNncoAWOYeA2CZftkAZgDsEo4BsOIaA2DPkWUA7D3VGAANgMXDz4BnYAAs6wAGwDL9DIDn9HME0BHAPKoXEUFHAFc6xRgA1wDAN8Ce+H9yg+fp6dPQ+TlHADco36ZPZgAsc5EBsEw/A6AB0JeAq8ivLwGv91RiAFwDAP9LQwD4MwbA9fbNfrE3AJZ5ygBYpp8B0ABoADQAbugsYgA0AG6o4zjROQUMgGWdwQBYpp8B0ABoADQAbugsYgBcAwB/o6EI4E86ArihDtoHiQyAZU4yAJbpZwA0ABoADYAbOosYANcAwLc0BID/wQC4oQ7aB4kMgGVOMgCW6WcANAAaAA2AGzqLGAANgBvqOE7kS8AN9QEDYLGQfg9gV0I/BeyngNUP4v19fgp4rVOLAXANAPythiKAz3UEcK2+2K/7HQEs85wBsEw/RwAdAXQE0BHADZ1FDIBrAODvwB69B75gmp6Fzo8ZAAsk3NRJDYBl7jEAlulnADQAGgANgBs6ixgADYAb6jhO5EvADfUBA2CxkL4E7EvA/hSc3wO4/hOJAXANAPxvDUUAn+UI4Pp7Z3+kcASwzE8GwDL9HAF0BNARQEcAN3QWMQCuAYB/0BAA/rABcEMdtA8SGQDLnGQALNPPAGgANAAaADd0FjEAGgA31HGcyJeAG+oDFQDCLwA71pXlWievbdtesa78eo1786+f3/Ay5cZl2t5dy9mut9q+BOxLwL4E7EvA6z1vrGy/0rlvrXNuczVYOafp6Wk6nY4M9Gf6PJXZDUy8u6EI4A+sGAF8OvAC4AHA7vQqg7MrtOkQ8HlgKds/AgwBB4FbgUcAH6X7L6EmBahuA+5xnnQa+GwdASzrAlsEAHOQ0+tIYsq3LweB2p8DXG6j11nkU56vtq+0fyWo7E0feUc+sb+3LbGe71c98+165cZyGkTeYd/7SF3ehuXak29T+Xld8/Zoe6zH19tjHq8Fifppu7qd1rWcp8t9ke3T4n5gFBgHJtLvkrRd+7RNp/C0vP3gKcb3zjA6Ml8JPcI8iwwxf3aEhYUhFheGmT8zwvzkOCxsgx0LcHKY6lUuk1m30LqqpbkmzfXT/0v6qUzJP5XsogkaAiJNb5eYARaBmCuN8gjJ1M6YJI/Sq32ax765tK7yVQ+VJ6lVD60rb7Xji2mf1nOXyT66RV5nlatyNKxF2b3tD5eueO5RY/SbTZloWYk0j8Ly18XEvsgw7GUb+3r7W36s9h4jecXyY7z3WMuPwV4xVmpcfkyVnXxXSm0ArE4QU1PvbQgAn1opvRywPhq4GNCJ8e1rAOBy7vqTdLQ8Je0UAH4kbctB8fx0FOdaEbanjStgADQApt5jADQApq5gAFzljGoAbC0CeP4BMPwc4Kb/MlaKAPb2ia8CrgeeCPzPHgBUZFD/4nk6zwoYAMsENgAaAA2AjgB2+0BcuDIAGgBXUKDVS8B/Cnt2lQ1w06eg870rRgBLAPCVwDOA+2Q1DJC8McXqPwO8CvhEWSuceiUFDIBlfcMAaAA0ABoADYB3uR6+0onVEcDWIoAfaAgAn9w4AOpfpK8AbwTelPWUA+l+wM+lGzqeD7waeDDwj2VDtVMvp4ABsKxfGAANgAZAA6AB0ABYayRpNQK4QQA8/GnQT9P8Alz1Z40D4PcDvwfcPT3ksZp2uifwk8BLawlso3UpYABcl1x3MTYAGgANgAZAA6ABsNZI0ioA/mVDEcDvahwA9aSv7v97dg3RPgx8CnhJDVubrFMBA+A6BesxNwAaAA2ABkADoAGw1kjSKgAebggAH7siAG5PT/7q3r2/Su8Y0MMbenXASk/x3g/4LPDvgL/rEe0xwLXp8rCexX8e8FrgYcA1tQS20boUMACuSy5HAO+sgF8Dc4cefgrYTwGn3uCHQFY5q/oewNbuATz/AKiPxOlSbsCeeELL35EifHrv3+OAq7MO8V+BhyQA7O0nusz7nHmuA7QAACAASURBVPRqGb3LSA+B6GGRj5cN0069kgIGwLK+4QigI4COADoC6AigI4C1RpJWI4AfbigC+Mg1LwHXaruNNp8CBsAynxgADYAGQAOgAdAAWGskaRUAPwZ79KLzgmn6JHT+vQGwQMJNndQAWOYeA6AB0ABoADQAGgBrjSQGwFoy2aglBQyAZUIbAA2ABkADoAHQAFhrJGkVAP+moQjgwx0BrOXcPjQyAJY5zQBoADQAGgANgAbAWiNJqwD4vxsCwIcaAGs5tw+NDIBlTjMAGgANgAZAA6ABsNZIYgCsJZONWlLAAFgm9IYBsLfYpaWX32nTtm2vKKvZulLrvRUx6XveMeXbF5bJ0a+BuUMUvwbGr4FJvcGvgVnl7HP+XwOzUuG959h1nSIbMm4VAD/dUATwQY4ANuT+TZeNAbDMJRkArudxq7vC1PkHwIC5HPCWAz1tk01un4PeWoLlJ3gtLweOy+WRw2ZvHWJ9tbJ3ZqOvFlWu8lRblmtz3na9ckrp8zRKq+35JJsoJ2yjfcsNbFF+no/s8jyUf5724lTgeFYfpc99om6nvDXXdqXXcjZpt4pRt9QH4e+WkhwEDgH7gcuBS9Ly7iXYMU9n/22Mj5xkP8eZYJKdzDLKPCPMschwWupuWWCI04yxyBCzaa5tk/MTzNym+sPZU2OgCOHx1BzJdSbrXnrd6yhwDDiZbLR+Ki3LNg4t2YZrI59osvbF/pjnh6S2RRfTsvLNp8hP9dT+yVSff0pzrctGWkZ3ivpGPrLJu5XaoUmvxlVa7Zcr1bYoP7pRtGsmMsj7jJbl4/x40rb4aXtur0zDXtuXOwbzf1jy40Pbo09Fn+vt46pjpFkN6Hr/KcrzqXte6PHTiqvh3JXzHTgAvAb2dA/DDU/TM9D55ip5B5jecEZOuCkVMACWucUAeBf9DIB3pgANlAZAA6AB8A4QXc8/hnVP0AbATKlqXJoyANbtPANrZwAsc70B0ACYFHAEUJFARwBTdxCPOALYc3ZwBLDVS8D/r6EI4AMcASzDhM2b2gBY5hsDoAHQAOhLwN0+4EvAPWcDXwLuPT22CoCfawgAv94AWIYJmze1AbDMNwZAA6AB0ABoAFzzPkNJ5AigAbBswHXqZhUwAJbpaQA0ABoADYAGQANgrZGkVQD8YkMRwK91BLCWc/vQyABY5jQDoAHQAGgANAAaAGuNJK0C4JcaAsCvNgDWcm4fGhkAy5xmADQAGgANgAZAA2CtkcQAWEsmG7WkgAGwTGgDoAHQAGgANAAaAGuNJG0C4OT1sKfn9aC1KpkZTU/DxL0cAVyvbv1ibwAs85QB0ABoADQAGgANgLVGkjYB8MSRZgBwn14c7xdB1/JvvxkZAMs8ZgA0ABoADYAGQANgrZHEAFhLJhu1pIABsExoA6AB0ABoADQAGgBrjSRtAuAtNzUTATxwqSOAtZzbh0YGwDKnGQANgAZAA6AB0ABYayRpEwCP3rKNPXvKhvjp6SUuO7DkS8C1vNt/RmW9o//a23SNDYAGQAOgAdAAaACsNbYYAGvJZKOWFDAAlgltADQAGgANgAZAA2CtkaRNADxy60WNRAAPXVx90q8DTNdqpI36RgEDYJmrGgPAbjX0BfmY4rNJZRW8I3XkfdEyGeblardscnst161PfP9T9lqum261OkSdVtNiZ9oZ+ahcLasty7U58pLNaUDp8zSxPS9TNlFO2Oaft5Jt3t4oX/nHJE3yPLQ91+ziZDie1Ufpc5+o2ylvzbVd6Xve96DdKmY3sAu4W0pyENBTffuBy4FL0vLuJdgxT2f/bYyPnGQ/x5lgkp3MMso8I8yxyHBa6m5ZYIjTjLFoADQAGgBrnajbBMAbpnY0AoCXd84YAGt5t/+MDIBlPtsgAOaF9gLDchVaC9oCApZLG/kHjOTgk+/L00Z5kSaHqjx9DmZ5HVb7CHwvJIVtpM9Bp7dOvZC4XJ1XArZeiIy8AsYCAGM9b0PopH15HfI2937nNAdoAVyAXy8whtaxf18Gdz3+1NG6I+1W9kMZFwr0tE3ziQR1o4A4UusBe1+zxL7Lj3Jo+xEOcYR9nGCMWS7lKOPMMM7JNNfyTAV+wyzeCfa0Z47RCvz0O8ZBTrCPWcaqHCYrdFShcIwDzDNarc/OjjF/ZoT5M6OM7Jhj/vgeGK7uL4Kbt8HJ1F61YzJJdiq1OaSQzW1pZbFHD2kjeaWBtNkLaOwKXaSH0qhq2qZ9OZvLFSpXwKxytf86uvXKbdUNlIdcJrv4X0dp54AoR/WIZaXXTzb6zaS6RHeI/ye0T3VSvnnd7vQPhFa0M/pVLOfnktinglTR/DjL+1zex3r/Yes97vMKrXQshhjLnZN6j5He80TveXG1432tE/dd//FcWnr5WonO+34D4HmX2AWsQwED4DrEWsbUAHguumYA7PaPGEQNgAZAA+Cdz5kGwDYB8MtTuxgvfAhkZnqJe3f034gvAZehwuZMbQAs84sB0ACYepAjgI4AOgJ4RzjUEcDlhpY2AfDaqT2NAOAVnerWP98DWMYKmzK1AbDMLQZAA6AB0JeAu4FfXwLOrocbAA2AZYOrU59/BQyAZRobAA2ABkADoAHw3HnU9wCuNqS0GQH8wtRexvdsLxrhZqbPct9OddOtI4BFSm7OxAbAMr8YAA2ABkADoAHQAFhrJGkTAD87tb8RALx/57gBsJZ3+8/IAFjmMwOgAdAAaAA0ABoAa40kBsBaMtmoJQUMgGVCGwANgAZAA6AB0ABYayRpEwD/39TBRiKAD+gcWykC+HTgBcAD0suXdOPn2VWE+Bjw0PSolNhD74B6EfCbWZqnAq8C7gHcALwUeF8tcW20bgUMgOuW7E4JDIAGQAOgAdAAaACsNZK0CYDXTF3WCAB+c+foSgD4aEBvrh8D3p7eSroaAH4U+ASw0gsZHwwIEp8B/DnwFOAPgIcD19QS2EbrUmAQAXCt/1rUgfXKVt3NHP+l6L+Wzy2jrAHQAGgANAAaAA2AtQbeLQaA0eZHAB+pCYD/C3jZCmL9bnrY5Puy/X8KnACeU0tgG61LgUEEwLX+axEAPhLQfytrTQZAA6AB0ABoADQArjVWVPvbBMBPTx1id+FTwCenz/KgzpGVIoAbAcD7A3o0WdeV/wx4dfrujfJSlO89wOsyMX8REBB+ay2BbbQuBQYRANfqtALAR6X/aNYS0wBoADQAGgANgAbAtcaK1gHwU1P3agQAH9K5vikAfAjwT+mNmd8A/DfgWuAHknj68OIbgLdlYj4fuBK4opbANlqXAgZAKoLJ71vQ8s0pnP2VdIOq7m9YbjIAGgANgAZAA6ABsNbA22YEcKMAePXhWfTTND+/xLuvWvNLIHUvAfdq9O3Ah9IXy/UVbUcAa/Wi5owMgHcFwO8APpk+1a7LxX8IKAyd/1cSHjAAGgANgAZAA6ABsNao3CYAXj311ezeM1SrXisZnZxe5GGdLzUVAVwJADWO6r573QOoZT0JHJPvASzy4OqJDYB3BcBexXTD6mPSk0i9+wyABkADoAHQAGgArDVMtwmAn5i6ohEA/PaOrtIu+yUQ3cunK2iKAP5ViuQtKnCYXvGSa3IAeCCgh0AUXvx64PfTq16elgz1FLDuvddTwB9MTwG/A/g2PwVcq3ut28gAWA8AHws8bBl1EwBql44DTfdNP30cdKVJDxjnk2y1bbhne+TZu11msU/L2h/p83x6v8e5Wn550bmdmtg79Zat9Z2ZUbT9dE89l9NENr16RJsiyzyf0CLSaL23XbLXpLpHvWIe6bWe1zu1M8z0/6jeUpVPUVTsC9krmyhTy2qn8tY2GcW+yFzroYXKjc9naX+qR5QVybVLy7vTXP/Y66KJ1nelZo4CE+k0fHdgf9q/A9D618xxz8uuZz8nuJwbOMQRDnKMA9xSzceZ4QDHGOM0I8wxPjdTtWx07ixzo9tZHB5m7NQ8s7tGGJ2bZ2ZsN7OMcZTLOM1ObmOiymmGcU6wjwWGqnXZqJSbuKzK7+h0d3765r1ddU+mYUCSaPhQ2zSPfdEdJJO0j+4gPQ4CykZt1HrMlWZ3cuDCNtixAMOLsDDEyO5ZFm4f4uyJXdBt4h3lnkp6apuW9SYy1S/vouGHqE/sV17yifygaTzpr3YpL9nJRsv6bnDeVdQuTZrLTj/ZRLla13SnPqmM4/jRZTr98mNT+1Vg9LW8v/cei3GuyPtxKvNcGRq3lztWoyGa5/vjOIj9sS/Od3m6KGu1ed62Xrv8PNBbn+4xv7SkW8nanw4fPox+mubn57nqqqu0eD4/rVaNSy0A4LOA38t6Zbw1Q1fRdOPg54HHAVen9/r9cbqXT0e4brP6k56HQKSLHvjQgyH3TEffi4H3t++1wShxEAFwtf9avim9+uUz6b5APQ38rvTYenXU9kwJAF+Tzvb5XgPgHSN16GIA7CphAKxUMADe8b+CAbDGiNufAJg3rM0I4Eem7ttIBPA7O18438Baw/c2OR8KDCIArvZfi4Du9Sleon9X9RDIW4DfXkF8A2AVdnEE0BFARwCrYJsjgNmp0hHA3nGjTQD80NQ3sKvwHsBT04s8qqN4yHmNWJ4PtnGeNRQYRACsIUttEwOgATB1Fl8C9iVgA6AvAa8+dhgAa4+tNmxBAQNgmcgGQAOgAdD3AHa5xxFA3wO4xnjSJgAenvqmRiKAj+38gyOAZZywaVMbAMtcYwA0ABoADYAGwHMPmvghkNWGlDYB8C+nvqURAPyuzt8bAMs4YdOmNgCWucYAaAA0ABoADYAGwFojiQGwlkw2akkBA2CZ0AZAA6AB0ABoADQA1hpJ2gTAP596ELv2LPcKsVpVrYxOTS/wpM6nHQGsL1lfWRoAy9xlADQAGgANgAZAA2CtkaRNAHzf1EMaAcDv6XzKAFjLu/1nZAAs85kB0ABoADQAGgANgLVGEgNgLZls1JICBsAyoQ2ABkADoAHQAGgArDWStAmA7516eCMRwKd2/sYRwFre7T8jA2CZzwyABkADoAHQAGgArDWStAmA/33qEYwV3gM4O73A93c+bgCs5d3+MzIAlvnMAGgANAAaAA2ABsBaI4kBsJZMNmpJAQNgmdAGQAOgAdAAaAA0ANYaSdoEwHdNfWcjEcBndD7iCGAt7/afkQGwzGcGQAOgAdAAaAA0ANYaSdoEwHdOPYqxPfpO+8an2enb+ZHOhwyAG5dwU6c0AJa5JwHgW4CLU07xJnx9CF3fh9UBmL+LSeu9H0mXndLtXKY2sl1uym0jT22LfHr3Kw/VI7fVtuXqE9s1Vz6RLvJQGVEv7VceUV6syzZs1joJhR7SIb6pq3wiT23L2xPa5l8dCI1V1nTKR9vGknjj3XqG2UTKcnfavR/Q8q60PgeMpm+7alMUdSbt13xH2t5bjZNpXy6XlmUf6WWT22m7yo/6yDZPo+2qu7ar7pqUZgi4G3BJ2qZ9asvXLXHg3ke4jKPs5DSH+BcOcgsHOMYhutur+eJNDC0ssOufzsIJ4NbUbpUVbVQZWj+VyggN1fs7cOo+2zkxup9jHOAWDnK0yv0yZtjNJHtZZIgbuLzaNscIN916GfNnRmFytJuvdLg56+Papp90PZ7pEpqFP2SjtnYAuXfHEgwvMrJ7loXbhxgbP838mZFqeXzvDMPDi1UhO7efruo0OT3B6clxODncLW/3HBwfhYnk0IUhuG5bVwd1QSXXcmWburf6ieqoSXPVMfwu3VSvSBNdW2m0rHlMKlL55uVMZmVov/LODz35Q4eC6rWkjJRYG/JzkI6FOP7ieAhxs/KrRaXNj9nIL47LmC933OVl5sdzXkbULcrIz49x7srPlSFKfv7o3b/aOSbXApaWfry3wa2vGwBbl9wFrqKAAbCsexgAKzAzAJ7rRgE0AXD5WGsANACqoxgA0+FiACwbflZMXY1Lvz/1uEYigD/a+R8qSP9m6b8JT1tIAQNgmTMNgAbAOwIuEflxBNARQIG/I4A9Z1dHANuMAP7O1BMaAcAf63zQAFjGCZs2tQGwzDUGQAOgAdCXgLtnEV8Czi4j66pAXPOOk6wB0ABYNuA6dbMKGADL9DQAGgANgAZAA+Cd7jv0PYArDSttAuBvTT25kQjgczsfcASwjBM2bWoDYJlrDIAGQAOgAdAAaACsNZK0CYBvnfoedhY+BXx6+nZ+ovM+A2At7/afkQGwzGcGQAOgAdAAaAA0ANYaSQyAtWSyUUsKGADLhDYAGgANgAZAA6ABsNZI0iYAvnnqaY1EAF/Y+WNHAGt5t/+MDIBlPjMAGgANgAZAA6ABsNZI0iYA/vrU09m5Z6RWvVYyOj09z0913mMALFJx8yY2AJb5xgBoADQAGgANgAbAWiOJAbCWTDZqSQEDYJnQBkADoAHQAGgANADWGknaBMA3Tf1gIxHAKzt/5AhgLe/2n5EBsMxnBkADoAHQAGgANADWGknaBMA3TP1IIwD4c513GgBrebf/jAyAZT4zABoADYAGQAOgAbDWSGIArCWTjVpSwABYJrQB0ABoADQAGgANgLVGkjYB8FenfpQdhQ+BnJme5xc6v+8IYC3v9p+RAbDMZwZAA6AB0ABoADQA1hpJ2gTA10z9eCMA+JLO2w2Atbzbf0YGwDKfGQANgAZAA6AB0ABYayQxANaSyUYtKWAALBM6AaBukr0kfQhdGS6kXOObmDvTPs016UPpmobTB9NnsvW8Qson8tA8JqVXXkofeWlf5D+eLWt/Xp6WI8/e9JG/PuKufSo/8kz71GOG0vIisDsta64ko6m4HVl1Za9slJ1s9NM2VeVk+sn+TLIJ+WSnaTLLV1Wfy/KKcjVXHlE35a31+EV1wv5gJqHK6QCSTftVJ9VFy5qHq8K1kkd1iPV8v/JSmqh7rMe2kPV4Kif0CPvQM9dA+6Id2j+R9JD+au9eOHDFvzDEIgc5xv34PJdylK/hSxzgGOPMsJ8TTDDJZbM3MXoMmAY0vzW1/cuA8rsMUK/eldY1P5D5LfTYAaf2bOf06BhHOMRx9jHJXk6wj2Mc5AYuZ5ad1TbNb+Eg1x37as7OjcKNw12f6qf8Yq62qA6a8v4R7ZdfZCsttE2HnH7SQ/Ws/LfA9tGuc/YfPMHs7BijO+aYPTnG2O5ZFhaGGB2Zr7afPD4Bk8MgX4QPw+/RV1WmlrU9flqP/q1tSq9tNyefqn/EFIdP+FPptKyf8tbxojb39vnQOfp79DnNw+9qZvQv2Z1KdZR2kd+5Y72nTtEP4xQUp5M4Fen4UR2WIp0K1u9EyjwaGRmpQ8W5JT//KX2cf1TJSKd573ko1qPM/JwX+Sx3PszL6C3vjjyWlp6eiXBhFtsEwFdOPa+RCODLOm+TWDpDysn5JEFfADwgnS3l6LMrKHs34A3At+vQTEfNu4FfBuZTmkcAH01HhjZptLkNuMeF8dbWL9UAWOZjA6AB8A54yCHCAGgAVH8wABoAszGmTQB8+dRPsGOP/sPY+HRmeo5XdN66EgA+GrgYGAN0nXg1ALwXIGDUW6WvB+4N6CPDHwauzADwI+nfmnP/emy89k65lgIGwLUUWn2/AdAAaAB0BNARwCo45AjgWsPJFgPAaK4idwK31QBwOWl+CvhR4IE9AKjPl8S1gLUk9f4CBQyABeKli2VT4EvA1eVSXwK+45KcI4COADoCmM6uvgQcw0ybAPiSqRc2EgF8TefNK0UASwHwg+lGlGf3AOCN6eaIzwCvAj5RNkw79UoKGADL+oYjgI4AOgLoCKAjgI4A1hpJ2gTAX5j66UYA8Fc7v3Y+APCXgOcDDwKOJvF0t7Huzv5cunlC+18NPBj4x1oC22hdChgA1yXXXYwNgAZAA6AB0ABoAKw1kvQDAF53+Mtcd1i36cHi/CKfvur/Ng2Aiuo9E3gkcN0awunS8ieBl9YS2EbrUsAAuC65DIDVc1l+CthPAfsp4O7JwE8Bp6eCfQ9gnaGkTQB80dTPMlr4EMjc9Byv77yxSQC8CtDDI4K/IzU000MinwJeUsPWJutUwAC4TsF6zB0BdATQEUBHAB0BdASw1kjSJgD+7NSLGgHAN3ZevxIAbk8PfughkL9KL2LSwxt6rUvvU7wKG7wD+EbgUenev17NHgNcC3wlvWjpecBrgYcB19QS2EbrUsAAuC65HAF0BDC9yszvAfR7AB0BzN4X6ghgnaFkiwHgs4Dfy2BPPCHw+470qpfPA48Drk7v/9M7/uItrpIr7BVI0aTLvM9Jr5bRiyL1EMgrgY/X0dY261fAALh+zfIUjgA6AugIoCOAjgA6AlhrJGkTAH966hcY3ZO/kb9WFe9kNDd9hl/r/Opal4DXn7FTbAoFDIBlbjAAGgANgAZAA6ABsNZI0iYAvnDqJY0A4Js7rzEA1vJu/xkZAMt8ZgA0ABoADYAGQANgrZHEAFhLJhu1pIABsExoA6AB0ABoADQAGgBrjSRtAuBPTL28kQjgWzuvcASwlnf7z8gAWOYzA6AB0ABoADQAGgBrjSRtAuDzpl7JSOE9gPPTZ3hb52UGwFre7T8jA2CZzwyABkADoAHQAGgArDWSGABryWSjlhQwAJYJnQDw/cAlWU4L6fUIsUlm+jisvokZc+3T+mx6r0jYxgfVY11Pw8e2eFo+8tD2nVm52XJ8izaKVJU05d+o1ZuZtD8eFNO6HtLPi4xqxnd+VYSW87RaHk/5aJ+gcH9aV97apnkuS5SpF+nqp3XNT6a5bG9PnwSfSfnH94ajLUqjslUnlaf9UYaW8wfglHfUbTds75xibPw0Q8MLLC4MMzY2yz5OMMYsiwwxxyj7OM5wtdb9zTPCcfZV27Rf6yPMM8QCi8lyhDlOsL9KoTSa5s+OMLR9kcWzQywsDFXlLS4McXpyHBaGGJlQA7vTzt2z5/aP7Ih3zcDEyCR7mazKU726ruzWS/W+nOsZ4zQHOcY38BkOcYS9nz0NU8CtwAmoxmjpoG2q2qm0vi/tV6b3Bg4B+iiTpnvA3D1gZqxL+iov9DnBPmYY5yiXoeVjHGSSiWpZ27RP67OMcWJ6H6e/sLdbpr70KT9PAreleQiQv15nV+qPo8n/sg//yrc65PRT1aIPalk/9R31bZWj9Gpv9Af1A/VxaaI8pUnepyYy27w/qo43pbziZegqQ3noJ7doHlPUQeUpf+WrXxyDylt1VV1666c6a5Jeykdz2Ud9Is+Yqw3qRvGSjWhT6JW3L9dB22UbuuVaqF5RvupXvdlNGzWpM1U9Nm2L9Thw8/NYFJ6n07bIKz8n5g3LtDy3GPnmeXaPhuWnKAOWlh6/gk17m9sEwB+fek0jEcC3d6p3MHcyp7cnmEs6rwoYAMvkNQAGhBkAKzAyABoADYBBvpprMgBGn2gTAH906lcbAcDf7/yCAbCMEzZtagNgmWsMgAZARwAdAXQE8FxE0BHA1YYUA2DZgOvUzSpgACzT0wBoADQAGgANgAbAWiNJmwD4I1NvYGRPfotQrSreyWh++jTv7PycI4Drl64vUhgAy9xkADQAGgANgAZAA2CtkaRNAPzBqTc1AoB/1LnSAFjLu/1nZAAs85kB0ABoADQAGgANgLVGEgNgLZls1JICBsAyoQ2ABkADoAHQAGgArDWStAmAT5/69UYigO/p/JQjgLW8239GBsAynxkADYAGQAOgAdAAWGskaRMAnzb1Zi4qvAfw9unT/HHnhQbAWt7tPyMDYJnPDIAGQAOgAdAAaACsNZIYAGvJZKOWFDAAlgltADQAGgANgAZAA2CtkaRNAPyeqbc2EgF8X+cnHAGs5d3+MzIAlvnMAGgANAAaAA2ABsBaI0mbAPjkqd/ioj1jteq1ktHt07N8oPNcA2CRips3sQGwzDcGQAOgAdAAaAA0ANYaSQyAtWSyUUsKGADLhDYAGgANgAZAA6ABsNZI0iYAPmHqdxqJAH6w82OOANbybv8ZGQDLfGYANAAaAA2ABkADYK2RpE0AfNzU7zcCgP+j86MGwFre7T8jA2CZzwyABkADoAHQAGgArDWSGABryWSjlhQwAJYJbQA0ABoADYAGQANgrZGkTQB81NQ7G4kAfqjzI44A1vJu/xkZAMt8lgDw/wATwEXAcMpRy5o0H+/ZllZD/aXeStyebVhYJk8g0g4B+t73aCom1lWd3SmbqJKy0hTQpnTap3UVqapqfirN59L+PE/Z965rm8qKvFT2Jcp3CYYXq9/IjjlGdsxXxQ8NLzA8vMjCgjKC+TOjDA0vVtu1fHpyvEqzfWiB4YsWWbh9iNGd81Uemk6fHKu2jY2fZlT5bp9nhDnGOM0YsywyxBCL1XyWsWp5iAWGWWSCScaZqX6aRpiv9mqKfVUdWWQfx8/5QWnnGK1SRlotjzJ3l7KUaLLqD3LLXJVO5c+n9Kqj8ld61VzLylN1nWVnVc+YFhiq0qluBzmG1vdzorJX3ZViHye4nOur7fv/7SR8Afgy8FngKDCVdSf5Wv69OPlL/eZA6gMyuz9wX2AfnNqznZtGL+MIh6q6zTHCIsPn2nI0gZ/aeoyDSdVxbuEA1x+7nLMndlHJrObckH4ngRuB6NY3pf1dF0D0OdVzV+qTMdd+LetwUr0vTf1Otup7+u1PbTmT8judjg+Vp22ax/IkoPrEds1Vj70prey63aRbL02y0TESLlLZykfb1V1UXrQh0iqd6hxtimWVrfpqCj1Uvo7DTlZv5Stb5a1ylW8cq3GaUR6yUV2UV28dpY3Sa1Ka0Eflqb7qrvpFWm2TfX4u0D6VoUn7tK7zUHX+0gb94twVFdS2OAHJTpWXTZyMtKwTUUx5g7RNaUOc2Nd7fgy7LJtzAkVFYWnpm3ODC7LcJgB+59S7GC58CnhhepaPdJ4hrdQjpy+IaC70vClgACyT1gAYA68B8BxsqksZAA2A5+DRAJidZQ2AnY5Y6rwCVTUuGQDLBvdBSG0ALPOyAdAA6AigI4COADoCWGskpZ4C3wAAIABJREFUaTMC+Iip/95IBPDjne8/38BaSzsbNa+AAbBMUwOgAdAAaAA0ABoAa40kbQLgw6fey/Ae3TOx8Wlh+hR/03mqAXDjEm7qlAbAMvcYAA2ABkADoAHQAFhrJDEA1pLJRi0pYAAsE9oAaAA0ABoADYAGwFojSZsA+JCp9zUSAfxU53scAazl3f4zMgCW+cwAaAA0ABoADYAGwFojSZsA+KCpP28EAD/deZIBsJZ3+8/IAFjmMwOgAdAAaAA0ABoAa40kWwwAnw68AHhAeiGT3tNzdhUh9KKhq4AnJLsPAi/seVGVbjh8FXCP9OKolwLvqyWujdatgAFw3ZLdKYEB0ABoADQAGgANgLVGkjYB8Fum/pKhwodAFqdP8fed71opAvjo9EbRMeDt6Q2TqwGggE+QKHAUe7wnvWnyu5N4DwY+BujFg38OPAX4A+DhwDW1BLbRuhQYRABc67+WbwR+A/iW9FrV3wZesYKqBkADoAHQAGgANADWGnjbBMBvmjrcCAD+Q+exa10CfgTwkTUAMCJ6Gl/1inpNWv6HFO3T6+F/N70f8fsyMf8UOAE8p5bANlqXAoMIgKv916LXGV+bOuIrgSuAvwL+M/DryyhrADQAGgANgAZAA2CtgXeAAfDJwLsBRQvzSd+k0WXfv0hRPkUFX5cZ/CIgIPzWWgLbaF0KDCIAhkDL/dfyrNT5LsvuZfiPwE8C9zEArvBpOH8K7k6fnVM/8ZdA/CUQfwlEn2zzp+DycaNNAPyGqQ81EgH8TOdRTUQAfxh4Q/qAYy7JzcCVwB8B1yWbt2UGz0/7FYzx1LACBsDuPQlx38Kb0pdQH5/p/FDgb1JoOr6EGbsdAXQE0BFARwAdAXQEsNbQ3CYA3nfqIwztiQ/C16reXYwWp0/yhc53NgGAq0UAFeHT/YG6z88RwI25akOpDIB3BkDdyKpXp1dfv07T1wGfAw4BR3tUNgAaAA2ABkADoAGw1gDcDwA4c/h/c/Lwp6r2LM3fzq1XvbcJANQ9gNenJ4bjHkA9PSzouycQ9wBqTK0+PZIm3wNYq2dtzMgA2EgEULw4CgwB/z79FFjUpPl4tpw5KtRf6nWeLp3EtJCuu0ZeaXukVZE7U/EqJtb1wH388ydI06SsNO3I7LRP6ypSVdX8VJrPrXDJNy/DAGgANAAaAA2AK47Ahw8fRj9N8/PzXHWV3oRCB5je2LC9ZqoqMHHF1CcaiQBe2/n2leq7PY0aup1K98prBFpUM8WOy9RST/ZqlPmh9BTwuwBdVaveNA3oKeCPpgCMIoJ6CvgdwLf5KeA1fb4hAwPgnQHwmcDrgfwewJ9K7ypa5R7Af05PwwdFiYrUz2OeQZvgKQc1meg2WE0592ldh1JMkS7nwYA55SH+FPRpUv5aFtgF7AWPar/AMOBQc8U8VZaWdcuOwE+T6hXQGPlFPVWG7AMgtV3LEwtsH53j4MFbmGCS0ZTZCPOMMYvmQynT04xV8KSf7GJ5hnFm2clodR6BEeaqZc2HWWSBoWou+5hUlrZorjJi0r14Km+R4aqMnZxmHJUwcy4/5R6Ttque+u3jRGUzz2iVt8rVssrQ/uPsq+qlNEqhXGWneqlc2exkttqv+h5nf7Vd67KZYTeT7D1XnuoQZUR9Ik/V5TKOVukOLR5hz7/Md/+fviXFpb8A3Npdvv0LcMsUaNM/QnUr2sXAVwEHk5u+uQMX6X/yA8A3pefw1A8eDKfus53J0b3cwgGOcKj6HeNg1T7VX23UuvZ3fTXGsWMHOHtil25+7N7J80/Av6VWyE3aropoHv9oaK71uIc0+r/+odEkt0Q/Uz9UX44+l/+Do3SRRjYaWlWm+q9+2hb/BIWwyk+/6OOyiWMr/8dJw5O2q+76xyiOh6ir2iAbteN4Wo4bRZRG9VDecfzoOL40tU3HmdoYx6b013LYK43KU95Rz9BR7Yh//KIOslXZOoZVbtRVbddyHNtqX36uiWN/f7KLdCoj2hJaqfwoJ+yi/eeGfBUUJwotq4KqbNwTqHkuZDhPafITVTgrzqvnDtPs/BplRaeKPNRo7bujoUtLd88zuCDLbUYAv3rq6kYA8Eudh0mr5YBV98z/XgZ74gn1gu9I0b7PA48Drk5iq0e/GXhishMQ6j2AOQjrcvCrU1TwBuDFwPsviLMGoNBBBMDV/mvRKfiL6Sng16QHP/SfyBtXfwrYAFgdKwZAA6AB8A5oMgAmCDMABktsMQAcAETa2k0cRABc7b+WTwD3B96S3gM4Bbw1vZl8uZ6Q7gE0ABoAHQF0BNARwO5J0hHAlbChTQC819Sn2F74EMjZ6ZNc33nIShHArU1HA9C6QQTAJt1qAPQlYF8C9iVgXwJW5NeXgNccW9oEwENTn24EAI90HmQAXNOz/WlgACzzmwHQAGgANAAaAA2AtUYSA2AtmWzUkgIGwDKhDYAGQAOgAdAAaACsNZK0CYCXTV3D9j3xBopa1buL0dnpGY52vtkRwI3Jt+lTGQDLXGQANAAaAA2ABkADYK2RpE0APDj1/xoBwGMdva7vvL62ppZ2NmpeAQNgmaYGQAOgAdAAaAA0ANYaSQyAtWSyUUsKGADLhDYAGgANgAZAA6ABsNZI0iYA7p/6bCMRwOMdvRjDEcBaDu4zIwNgmcMMgAZAA6AB0ABoAKw1krQJgHtv+0IjAHjb3vsaAGt5t/+MDIBlPjMAGgANgAZAA6ABsNZIYgCsJZONWlLAAFgmtAHQAGgANAAaAA2AtUaSNgFwz79dy7bCp4CXpmeYvtsVjgDW8m7/GRkAy3xmADQAGgANgAZAA2CtkaRNANx105cbAcBTl97bAFjLu/1nZAAs85kB0ABoADQAGgANgLVGEgNgLZls1JICBsAyoQ2ABkADoAHQAGgArDWStAmAO47cwLY9GqI2Pi1NT3Pm0OWOAG5cwk2d0gBY5p4EgJ8H4o3rY8DO9EuZS2V9J1PzobRL4LQbGAUuSna3A3M9FZJdpJG9bJS9pgVgR9qvPCYA2StPVUfptD2HNKWTndJp2gXsBc6k+miufCN/lbeY7FW+9kXdlY+mHUuwsI2R/dOM7JhnfGyGQxxhgknGmGWIRYZYqNYXGWY0NXKGcUaYrwBK27Ss6QT7UvFKOVztVw5Kry13LC9U+zXt4zjDLLKPE9X+EeYYZZ7j7Kvms4xV2zVpn2yVV7fU0VTDIfYyWW1RfkqnsjXJNvIcnz7N8C1JOxW/ANMHRpgf6lqML84wOTRRtX1i+iTDc7A0BJMX76zKU9njsycZWoDZXSPsuXW+q/8wzHXgtrFOVabaNs/Iufbuv/5k1xcngCngCKB6/Euafwamj8B1Z7rdRK46BpwG9C7/y3fBzoPApYCu6nwTcB/gHnDzAzrcwkG+yNdyA5czyUSl3REOVduPcZAbjx6CMyNwchucBI4DN6a+IwDQT4WdAv4t9WWtq86qjJbzvq5+rG2aq8/KLulQCbA/9dU4JOQ+9bmwV17qx3EcaK4+mnxS1TH6u/ZJlJiUVmVFeUoje0mv4yH6+c2pXVrX/jg2lE+0R3mEXRzDYad9qqOOM5Wp/XdLdVR6HadxLMZxpXW1I46v0FTb1R3zc0R+rGt7tDM0j/Yqj0gXWstWv6ij9M6P+24n7J5HNCl9aKo0sRx9QfNz3wOOgkP0cPx05oTcKdFjs93nFvM8wk5ixYlQhspfjdR2LWvfnZ2+tKTOf2GnNgHwouuPNAKAt9/rkETTkZE778IK6dIbUcAAWCajAVD6GQANgAbALgAbAHvOqAbAXBADYNmA69TNKmAALNPTAGgAdARQ0SADoAHQEcA1R5M2AXDbPx9l23jhJeCZaZbuc5kjgGt6tj8NDIBlfjMAGgANgAZAXwLW7QAGwDVHkzYBkH+6BQoBkJlp+LoDBsA1PdufBgbAMr8ZAA2ABkADoAHQAFhrJDEA1pLJRi0pYAAsE9oAaAA0ABoADYAGwFojSasA+NkTzUQA7189lOeHQGp5uL+MDIBl/jIAGgANgAZAA6ABsNZI0ioA/sNkMwD4TdXj6AbAWh7uLyMDYJm/DIAGQAOgAdAAaACsNZIYAGvJZKOWFDAAlgltADQAGgANgAZAA2CtkaRVAPz7Kdhd9hQwJ6fhW6r3kjoCWMvD/WVkACzzlwHQAGgANAAaAA2AtUaSVgHwbxsCwAcbAGs5tw+NDIBlTjMAGgANgAZAA6ABsNZIYgCsJZONWlLAAFgmtAHQAGgANAAaAA2AtUaSVgHw6oYigA9zBLCWc/vQyABY5jQDoAHQAGgANAAaAGuNJK0C4McbAsBHGABrObcPjQyAZU4zABoADYAGQAOgAbDWSGIArCWTjVpSwABYJrQB0ABoADQAGgANgLVGklYB8MNTsKvwKeBT0/BIRwBrObcPjQyAZU4zABoADYAGQAOgAbDWSNIqAP7PhgDwMQbAWs7tQyMDYJnTEgDq7DeWctrZnUvZHWnTMBUkoIEyJu3bBch8NG1cBGaAi9I+zZVWk17GHnbafjsge01DKc3+tK6yIm9tUjqVo+3KT3ntTsuqx+XpQ+7aJhtNUddY177IS3Plp7xU370p3e4Fto/Osfduk+zbfoIJJhlhjlHmGWem+o0wzxALDLPILGPVNtlpPsYsCwwxyd7KZozTqci5KvVOZqt0shnnJIsMVdvnGWGIxaqs/Zyo8o1tx9lX5a1J+cmuK9lCVS9Nyke17Mp8W2WvfFRe1EH5xXTp3FF2HTnb1Woutf3irh5Lu2DbqTRXUSfST8sqQl9VUrop4FTy063Jh7Eu3+W+1fItwL8kzcM3Sqd8tF02X0h5Rj+pztspzdcB90i/y4BDsHBfuG7PPbmFg/wtD65ae4J9HOWyStdJJriBy5mdHePk5DgcGe22dzL9Tqa5yrgZKnepHuqPMan9qq/6a/Qb1S/aJy3UH5VW7VeeXbfD3dIxpD4a6bVN/S7amPdF5SVbla/8lZfWqw8ZZPWLuoWd5nF8qJ8rjfLS74asTqqfpqiLtNAku3/LtFDZkV62cTzLNvRReapfHIeyj+NU7dP2qIu6b6RTmdInzgGaX5LVQ4vKN47baGvvttBYbYo6qt36qfzYH+l7zznKX3lqkn+13D3MsgyVcW9FtC0aECfGOMmFbQisvPKGaj1s831RbpzgtH7XyNfSUthduPkWBMBXAD+eBP974AXA55ZR+BDweSD3gk6qOvoOAjqbPQL4aNazNIrels5aF85pW7hkA2CZcw2ABsDuGGcANAAKjAPyDIAJBA2A+RDTKgD+ZUMRwO9aMQL4c8ALgccDXwJeDjwTuAKYrTG0/kki+qckWwHgR9K2TYDrNVrQ5yYGwDIHGgANgAZAHUOOAHYjowZARwBXGVNaBcAPNASAT14RAL8MvAl4c2qyonlHgSuBP1xjaP0q4HrgicD/7AFARQbjGkHZCO3UqypgACzrIAZAA6AB0ADoS8C+BFxrJNlCAKixT//yPBT426zxh4HPAP9pDUFeCTwDuE9mFxHAG9MNEcrnVcAnaolro3UrYABct2R3SmAANAAaAA2ABkADYK2RpFUA/NOGIoDfu2wE8O7p7uP7Al/MGv9uYBp47iqC6EbOrwBvTBHEMD2Q7gfUPYS6S/b5wKuBBwP/WEtgG61LAQPguuS6i7EB0ABoADQAGgANgLVGklYB8I+nYKzwNTCz0/C0ZQGwJAL4/cDvAYJIPeSx2qR7Aj8JvLSWwDZalwIGwHXJZQA89/SlnwL2U8B+Crh7QvBTwH4KuOY40hcA+PeH4RpdydXtnPPwwau0JApUZC+flrsH8CbgZ9a4B1BP+ur+v2fXkO3DwKeAl9Swtck6FTAArlOwHnNHAB0BdATQEUBHAB0BrDWStAqA72ooAviMFR8C0X1+egr4CYBg8GXADwNfu8pTwPcDPgv8O+DvekR7DHBtujysFxI9D3gt8DDgmloC22hdChgA1yWXI4COAILfA+j3AFZnAr8HsHtC9HsAa48irQLgOxsCwB9Z9UXQv5xAbTwBXbwHMN779zjg6kyg/wo8JAFgr266zPuc9FItvSxSD4HoYZGP1xbYhutSwAC4LrkMgAZAA6BfBJ3OAwZAA+A6x48tCIDrVMDmm0kBA2CZN3wJ2JeAfQnYl4B9CdiXgGuNJK0C4O82FAF8tj8FV8u5fWhkACxzmgHQAGgANAAaAA2AtUaSVgHwtxsCwOcYAGs5tw+NDIBlTjMAGgANgAZAA6ABsNZIYgCsJZONWlLAAFgmtAHQAGgANAAaAA2AtUaSVgHwN6dgZ+F7AE9Pw/MdAazl3D40MgCWOc0AaAA0ABoADYAGwFojSasA+OaGAPCFBsBazu1DIwNgmdMyALwo5aS53pLcO90O6Ml2TXpzrGx2gjywO22SiXZp0m69CUn7NE1k+Y2mT2Xr9Qv6qI7WVaxstK6fPssdL2tWPrtSGm3XA/vKV3bad0nKW+sqXz/VRZ/jjnXZ5VOsqw7KU7aa7015T8DIJdMMDS8ysmOO8ZGTjDPDKHOMMM8QiywyxAhz7OcEE0wyxiw7mWWe0Wq/7DXJ5jRj1bLSaNrNDIsMM8NuJtnLPCPVvks5Wm0fYoFR5plhvMpbeexlsspXk+b6RX20TaUrjSaVrTSqi5ZV6+72k0xwG4fmjrBr+uw5fy3sguFTsDQK225JQskvp9Ln0VWs1sO/0lr7NJ9K2zWXrhdnn0KXvkqj16tqv9aVTpPWtaynUbX91rQsifRRpX3J75rfP3118zJYuAyO7jnAEQ5xHV/NDdyLa7mi0uo2Jiqtj7OPE7P7OXnD/m75GuCPpPxjsFf/0nK8EFptURfXdvUvtUVfC402q/3aJ5v41Lvs1G/DTnPlF30++rq6gvr4/iw/2aifqj9rUjl5/5c2yltpVE7km9tpOT/G4phSfrLXV0mjjSpLU+ivdbVHNqqf6q46yj4O9TgutF156xf1lr3SRxvzuki7+KUkVdePJ49jm9LeLdlGu/IyQud8m+qiuqsdUSeta7vOE6pT+CzarPKUV56P6pLbVctqeCSKndoux4STovL5XGnydHEyinnkofNA5JU7PfIKEe4a+VpaWq7cdrcZANvV26WtroABsKyHGAANgN1hzQDYHZcNgAZAA+CKo0qrAPhrDUUAf9oRwDJM2LypDYBlvjEAGgANgI4AOgJ47uqFI4CrDSmtAuB/bggA/5MBsAwTNm9qA2CZbwyABkADoAHQAGgArDWSGABryWSjlhQwAJYJbQA0ABoADYAGQANgrZGkVQB83RTsKHwK+Mw0/LwjgLWc24dGBsAypxkADYAGQAOgAdAAWGskaRUAX9MQAL7EAFjLuX1oZAAsc5oB0ABoADQAGgANgLVGEgNgLZls1JICBsAyoQ2ABkADoAHQAGgArDWStAqAv9xQBPCXHQGs5dw+NDIAljnNAGgANAAaAA2ABsBaI0mrAPhLDQHgqwyAtZzbh0YGwDKnGQANgAZAA6AB0ABYayQxANaSyUYtKWAALBPaAGgANAAaAA2ABsBaI0mrAPiLDUUAX+sIYC3n9qGRAbDMaQZAA6AB0ABoADQA1hpJWgXAF03BaOFrYOam4fUGwFrO7UMjA2CZ0wyABkADoAHQAGgArDWSGABryWSjlhQwAJYJbQA0ABoADYAGQANgrZGkVQC8sqEI4JscAazl3D40MgCWOc0AaAA0ABoADYAGwFojSasA+JMNAeBvGABrObcPjQyAZU5LAPjPwA5gIeW2E9AZ8aI012Z9JF3740w5DCi5bGWn9bDT+njanlUwvKWilM1itk/b9NM0lG1XNrEvtqvI3anI/cBESqNiNUVVoroqR8varjyi3BlgMmtm5LsrVV/5atJ8NG1TuVqWrfIbhp2X3MbEnkmGWGSM2Wo+xAJjlWYqrrtFf7vNW+Qgx5hjlEkmOME+Rphnlp1VmlHmKptxVEGqfcpP64sMV3lpUlk7OV2ta3mhSjVcLXfTzVX5TTDJPCPVXm0b5yT7OV5tV96y1/audItVvVRerO/jBBO3nmZ+BwwtwOIwzI2OMDo3nzkKRqcA/aSftFEWp9L8euAEcFOyURO0L3yjnJRW2up8fRA4kPy1D7g//Ou99nGUSznNGNfxNZVu0u84+/gSX8MsY1U7pccM40yeneDEjQdheBHOjMDktjv6RvQRSbwTth88xfjeGcZGTlda7OfEXfRXueEP+SiWVVbUS9vkh/mqIV3fz82PsLgwzOiOOebOjHJyUp26O+3cPVv1na4P5qv2RP/R/r3Vlq6fIj/tj/4Uc/UE2apsTapTpdP8VzN1fC/bhxYYvmiRhduHODu1q+sb/fOj6Xh26KuvC4jjeNFcv5uTvfbFT+l0/MRP+clW3V5+jeM1P43I5zHpHBDHd/gj+k3UL/LTPM4fyld5qrylO3XBJGp2Lsl3R1vi3FOlVWWjUO24NW3TvvwcmDciMo0TTtjG+TNvYG5TeXyZE1R+/gzh49437VMaWFrqzi/kZAC8kOq77F4FDIBlfcIAaAA0ABoA7/jfzwBoAFxlTGkVAP9DQxHAtzgCWIYJmze1AbDMNwZAA6AB0ABoADx3dcMRwNWGlFYB8HlTMFL4FPD8NLzNAFiGCZs3tQGwzDcGQAOgAdAAaAA0ANYaSQyAtWSyUUsKGADLhDYAGgANgAZAA6ABsNZI0ioA/lhDEcDfcQSwlnP70MgAWOY0A6AB0ABoADQAGgBrjSStAuAzGwLAdxgAazm3D40MgGVOMwAaAA2ABkADoAGw1khiAKwlk41aUsAAWCa0AdAAaAA0ABoADYC1RpJWAfCHGooA/qEjgLWc24dGBsAypxkADYAGQAOgAdAAWGskaRUAn94QAL7HAFjLuX1oZAAsc5oB0ABoADQAGgANgLVGki0IgK8Afjx91eDvgRcAn1tBjI8BD4XqLfBiD71K/EXAb2b2TwVeBdwDuAF4KfC+WuLaaN0KGADvKtnLgV+C6nMQ0Un/HPihZdQ1ABoADYAGQAOgAbDW4NsqAD51Ci4qfA/g7dPw3hUjgD8HvBB4PPAlQGPnM4Er0vjZq8lHgU8ku+X0ejAgSHwGoDH3KcAfAA8HrqklsI3WpYABcHkAfCTw7TWUNAAaAA2ABkADoAGwxnABrQLgUxoCwD9bEQC/DLwJeHNqvD4yeBS4EvjDZQQRAP4v4GUriPW76UOW35ft/9P0Eczn1BLYRutSwABoAPS3gP0tYH8L2N8CvuNMqE/mxnfF8/OjvwW8rsF1OeMtBIAKfujff13S/dusrYeBzwD/aQUAvD+wHTgG/Bnw6vRVc5kryvce4HVZ2l8EBITfWiy+M7iLAgbA5QFQnVeXgPX7JPCSdD9Cr7UjgI4AOgLoCKAjgI4A1sKLVgHwiQ1FAP9i2Qjg3YF/Ae4LfDFr/LuBaeC5ywjyEOCfEjh+A/DfgGuBH0i21wFvAN6WpX1+iijqsrKnhhUwAN5V0PsBM8AR4NLUIdVxv3GZ+xoMgAZAA6AB0ABoAKw1NLcKgI9rCAD/x7IAuJEIYK9Gus3qQ8B4ejDEEcBavag5IwPg2lqOAFPAk1JnzVMkAPz/AN3+oIeavg14HHA7cFGaK8lpuqOEtms+nB6c0vUW2Wk97LSuY0LzbApv6fKMslnM9mlbXLZRVWJSNrEvtqvI3anI/fgSsC8B+xKwLwHfcc7wJeC1R4V1WBw+fBj9NM3Pz3PVVVdpUVSlSNn5mLrj0kYB8JbD8G/d+rI4D19Zsb7L3QN4E/AzK9wDuBIAqr5nAN0DqGU9CRyT7wE8Hz0k5WkAXFvcAMAnA3/dY54A8MPArrQrYE6QpylAUMd6AKBgMOxyANR+2WmbshYUBgTGsuZjWd7KK/LQPCaVqyn29dRc2QQIquo9rHmu6oJMgeOpVB2BY0wnUzA/1kdTPhMJOmWrJsW68lE5AtLYpmXZCUSVPqa96RSpOqocTbLVaULTJUuwsK2773japjYFR2uT8g0oHobtnVMMX7TI/JlRWBhi++gcozvnq8Tje2YYYZ7TZ7saDm1frNZHmeMAtzDEAosMV/NR5quon7Zrrt8Ek5W99ndTjZyzV357q9tlJPlp5qt9XRpXmlgeQmXqDQli+25ZB1PZmmvffk5U+7ulKrduDZSH6qByJ5moth/jIDOMV7+jXMotHOQ4+6q8ZaM0sld9bzp6GZwZgTPbunpKZ/nu5qSpdNU2zeP/l1iWzl8FHJpj98QM42MzHOQYu5mpypIaapvqoXloWOmctqt+qkdosXh2qPLB3PwIw8Pd/3S0ru2zJ7v9f3FhiPGJGXZvn2G40m6+apfy1JT7JrZ1S+9qL420rn3SUPaaa1I+k9MTnP6Hvd07ndRWtTP6X+gh4+h/2q9+HffKRd/T+o2pH6u/all56hf9VzrruoOqrv8jq38YNek4Vp+Mc0faXK3rDhXVV//lxRQHcqTXXHb5eUK2sgv+UBlhr7L2pXJ7/4mNdFEvVThsos6RZ35CiXNRnBNlmx+oeZvCNj/vZc2rtIgp/qnOz5O957079i0t6WrlhZ1ajQA+qqEI4IdWfAhEt0rpKeAnAIJBPdzxw8DXLnO17ADwwPQQiDrk1wO/n26telryip4C1oMiegr4g+kp4HekqIqfAj4PXdcAeFdR1Rk/kp48OpguAesxdN2zIAzKJwNgl2u6kwHQAGgANABWJwMD4HLjdasA+O+nYLjwNTAL0/CxVV8E/cvA89J/IX+XvQfwEPD5dDns6vRevz9Or4jRf2D6l+dPeh4CkWR64EMPhtwzweGLgfefB/Zxluk9dxbizgroySTd86e42G3pvUV6L6D+w+mdDIAGQEcAHQHsnhccAcyiigbAAQFA80MfK+AIYJnzDIAGQAOgAdAAWCmgS7K+BLzakNJqBPDbGooA/i9/Cq4MEzZvagNgmW8MgAZAA6AB0ABoAKw1krQKgA/bvLmzAAAbqklEQVRtCAD/twGwlnP70MgAWOY0A6AB0ABoADQAGgBrjSQGwFoy2aglBQyAZUIbAA2ABkADoAHQAFhrJGkVAB/UUATw044A1nJuHxoZAMucZgA0ABoADYAGQANgrZGkVQB84BQMFT4FvDgN/9cAWMu5fWhkACxzmgHQAGgANAAaAA2AtUYSA2AtmWzUkgIGwDKhDYAGQAOgAdAAaACsNZK0CoDf2FAE8B8dAazl3D40MgCWOc0AaAA0ABoADYAGwFojSasAeL+GAPDzBsBazu1DIwNgmdMMgAZAA6AB0ABoAKw1khgAa8lko5YUMACWCW0ANAAaAA2ABkADYK2RpFUAvKKhCOC1jgDWcm4fGhkAy5xmADQAGgANgAZAA2CtkaRVALz3FGwvfAr47DR82QBYy7l9aGQALHOaAdAAaAA0ABoADYC1RhIDYC2ZbNSSAgbAMqETAL4fGAFW+gD67akUfStTywtpfTjNlS5stE/blbW2h01UdGdPjWUf23LbKGO8p155mUqn32r/Jeb1UtF5+WrPiaw+UXdtUt3DVvNMG/U6QcNo1kSt7wJU3R0py4m0PJOyUvZRHe3TtAhov+Zn0jblE1XNy9jfI6fKUbmSRGlUTQGt1pVXbDuYlaXtSjcE7IXtnVOM751h/0hXhxHmmWUn84wyOT3B4sIQIzvm2Td2vNo/XFUUZhljscoEZmfHGN0xx8j2ecaZYYQ5RplniMVqXfP9HK/yjv2LqV/MM8IkEyykvE4zxgzj1frJyrr7m7x1gvmb93S10nQs00LbTqU2q33xk94nMz9Je7lR7Y9upCbIZ5cDX9PVY2z8dNXeserbsJJqoWq32qyf1mOf9p9gH7dce4/uYRG/8GXore3y5VzmIyVWP4i+MDEHx0e7dduxxPZds1VdpO3O7d26dLuMFB1i8WxX/8oHJ8cY2z1bzU9PjsPkMMhlN6TyZBT9Mg4tZSk9ZKcyoy7RhtxeNtqu/qVlzaVtrGvbuSrGcaXMtazM5aQ72nDHcjgiP/ZU2bDVPA6c/Pwkm3BoLGsuGzk0psi395yldeU9C4wlY+UX3wKOvOKAzc9xsU/zvP5RZn4eU96Rh/bn55K8rDjX5PWM82e33UtLD8vadWEWWwXAezYUAfyKI4AXprec/1INgGUaGwANgAZAA6ABsDqPGgDXGk5aBcCvaggA/9UAuJZf+3W/AbDMcwZAA6AB0ABoADQA1hpJDIC1ZLJRSwoYAMuENgAaAA2ABkADoAGw1kjSKgAebCgCeMwRwFrO7UMjA2CZ0wyABkADoAHQAGgArDWStAqAFzcEgLcaAGs5tw+NDIBlTjMAGgANgAZAA6ABsNZIYgCsJZONWlLAAFgmtAHQAGgANAAaAA2AtUaSVgFw4v9v70xgLDnqO/xbbIMBx+vYGAPBBJDMJe4kHCLBigwsN1FkBAmnELc55A2nDRgDMSJgR4C5L8mAYu4j4VhxCRRQbC6BBcicgSWAE3zM2hAb25noN1N/u9zbb169Vz09XeuvpdHM61dVXfX93+v6prq7akXaVjkP4Ooe6SJGAIuC22AiBLAuaAggAogAIoAIIAJY1JOMKoAHDSSAlyCARcFtMBECWBc0BBABRAARQAQQASzqSRDAIkwkGokAAlgHGgFEABFABBABRACLepJRBfDAgUYAL2UEsCi4DSZCAOuChgAigAggAogAIoBFPcmoArjfQAJ4JQJYFNwGEyGAdUFDABFABBABRAARwKKeBAEswkSikQgggHWgEUAEEAFEABFABLCoJxlVADXQCKAYASwKboOJEMC6oCGACCACiAAigAhgUU8yugCqchoY7ZEQwKLYtpgIAayLGgKIACKACCACiAAW9SQIYBEmEo1EAAGsA50E8FRJ189KukLS5T2vD+js9+tZ2/7pjUjj8vy39/eV000f5XbTR3nxO/Ll9Y+8/9tpQ7esaKfT+T2/jjQuL29f8MiPG+3xYSK9bSLSGG+UJyn/tB4k6Xqpepel35dIWk1FOVv3fb+OQ7qq3g6UtJ+kQ9Lr36VDuizv83FulNI4yZXpb/92Xh/H6ZzGm187b/wEwcNS3shjZD6WX1+UyjQC19G/XU6UfWk6htP6OL4ik38Ezk9pfSy3y+mNO+pqPj7GbyX5uC7fr12+2+50sTm/6+4y/OPXrlMwdh7XI7i5DTeUdCtJN0/1jDTe73oe4qD4mNvWy/ZmZi7D5f8y/fhv74vjRYziIxr73Y74uEU5bofbFe9FLFxexDk+ntHW4OTXK4mr2+n9Eb/fpNdRRuTNP+5xTMfGn5f8sxXpXI7/Nvf48TEultb+h7owCnaA/HNe2uERGEP0Pgc1Cs+/T/FeFser/ozyAkjk92uXEd/7HHZ+Luiea/JjRNl5+jhn9NWzW79oQ16nqKd/952TZp1T4kPSLevq8/Dq6jP6AI26DwEcFTcHm0MAAaz7iCCAa50IAogAIoBrnwEEMJ1REcC+rgUBrOtwyT0sAQSwjicCiAAyAsgI4ProHgKY/hn0SRUBvJYI4MmSniyt3Wz4TUnHSfpeT9sPl/RaSfdN3xSPm58p6eWS/pDSHy3pS2n83bvsJx4bv0VdN03uWQQQwLrPBgKIACKACCACeI1LzgjgrG5lHxsBfL6kZ0l6kKSfSDpJ0uMl3UbS7zsMfJPIoyR9QNLPJN1a0sckfUHSzkwAv5hulkj3jdR10OTemAACWPcJQQARQAQQAUQAEcCinmRcAfQg2xBPAa/d4Ow7j31Dar79VNJpkk5PO30n7q+S0L2/AMhzJT1R0t06Anjdzp3JBUWRZBkCCOAy1K7OgwAigAggAogAIoBFPck+JIDu+/w4070lnZU1fpekcyQ9rwDIp9LTTk/qCKAfC/MjXy7nlZK+UlAWSZYggAAuAS3LggAigAggAogAIoBFPcm4AuhHz4cYAbxJ3wign/n/haTbSzo3a7zv6/NI4VPnAHmppKdL+os0aujkN5Z0RLqH0I+A+/1XSbqnpO8WASbRQgQQwIVw7ZUYAUQAEUAEEAFEAIt6knEFcPeSAvh5Sb4Vz5vnRXpnnwDWjAB6VM/3Ch4j6cdzwLkiX5P0kiLAJFqIAAK4EC4E8JoT0DEPIPMAMg/gVXMP8hQwTwHP6U/aEMC8ER7MO7JPAL2v7x7AX0s6XtKsewDfJOn+Sf5sqPM2PyTyH5JOnJeQ9xcngAAuzizPwQggI4CMADICyAggI4BFPcm4AuiHbYe4BOwHeHsfAvF9fn4K+CFJBl8m6bGSbtvzFLAfEDlD0p0l3S+b6Tzn9gBJP5T08zQl/NMkvVrSfSR9qwgwiRYigAAuhIsRQEYAWQmElUA6C9+wEkh2YmQlkI26lHEF8EeSvLJSzealao6aJYDe73n8LGo+0DeyeQA9bPh9SQ+U9NU0/5/n+PM15Zgk0v7h6V7CUn2Z9ymSDk1DyX4I5BWSvlzTAvLOJoAA1n06GAFkBJARQEYAGQFkBLCoJ9kHBbCo3SSaJgEEsC4uCCACiAAigAggAljUk4wrgB6AG2IE8A4bjQAWtZtE0ySAANbFBQFEABFABBABRACLepJxBdAzpwwhgL5tr/cewKI2k2i6BBDAutgggAggAogAIoAIYFFPggAWYSLRSAQQwDrQCCACiAAigAggAljUk4wrgH5w9qCies1OdImkuzMCWElxqtkRwLrIJAH0w0sHdkqKB52WOcD+PZlKy+vLG8Ud0CnXaT3her5FGr8Xx4y//Tvev7xnzq8oJy83ynD6/PhRVvc4LiPS+WEw189lxKWMeM/78zKjHM9b5TSxbKXT+Um2vA0uLzg5rctxKCONjxv1SH/7m+KJDK5Mb+VLlbtqnhffCyO5KD/n5kPmm4uOvC7HW1TjCmfK2p1XzX/7wUo345D040WSfEwfI/C4iNjn4zuNj+dJ/KMM18/vuZwIiz+28drpXDenM5ZLU7Wi3c7jcn1M9yvuG/xzmKRbpunCPBee0/nHYug6uXyXcX6qs/fl9fIy8jEjmN9z2c7vsv3br10319VP3Lp+Eaoox1X18VZSnf2+X3fbkX98zSLi4GPFA6ze7+P49XmpHj52fMXdlvjY/C7VxVxumvW3Tut9cfwL09//k8qMY5mJF2zw7GlrW3yvLkh/O8jeF78jXX4++L0kp8+3/HsXB+skWXvpYHbPAfnraGjkzb/DLre7PGx+rL3PWaurJ60VtG3byVmg+uoV34eNzoXxvcm+O9coKn9fWl3d+qnkxhXAswcSwHuYat9awLMCx/5GCCCAdYFCAK/qQXOQCCACiAAigAhgt3tBAOs6XHIPSwABrOOJACKA66NcjACuLxjACOD6GYURwGyY+eqTLCOAe7R9uwfTNnVELfVLnn5viEvAnod5U+tb1wuTe2kCCODS6NYyIoAIIALIJeD1y71cAs6uo/vEwAjg1o4Aev7kIQTwaASwzhMmmxsBrAsNAogAIoAIIAJ41XmUewA36lLGvQSMANZ17/t+bgSwLsYIIAKIACKACCACWNSTjCuAX0g34hZVbUYiP+V0DCOANQgnnBcBrAsOAogAIoAIIAKIABb1JOMK4K6BBHAHAlgU3fYSIYB1MUMAEUAEEAFEABHAop4EASzCRKKRCCCAdaARQAQQAUQAEUAEsKgnGVcA/22gEcCHMgJYFN32EiGAdTFDABFABBABRAARwKKeZFwB/PhAAvg3CGBRdNtLhADWxQwBRAARQAQQAUQAi3oSBLAIE4lGIoAA1oFGABFABBABRAARwKKeZFwB/LCkGxTVa3YiLzN4LCOAlRSnmh0BrIsMAogAIoAIIAKIABb1JOMK4JkDCeCjEcCi6LaXCAGsixkCiAAigAggAogAFvUkCGARJhKNRAABrAONACKACCACiAAigEU9ybgC+N6BRgAfxwhgUXTbS4QA1sVsTQBXVlZ08MH+ky0IbNv2jwUwDpiRJpaTur6k/VOaWXz9vsuJdFGky4i1SPPlqbrpnMb7fKx8i9c+7qHpjSjHeS5P5Ud50ZY9qT5+3+nz4+Vro3q/33c+p40t6uPX/rtbL+93+m4+v87Tut5Rxy5m39dja4vN+VyXSH9+9l5fjNzG2C7WugG7jG6M+vLmbXUZPm78uBxvEc9oU7QrLy/2xe8LUt788+D0EbPY7/RRh8jjrMHA7Yk85hR/e3/kOyJjnX+28rJ971W0J9rk1z6my/KWvx8cg0l8joJzfpxIk8c1YhKfsajr+uvV1ZO6H4Ki19u2nVyUblaiZY/r8srOIbOrt7p6YlXdNyPzuAL47oEE8ElGsV1S/sXfDDyUOTIBBLAOOAI4g1/ZyRsBRAARwPWvEALYPZWUnUMQwB4C6coUAljXve/7uRHAuhgjgAhgIsAIICOAjAB2TweMAF6TyLgjgG8faATwqYwA1nnCZHMjgHWhQQARQASQS8DpM4AAIoAbdyjjCuCbZ9xCskin5xH6ZyKAiyBrKC0CWBcsBBABRAARQARwxnmAEcCtHAFEAOu6930/NwJYF2MEEAFEABFABBABLOpJxh0BfMNAI4DPYQSwKLrtJUIA62KGACKACCACiAAigEU9ybgCeOpAAvgPCGBRdNtLhADWxQwBRAARQAQQAUQAi3oSBLAIE4lGIoAA1oFGABFABBABRAARwKKeZFwBfM1AI4AvZASwKLrtJUIA62KGACKACCACiAAigEU9ybgC6Mn4Dyyq1+xEl0pam1CbiaArSU4xOwI4OyqeAv/JaXmDb0o6TtL3OskRQAQQAUQAEUAEsKh/RwCLMJFoJAIIYD/o50t6lqQHSfqJJK+j9HhJt5Hk9aFi22cEcNeuXdqxY8dgH7uyWfyHXgnkO5LuktrQ8lJwP5B0p8aXgvu2pHv0fJ5aWwru7BSLtpeC++xn77XU93taS8H9SNJRxecoloLzGMYQI4BrywgyAlj8yWsnIQLYH6ufSjpN0unp7f0k/UrSTknv3xcFcOfOnTrtNDd5mG1rBPCM5OluQ8sC+CFJf9e4AL5H0toKAp2tNQH0clqORdsCePzxK0t9v6clgJ+W9ODiExQC+JKBBPBVCGDxp66thAjg3vHyqN5Fku4t6azs7V2SzpH0PARw/occAbToWEL3z2Bdkf3t/X7fo6C5FDkNArjOzWy8uoZ/vOUjxrEvfl+Q0kS++G3G3vLygnfk8fsux2ktepHHg/0I4Pxv++wUw04EjQAWxiKtBYwAFvK61iZDAPcO/c0l/ULS7SWdm719pqQ9nWGNtS/a7t27dfDB/rPd7YQTTtApp5wyWAO2b39tQVmLXAL+oxnluWOPH4fo0Smdb14O4YoOPUQgL8ppnD9EIt6LSyeO6x+nnVGO88RPCF78tkCE1M0TQNfR+XIx9N8fl3Rs2t93CSfq3BXKvA3mFXXsorPYHJrtDPmJ9LkY5QIbWdzG2C6RdFDil8fIA+VP7IlZ3wigOZhV5I94hgAGg/zzEvuizRdmopfn7wqg8wW3yOOs3u96uG3+7c2c3pdi4f2R7/Ds85J/trpLwZlLbK67X/uYZuYtfz84er+P71NN3t78OJEmxxsxiToG5/XXxx138VLf7+3bX90Tw/JdKysvLk/cSbn3OeRzku5fXN7Kiu/kmdbmewCPPPJIV2ozL6kmATxe0vUqAVwm6Z83u76VdST7sgQQwL3JLTIC+CeSfrksfPJBAAIQgMC1koAHGv5rk1ru/2Z+JukmA5X/G0m3yv4zGqhYitlqAghgfwT67gH8tST/S5XfA2h+N0vDBlsdS44PAQhAAALTJ+Chbt9TvrqJVbUEXneg8v+A/A1EcmLFIID9AfF9fn4K+CGSLIMvk/RYSbftPAU8sXBSHQhAAAIQgAAEIDCfAAI4m9HLJT0t3Zj0jRnzAM4nTAoIQAACEIAABCAwMQII4MQCQnUgAAEIQAACEIDAZhNAAOsIl6wWUneE4XJ7Ns+XpkvYjrvvP/lXSY9Jh7izpDdK+rM0Dc47JLl9W709Ko2+eoZnPzbpRyP/L6tUSb23Ok7z2uD2xGPLERtPQ5SvPLPVbTByPxLq2yL+ND3K+mVJL+g8CDX1eJS0oYV4+LYUT05/I0m+R8urFb1IkmdDj23qsShpQwux6J4jPybpEZLuJ+mLC5xfp/Ad3+rzPccfkQACuDzs0tVClj/CsDktgMdIum9PsRarH6ZJz16RVjz5jKTXSXr9sNVYuDTP++A5S24g6Z0dASyp9xTitFEbDMSdnGPzpRl0ptAGV82Li344zYfpeLxF0h0k3S3Vu4V4zGtDK/Hwkhj/7Wmo0lxCz0kyftP0z10LsZjXhlZikX9tLeV/n+ar8ffeAthCLBY+MZOhfQII4PIxLF0tZPkjDJtzIwF8gqTXpCeaY3TNHcqzF1p7adj6dks7Op1M8xHAknpPKU59bYhOLh8t6LZ9Sm3I6+ZR2W8lQbeItBYPt6Xbhhbj4cneniHpVEk3lnR+g7Hoa0NrsfDULv8u6S/TXLLxnW7xe7G5Z3NKnwQBBHC5MCwyV+ByRxg+lwXQTzd7dlv/fE3SiZL+My1754mvvfZxbL4E6ZOZJyyNmWuHr1V5iX3y5LXrNqr3dRZY1aW8Jsun3EgAPdeW5fbnkt6aRjt9pCl/1nz59+mSbp2QtBYPV7vbhpCOFuLhddE8LZW/o/7HzTP2xuzHrcRioza0FAvX1atFfVDSu1I8QgBbicXyZzZyNkkAAVwubIusFrLcEYbP5Ut1Xi5gtyRfJvJSHfdMIyBvkHTDtOhpHPl26R40T1vvOau2euuTJ18S3qjeFsDSVV3GaN8sAfzrJORXpktH7tS9hMLbJE31s+bOzfc6/a0kL9HgrbV49LXB7WgtHoekET9PSv+RRmPR14aWYvHMdN/fjsTfQh4C2Nr3YoxzIceYAAEEcLkgTHlUprRFniTUl+0ellZYZwSwlNzy6WYJYLdE3xz/gHQpaYqftYdKem+Sjk9mlW9ppGNWG/qiO/V4uM4+l3udub9K92i2FItg3m1DK7HwCLivlvgfav+D7S0XwBZjsfxZjpzNEEAAlw9V6Wohyx9hc3OGAD48jQj+U+cewOemybB9o/YUtj558g3X8+o9pTgtIoAeSbhPAj+lNvip8dMlPVLS5zsfjFbisVEbZknHVOMR9fUCyP6H7nGSPpqeEG7pu+F2dNvQSix8j59H672Ac/Sph6V4fCCN7vuKi6+8xD3W3fPrlL7jUzjfU4cRCCCAy0NubbUQd9h+Is03iB+RLgH7ZuU7pZPWuekpYD8laen7VLqpfKufAvZlXN8bZ3nyk8leRsmXSj31hS//zqv3FOK0URvumvifkz0N/C9p9Zk3pY/nFNrgqnh1HD8l7lHjr/Z8dfy049TjMa8NfqLZ58Wpx8MPaZ2ZngQ+PD2hfWy6J/a89OTp1GMxrw2txMLLrnmmgnzz5XhP/+TbI65o4HuxfE9IzmYJIIB1oWtptZBPSLpXkiZfKvpKmhfQ/3l6u6OkN6d5AD2S4Ck+XlmHZ5Dc/u/6Pdm6mTFPnu/TchtK6r3VcdqoDb7E65Ea3+vnjsIPgTgOnocx37a6Da6LRy8ul3RZqljEwg8PhRBOPR7z2uBLwy3Ew3N4/nkSPY88fT3JuZ/Kjm3qsZjXhlZi0Xeii/t5Yx7AqcdikJM1hbRFAAFsK17UFgIQgAAEIAABCFQTQACrEVIABCAAAQhAAAIQaIsAAthWvKgtBCAAAQhAAAIQqCaAAFYjpAAIQAACEIAABCDQFgEEsK14UVsIQAACEIAABCBQTQABrEZIARCAAAQgAAEIQKAtAghgW/GithCAAAQgAAEIQKCaAAJYjZACIAABCEAAAhCAQFsEEMC24kVtIQABCEAAAhCAQDUBBLAaIQVAAAIQgAAEIACBtggggG3Fi9pCAAIQgAAEIACBagIIYDVCCoAABCAAAQhAAAJtEUAA24oXtYUABCAAAQhAAALVBBDAaoQUAAEIQAACEIAABNoigAC2FS9qCwEIQAACEIAABKoJIIDVCCkAAhCAAAQgAAEItEUAAWwrXtQWAhCAAAQgAAEIVBNAAKsRUgAEIAABCEAAAhBoiwAC2Fa8qC0EIAABCEAAAhCoJoAAViOkAAhAAAIQgAAEINAWAQSwrXhRWwhAAAIQgAAEIFBNAAGsRkgBEIAABCAAAQhAoC0CCGBb8aK2EIAABCAAAQhAoJoAAliNkAIgAAEIQAACEIBAWwQQwLbiRW0hAAEIQAACEIBANQEEsBohBUAAAhCAAAQgAIG2CCCAbcWL2kIAAhCAAAQgAIFqAghgNUIKgAAEIAABCEAAAm0RQADbihe1hQAEIAABCEAAAtUEEMBqhBQAAQhAAAIQgAAE2iKAALYVL2oLAQhAAAIQgAAEqgkggNUIKQACEIAABCAAAQi0RQABbCte1BYCEIAABCAAAQhUE0AAqxFSAAQgAAEIQAACEGiLwP8DD7e7iYTBJ0YAAAAASUVORK5CYII=\">"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x11519d5c0>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fig = plt.figure()\n",
    "ax = plt.pcolormesh(rainfall_eff[0].T)\n",
    "plt.colorbar()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## DECADAL CHANGES IN EFFECTIVE RAINFALL\n",
    "\n",
    "#time limits for different time period (periods given below)\n",
    "p1_beg = [1951,1980]\n",
    "p1_end = [1979,1993]\n",
    "p2_beg = [1980,1994]\n",
    "p2_end = [2007,2007]\n",
    "\n",
    "rainfall_p1 = [rainfall_eff[1], rainfall_eff[3]]\n",
    "rainfall_p2 = [rainfall_eff[2], rainfall_eff[4]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "rain_diff = []\n",
    "\n",
    "for rain1, rain2 in zip(rainfall_p1, rainfall_p2):\n",
    "    rain_diff.append(rain2-rain1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "8007_5179\n",
      "9407_8093\n"
     ]
    }
   ],
   "source": [
    "## SAVE AS NETCDF\n",
    "## SAVE FREQUENCY DIFFERENCES AND P-VALUES AS NETCDF\n",
    "outfilename = \"RDA_effective_rainfall_diff.nc\"\n",
    "ferretpath = \"/Users/Siwen/Desktop/Ferret/bin/\"\n",
    "\n",
    "#netCDF output file formatting formatting\n",
    "fullname = ferretpath + outfilename\n",
    "\n",
    "##the following is only to make sure that we don't exceed the bounds of the number of rows of actual data.\n",
    "if os.path.isfile(fullname):\n",
    "    os.remove(fullname)\n",
    "    \n",
    "outfile = nc.Dataset(fullname, 'w')\n",
    "\n",
    "#use same years as chosen above:\n",
    "#p1_beg = [1951,1980]\n",
    "#p1_end = [1979,1993]\n",
    "#p2_beg = [1980,1994]\n",
    "#p2_end = [2007,2007]\n",
    "\n",
    "#create dimensions in NetCDF file\n",
    "days_dim = outfile.createDimension('time', None)\n",
    "lat_dim = outfile.createDimension('lat', 22)\n",
    "dim_tuple = (\"time\",\"lat\")\n",
    "\n",
    "#run through loop and save all differences of interest, as well as p-values\n",
    "for p1_b, p1_e, p2_b, p2_e, r_diff in zip(p1_beg, p1_end, p2_beg, p2_end, rain_diff):\n",
    "    \n",
    "    #convoluted line of code below guarantees that 2007 will become '07' and not '7'\n",
    "    yrs_p1 = str(p1_b-1900) + str(\"{0:0=2d}\".format((p1_e-1900)%100))\n",
    "    yrs_p2 = str(p2_b-1900) + str(\"{0:0=2d}\".format((p2_e-1900)%100))\n",
    "    years = yrs_p2 + '_' + yrs_p1\n",
    "    print(years)\n",
    "\n",
    "    #create variables inside netCDF file\n",
    "    rain_diff_out = outfile.createVariable(\"rain_effective_diff_\" + years, 'f8', dim_tuple)\n",
    "    \n",
    "    #fill with actual values\n",
    "    rain_diff_out[:] = r_diff\n",
    "    \n",
    "outfile.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#nice function lifted from RDA_intensity.ipynb\n",
    "\n",
    "##collect function designed to aggregate all observations from 3D storage matrix - takes every ensemble (3rd dimension)\n",
    "#then selections range of days and latitude range (I've mostly been using 15 day and 2 degree-latitude smoothing)\n",
    "\n",
    "#RETURNS every non-NaN observation, as a single numpy array\n",
    "#inputs: P is the full 3D matrix, dd is the day of interest(0-364), yy is the latitude bin of interest (likely 0-21)\n",
    "def collect(P,dd,latbin,dayrange,latrange):\n",
    "        \n",
    "    ll=int((dayrange-1)/2)\n",
    "    yy=int((latrange-1)/2)\n",
    "    \n",
    "    days = P.shape[0]\n",
    "    lats = P.shape[1]\n",
    "    \n",
    "    Psample = P.take(range(dd-ll,dd+ll+1), mode='wrap', axis=0)\n",
    "    P_collect = Psample[:,max(0,latbin-yy):min(days,latbin+yy+1),:].flatten()\n",
    "    P_return = P_collect[~np.isnan(P_collect)]\n",
    "    #print(len(P_collect))\n",
    "    #print(len(P_return))\n",
    "    \n",
    "    return P_collect[~np.isnan(P_collect)]"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}