pltWrfStageIVgrb.ipynb

{
 "metadata": {
  "name": "pltWrfStageIVgrb"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Plots WRF Fields and GRB data"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "from mpl_toolkits.basemap import Basemap\n",
      "import matplotlib.pyplot as plt\n",
      "#import pygrib\n",
      "import matplotlib as mt\n",
      "import ttuatmo.radar.ctables as ct\n",
      "import glob\n",
      "import netCDF4 as nc"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 14
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab osx"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n",
        "Welcome to pylab, a matplotlib-based Python environment [backend: MacOSX].\n",
        "For more information, type 'help(pylab)'.\n"
       ]
      }
     ],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "grbs = pygrib.open(\"/Users/aereinha/Downloads/stageIV/ST4.2010051912.24h\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 29
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#pcp = 0.\n",
      "#del(pcp)\n",
      "p=0\n",
      "for i in glob.glob(\"/Users/aereinha/Downloads/stageIV/ST4.20100518*.01h\"):\n",
      "    grbs = pygrib.open(i)\n",
      "    if (p == 0):\n",
      "        grb = grbs.select(name='Total Precipitation')[0]\n",
      "        pcp = grb.values\n",
      "    else:\n",
      "        grb = grbs.select(name='Total Precipitation')[0]\n",
      "        pcp = pcp + grb.values\n",
      "    p=p+1\n",
      "grbs = pygrib.open(\"/Users/aereinha/Downloads/stageIV/ST4.2010051900.01h\")\n",
      "grb = grbs.select(name='Total Precipitation')[0]\n",
      "pcp = pcp + grb.values\n",
      "grbs = pygrib.open(\"/Users/aereinha/Downloads/stageIV/ST4.2010051901.01h\")\n",
      "grb = grbs.select(name='Total Precipitation')[0]\n",
      "pcp = pcp + grb.values\n",
      "grbs = pygrib.open(\"/Users/aereinha/Downloads/stageIV/ST4.2010051902.01h\")\n",
      "grb = grbs.select(name='Total Precipitation')[0]\n",
      "pcp = pcp + grb.values"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "pcp = grbs.select(name='Total Precipitation')[0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 30
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "lats,lons = pcp.latlons()\n",
      "print(lats.max())"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "57.8475794356\n"
       ]
      }
     ],
     "prompt_number": 27
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "a = nc.Dataset(\"/Volumes/External/aereinha/lajunta1km/200906112136/Posterior_Diag.nc\",\"r\")\n",
      "b=nc.Dataset(\"/Volumes/External/aereinha/lajunta1km/200906120130/precip/wrf_precip_d01_2009-06-12_01:32:00_1\",\"r\")\n",
      "xlat = a.variables['XLAT_d01'][:]\n",
      "xlong = a.variables['XLONG_d01'][:]\n",
      "Tdfi = a.variables['TH2_d01'][0,0,:,:]\n",
      "shape = xlat.shape\n",
      "xdim = len(a.variables['XLAT_d01'][0,:])\n",
      "ydim = len(a.variables['XLONG_d01'][:,0])\n",
      "print(np.min(xlong[:,:]))\n",
      "print(xlong[:,:].min())\n",
      "print(b.dimensions)\n",
      "ref = b.variables['REFL_10CM'][0,0,:,:] \n",
      "\n",
      "\n",
      "#b=nc.Dataset(\"/Volumes/External/aereinha/DFItests/wrfoutnodfi/wrfout_d01_2012-06-17_00:00:00\",\"r\")\n",
      "pcpwrfnodfi = b.variables['RAINNC'][0,:,:]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "-105.803\n",
        "-105.803\n",
        "OrderedDict([(u'Time', <netCDF4.Dimension object at 0x13fd16ef0>), (u'south_north', <netCDF4.Dimension object at 0x10c89b638>), (u'west_east', <netCDF4.Dimension object at 0x10c89b5f0>), (u'bottom_top', <netCDF4.Dimension object at 0x10c89bd40>)])\n"
       ]
      }
     ],
     "prompt_number": 20
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "figure()\n",
      "m=Basemap(resolution='i',projection='lcc',llcrnrlon=xlong[70,100],llcrnrlat=xlat[70,100],urcrnrlon=xlong[180,260],urcrnrlat=xlat[180,260],lat_1=a.variables['TRUELAT1'][:],lat_2=a.variables['TRUELAT2'][:],lat_0=a.variables['CEN_LAT'][:],lon_0=a.variables['STAND_LON'][:])\n",
      "m.drawcoastlines()\n",
      "m.drawstates()\n",
      "m.drawcountries()\n",
      "shapefilelocation=\"/Volumes/External/aereinha/shapefiles\"\n",
      "m.readshapefile(shapefilelocation+'/counties/c_03ap12', name='counties')\n",
      "x,y = m(xlong[70:180,100:260],xlat[70:180,100:260])\n",
      "levs = np.arange(-15.,60.,5.)\n",
      "cs = m.contourf(x,y,ref[70:180,100:260],cmap=ct.Carbone42,levels=levs)\n",
      "#cs=m.contourf(x,y,pcpwrfnodfi,colors=pcp_cmap,levels=[0.01,0.1,0.5,1.,2.,3.,4.,5.,7.5,10.,15.,20.,30.,40.,50.])\n",
      "plt.colorbar()\n",
      "plt.show()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 21
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#m=Basemap(resolution='l',projection='stere',llcrnrlon=-125,llcrnrlat=25,urcrnrlon=-67,urcrnrlat=49,lat_ts=50,lat_0=50,lon_0=-107.)\n",
      "#m=Basemap(resolution='i',projection='lcc',llcrnrlon=xlong[0,0,0],llcrnrlat=xlat[0,0,0],urcrnrlon=xlong[0,-1,-1],urcrnrlat=xlat[0,-1,-1],lat_1=a.TRUELAT1,lat_2=a.TRUELAT2,lat_0=a.CEN_LAT,lon_0=a.STAND_LON)\n",
      "m.drawcoastlines()\n",
      "m.drawstates()\n",
      "m.drawcountries()\n",
      "#m.drawmeridians()\n",
      "#m.drawparallels()\n",
      "x,y = m(xlong[0,:,:],xlat[0,:,:])\n",
      "#x,y = m(lons,lats)\n",
      "cs = m.contourf(x,y,pcp,colors=pcp_cmap,levels=[0.01,0.1,0.5,1.,2.,3.,4.,5.,7.5,10.,15.,20.,30.,40.,50.])\n",
      "plt.colorbar()\n",
      "plt.show()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 50
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "pcptotal = 0.\n",
      "for i in glob.glob(\"/Users/aereinha/Downloads/stageIV/ST4.20120617*.01h\"):\n",
      "    grbs = pygrib.open(i)\n",
      "    pcp = grbs.select(name='Total Precipitation')[0]\n",
      "    pcptotal = pcp.values + pcptotal\n",
      "    grbs.close"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "red = [255./255.,214./255.,181./255.,142./255.,127./255.,99./255.,0./255.,0./255.,0./255.,99./255.,150./255.,198./255.,255./255.,255./255.,255./255.,255./255.,255./255.]\n",
      "blue = [255./255.,255./255.,255./255.,255./255.,255./255.,247./255.,255./255.,150./255.,51./255.,0./255.,0./255.,51./255.,0./255.,0./255.,0./255.,0./255.,0./255.]\n",
      "green = [255./255.,226./255.,201./255.,178./255.,150./255.,112./255.,99./255.,150./255.,198./255.,255./255.,255./255.,255./255.,255./255.,198./255.,160./255.,124./255.,25./255.]\n",
      "#red = red/255.\n",
      "#blue= blue/255.\n",
      "#green = green/255.\n",
      "cmaps = np.dstack((red,green,blue))\n",
      "\n",
      "\n",
      "\n",
      "pcp_cmap = [ #(1.,1.,1.),\\\n",
      "             (0.83921569,.88627451,1.),\\\n",
      "             (0.70980392,.78823529,1.),\\\n",
      "             (0.55686275,.69803922,1.),\\\n",
      "             (0.49803922,.58823529,1.),\\\n",
      "             #(0.38823529,.43921569,.96862745),\\\n",
      "             (0.,0.38823529,1.),\\\n",
      "             (0.,0.58823529,0.58823529),\\\n",
      "             (0.,0.77647059,0.2,),\\\n",
      "             ( 0.38823529,1.,0.),\\\n",
      "             ( 0.58823529,1.,0.),\\\n",
      "             ( 0.77647059,1.,0.2),\\\n",
      "             ( 1.,1.,0.),\\\n",
      "             ( 1.,0.77647059,0.),\\\n",
      "             ( 1.,0.62745098,0.),\\\n",
      "             ( 1.,0.48627451,0.),\\\n",
      "             ( 1.,0.09803922,0.)]\n",
      "print(pcp_cmap)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[(0.83921569, 0.88627451, 1.0), (0.70980392, 0.78823529, 1.0), (0.55686275, 0.69803922, 1.0), (0.49803922, 0.58823529, 1.0), (0.0, 0.38823529, 1.0), (0.0, 0.58823529, 0.58823529), (0.0, 0.77647059, 0.2), (0.38823529, 1.0, 0.0), (0.58823529, 1.0, 0.0), (0.77647059, 1.0, 0.2), (1.0, 1.0, 0.0), (1.0, 0.77647059, 0.0), (1.0, 0.62745098, 0.0), (1.0, 0.48627451, 0.0), (1.0, 0.09803922, 0.0)]\n"
       ]
      }
     ],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cmapp = {'red' : (red),'green': (green),'blue': (blue) }\n",
      "a = mt.colors.ListedColormap(cmaps)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#m=Basemap(resolution='l',projection='stere',llcrnrlon=lons.min(),llcrnrlat=lats.min(),urcrnrlon=lons.max(),urcrnrlat=lats.max(),lat_ts=50,lat_0=50,lon_0=-107.)\n",
      "#m=Basemap(resolution='l',projection='lcc',llcrnrlon=-117.,llcrnrlat=24.,urcrnrlon=-75.,urcrnrlat=49.0,lat_1=32.,lat_2=46.,lon_1=-101,lat_0=37.41502,lon_0=-96.46472)\n",
      "#For WRF LCC\n",
      "fig=plt.figure(figsize=(15,8.5))\n",
      "plt.subplot(1,3,1)\n",
      "plt.title(\"No-DFI\")\n",
      "m=Basemap(resolution='i',projection='lcc',llcrnrlon=xlong[0,0,0],llcrnrlat=xlat[0,0,0],urcrnrlon=xlong[0,-1,-1],urcrnrlat=xlat[0,-1,-1],lat_1=a.TRUELAT1,lat_2=a.TRUELAT2,lat_0=a.CEN_LAT,lon_0=a.STAND_LON)\n",
      "m.drawcoastlines()\n",
      "m.drawstates()\n",
      "m.drawcountries()\n",
      "x,y = m(xlong[0,:,:],xlat[0,:,:])\n",
      "cs2 = m.contourf(x,y,pcpwrfnodfi,colors=pcp_cmap,levels=[0.01,0.1,0.5,1.,2.,3.,4.,5.,7.5,10.,15.,20.,30.,40.,50.,np.max(pcpwrfdfi)])\n",
      "plt.colorbar(orientation='horizontal')\n",
      "\n",
      "plt.subplot(1,3,2)\n",
      "plt.title(\"Stage IV\")\n",
      "m=Basemap(resolution='i',projection='lcc',llcrnrlon=xlong[0,0,0],llcrnrlat=xlat[0,0,0],urcrnrlon=xlong[0,-1,-1],urcrnrlat=xlat[0,-1,-1],lat_1=a.TRUELAT1,lat_2=a.TRUELAT2,lat_0=a.CEN_LAT,lon_0=a.STAND_LON)\n",
      "\n",
      "m.drawcoastlines()\n",
      "m.drawstates()\n",
      "m.drawcountries()\n",
      "#m.drawmeridians()\n",
      "#m.drawparallels()\n",
      "x,y = m(lons,lats)\n",
      "\n",
      "cs = m.contourf(x,y,pcptotal,colors=pcp_cmap,levels=[0.01,0.1,0.5,1.,2.,3.,4.,5.,7.5,10.,15.,20.,30.,40.,50.,np.max(pcptotal)])\n",
      "plt.colorbar(orientation='horizontal')\n",
      "plt.subplot(1,3,3)\n",
      "plt.title(\"DFI\")\n",
      "m=Basemap(resolution='i',projection='lcc',llcrnrlon=xlong[0,0,0],llcrnrlat=xlat[0,0,0],urcrnrlon=xlong[0,-1,-1],urcrnrlat=xlat[0,-1,-1],lat_1=a.TRUELAT1,lat_2=a.TRUELAT2,lat_0=a.CEN_LAT,lon_0=a.STAND_LON)\n",
      "m.drawcoastlines()\n",
      "m.drawstates()\n",
      "m.drawcountries()\n",
      "x,y = m(xlong[0,:,:],xlat[0,:,:])\n",
      "cs1 = m.contourf(x,y,pcpwrfdfi,colors=pcp_cmap,levels=[0.01,0.1,0.5,1.,2.,3.,4.,5.,7.5,10.,15.,20.,30.,40.,50.,np.max(pcpwrfdfi)])\n",
      "plt.colorbar(orientation='horizontal')\n",
      "plt.savefig(\"runsVSStageIV1domainfeedback.png\",dpi=300)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "ename": "NameError",
       "evalue": "name 'Basemap' is not defined",
       "output_type": "pyerr",
       "traceback": [
        "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
        "\u001b[0;32m<ipython-input-1-9a93dd616971>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mplt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msubplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      6\u001b[0m \u001b[0mplt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtitle\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"No-DFI\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 7\u001b[0;31m \u001b[0mm\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mBasemap\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mresolution\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'i'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mprojection\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'lcc'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mllcrnrlon\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mxlong\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mllcrnrlat\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mxlat\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0murcrnrlon\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mxlong\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0murcrnrlat\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mxlat\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlat_1\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mTRUELAT1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlat_2\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mTRUELAT2\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlat_0\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mCEN_LAT\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlon_0\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mSTAND_LON\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      8\u001b[0m \u001b[0mm\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdrawcoastlines\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      9\u001b[0m \u001b[0mm\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdrawstates\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
        "\u001b[0;31mNameError\u001b[0m: name 'Basemap' is not defined"
       ]
      },
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAR4AAAICCAYAAAAZNaKrAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFoRJREFUeJzt3HFM1Pf9x/HXOfjDOqdYWyp3pFTuCqcUWAK1rml3dnE4\ntGy1XcayPzpHCCHrnFuWuGzJCv6hsvUfI3+UdtattiUs3RY6S28JxtuyKmWpDSbVOdqIO64pK1OL\nm0nB8/P7w18vUuQ44XgD5/ORNPG4j9/v+xP7fea+d6ce55wTABhaNNcDALj1EB4A5ggPAHOEB4A5\nwgPAHOEBYI7wQJJUUFCg3NxcXb58OfGzX//619qwYcNNH2tgYECLFi3S0qVLtXTpUt1111169NFH\n1d3dPeGct912W2LdF77wBX344YeJ33/16tUZ7wvzE+FBwtWrV7Vv3760He/jjz/WpUuXdPLkSW3c\nuFGPPfaYfvvb3yae93g8Onz4sC5duqRLly5pZGREd911V9rOj/mL8EDStQj85Cc/0TPPPKOPP/54\nwvPHjh1TZWWlli9frvvvv1/Hjx9P+dh33nmntm/frqamJu3cuTOdY2OBIjxIqKioUCgU0jPPPDPu\n5xcuXNDmzZu1Y8cOnT9/Xj/+8Y+1efNmnT9//qaO/9hjj+nf//63zpw5k/gZX5y/NREeJHg8Hu3a\ntUv79+/X8PBw4uevv/667r33Xn3nO9/RokWLVFtbq+LiYv3pT3+6qePn5eVJUiJYzjl94xvfUE5O\njnJycrR169b0bQbzWtZcD4D5Ze3atdqyZYv27t2rYDAoSfrggw909913j1t3991364MPPpAkff7z\nn5fH45HH49GpU6cmPXYsFpMkrVixQtK10HV2duqRRx6Zja1gHuMVDyZobm7W888/nwhFXl6ezp07\nN27NuXPn5PV6JUn//e9/E28O+3y+SY/7xz/+Ubm5uSoqKpq94bEgEB5MUFhYqG9961vat2+fPB6P\nqqur9c9//lPt7e26cuWKOjo69I9//ENbtmxJepxP378ZGhpSa2urdu3apT179lhsAfMct1q4oV/8\n4hc6dOiQpGu3RocPH9YPf/hDNTY2KhAI6PDhw4lbpsksX75czjktWbJElZWVevXVV/XVr341pfN7\nPJ4Z7wHzl4d/jweAtaS3Wt/73veUm5ur++67b9I127dvVyAQUFlZmd555520Dwgg8yQNz7Zt2xQO\nhyd9vqurS++99576+/v13HPPqbGxMe0DAsg8ScPz0EMPKScnZ9LnX3vtNT355JOSpHXr1unixYsa\nGhpK74QAMs6M3lyOxWLKz89PPPb5fBocHFRubu6EtbxZCGSm6bxNPOOP0z970mSBcc4t6P+efvrp\nOZ+BPWTOPjJhD9M1o/B4vV5Fo9HE48HBwcSXygBgMjMKT01NjV588UVJUk9Pj5YvX37D2ywAuF7S\n93i+/e1v6y9/+YuGh4eVn5+v5uZmjY2NSZIaGhpUXV2trq4u+f1+LVmyRAcPHjQZeq6EQqG5HmHG\nMmEPUmbsIxP2MF1mXyD0eDwzuicEMP9M97rm72oBMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEB\nYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnC\nA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBz\nhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A\n5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8\nAMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBH\neACYIzwAzBEeAOYIDwBzhAeAuSnDEw6HVVxcrEAgoJaWlgnPDw8Pa9OmTSovL1dJSYl+85vfzMac\nADKIxznnJnsyHo+rqKhI3d3d8nq9qqysVHt7u4LBYGJNU1OTPvnkE+3Zs0fDw8MqKirS0NCQsrKy\nxp/I41GSUwFYgKZ7XSd9xdPb2yu/36+CggJlZ2ertrZWnZ2d49asWrVKIyMjkqSRkRHdfvvtE6ID\nANdLWohYLKb8/PzEY5/Pp7feemvcmvr6ej3yyCPKy8vTpUuX9Lvf/W7S4zU1NSV+HQqFFAqFpjc1\ngDkRiUQUiURmfJyk4fF4PFMeYPfu3SovL1ckEtH777+vjRs3qq+vT0uXLp2w9vrwAFh4PvuCobm5\neVrHSXqr5fV6FY1GE4+j0ah8Pt+4NceOHdM3v/lNSVJhYaHuuecenTlzZlrDALg1JA1PRUWF+vv7\nNTAwoNHRUXV0dKimpmbcmuLiYnV3d0uShoaGdObMGa1evXr2Jgaw4CW91crKylJra6uqqqoUj8dV\nV1enYDCotrY2SVJDQ4N+9rOfadu2bSorK9PVq1f1y1/+UitWrDAZHsDClPTj9LSeiI/TgYwzKx+n\nA8BsIDwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFg\njvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcID\nwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOE\nB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDm\nCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwA\nzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgLkpwxMOh1VcXKxAIKCW\nlpYbrolEIvriF7+okpIShUKhdM8IIMN4nHNusifj8biKiorU3d0tr9eryspKtbe3KxgMJtZcvHhR\nDz74oP785z/L5/NpeHhYK1eunHgij0dJTgVgAZrudZ30FU9vb6/8fr8KCgqUnZ2t2tpadXZ2jlvz\nyiuv6PHHH5fP55OkG0YHAK6XNDyxWEz5+fmJxz6fT7FYbNya/v5+nT9/Xhs2bFBFRYUOHTo0O5MC\nyBhZyZ70eDxTHmBsbEwnTpzQkSNHdPnyZa1fv14PPPCAAoHAhLVNTU2JX4dCId4PAhaYSCSiSCQy\n4+MkDY/X61U0Gk08jkajiVuqT+Xn52vlypVavHixFi9erIcfflh9fX1ThgfAwvPZFwzNzc3TOk7S\nW62Kigr19/drYGBAo6Oj6ujoUE1Nzbg1X//61/W3v/1N8Xhcly9f1ltvvaU1a9ZMaxgAt4akr3iy\nsrLU2tqqqqoqxeNx1dXVKRgMqq2tTZLU0NCg4uJibdq0SaWlpVq0aJHq6+sJD4Ckkn6cntYT8XE6\nkHFm5eN0AJgNhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4\nAJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO\n8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPA\nHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QH\ngDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYI\nDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJibMjzh\ncFjFxcUKBAJqaWmZdN3f//53ZWVl6Q9/+ENaBwSQeZKGJx6P66mnnlI4HNapU6fU3t6u06dP33Dd\nzp07tWnTJjnnZm1YAJkhaXh6e3vl9/tVUFCg7Oxs1dbWqrOzc8K6/fv364knntAdd9wxa4MCyBxJ\nwxOLxZSfn5947PP5FIvFJqzp7OxUY2OjJMnj8czCmAAySVayJ1OJyI4dO7R37155PB4555LeajU1\nNSV+HQqFFAqFUh4UwNyLRCKKRCIzPo7HJSlFT0+PmpqaFA6HJUl79uzRokWLtHPnzsSa1atXJ2Iz\nPDys2267Tc8//7xqamrGn+j/wwQgc0z3uk4anitXrqioqEhHjhxRXl6e7r//frW3tysYDN5w/bZt\n2/Too49q69ataRsQwPw13es66a1WVlaWWltbVVVVpXg8rrq6OgWDQbW1tUmSGhoapjctgFta0lc8\naT0Rr3iAjDPd65pvLgMwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A\n5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8\nAMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBH\neACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFg\njvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcID\nwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOE\nB4C5KcMTDodVXFysQCCglpaWCc+//PLLKisrU2lpqR588EGdPHlyVgYFkDk8zjk32ZPxeFxFRUXq\n7u6W1+tVZWWl2tvbFQwGE2uOHz+uNWvWaNmyZQqHw2pqalJPT8/EE3k8SnIqAAvQdK/rpK94ent7\n5ff7VVBQoOzsbNXW1qqzs3PcmvXr12vZsmWSpHXr1mlwcPCmhwBwa0kanlgspvz8/MRjn8+nWCw2\n6foDBw6ouro6fdMByEhZyZ70eDwpH+jo0aN64YUX9Oabb066pqmpKfHrUCikUCiU8vEBzL1IJKJI\nJDLj4yQNj9frVTQaTTyORqPy+XwT1p08eVL19fUKh8PKycmZ9HjXhwfAwvPZFwzNzc3TOk7SW62K\nigr19/drYGBAo6Oj6ujoUE1Nzbg1//rXv7R161a99NJL8vv90xoCwK0l6SuerKwstba2qqqqSvF4\nXHV1dQoGg2pra5MkNTQ0aNeuXbpw4YIaGxslSdnZ2ert7Z39yQEsWEk/Tk/rifg4Hcg4s/JxOgDM\nBsIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggP\nAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwR\nHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACY\nIzwAzBEeAOYIDwBzhAeAOcIDwBzhAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAA\nMEd4AJgjPADMER4A5ggPAHOEB4A5wgPAHOEBYI7wADBHeACYIzwAzBEeAOYIDwBzhAeAOcIDwBzh\nAWCO8AAwR3gAmCM8AMwRHgDmCA8Ac4QHgDnCA8Ac4QFgjvAAMEd4AJgjPADMTRmecDis4uJiBQIB\ntbS03HDN9u3bFQgEVFZWpnfeeSftQ84XkUhkrkeYsUzYg5QZ+8iEPUxX0vDE43E99dRTCofDOnXq\nlNrb23X69Olxa7q6uvTee++pv79fzz33nBobG2d14LmUCf+jZMIepMzYRybsYbqShqe3t1d+v18F\nBQXKzs5WbW2tOjs7x6157bXX9OSTT0qS1q1bp4sXL2poaGj2Jgaw4CUNTywWU35+fuKxz+dTLBab\ncs3g4GCaxwSQSbKSPenxeFI6iHMupd+X6vHms+bm5rkeYcYyYQ9SZuwjE/YwHUnD4/V6FY1GE4+j\n0ah8Pl/SNYODg/J6vROO9dk4Abh1Jb3VqqioUH9/vwYGBjQ6OqqOjg7V1NSMW1NTU6MXX3xRktTT\n06Ply5crNzd39iYGsOAlfcWTlZWl1tZWVVVVKR6Pq66uTsFgUG1tbZKkhoYGVVdXq6urS36/X0uW\nLNHBgwdNBgewgLk0euONN1xRUZHz+/1u7969N1zzgx/8wPn9fldaWupOnDiRztOnzVT7eOmll1xp\naam777773Je+9CXX19c3B1Mml8qfhXPO9fb2us997nPu97//veF0qUllD0ePHnXl5eVu7dq17stf\n/rLtgCmYag8fffSRq6qqcmVlZW7t2rXu4MGD9kNOYdu2be7OO+90JSUlk6652es6beG5cuWKKyws\ndGfPnnWjo6OurKzMnTp1atya119/3X3ta19zzjnX09Pj1q1bl67Tp00q+zh27Ji7ePGic+7a/1jz\nbR+p7OHTdRs2bHCbN292r7766hxMOrlU9nDhwgW3Zs0aF41GnXPXLuL5JJU9PP300+6nP/2pc+7a\n/CtWrHBjY2NzMe6k/vrXv7oTJ05MGp7pXNdp+ysTmfKdn1T2sX79ei1btkzStX3Mt68PpLIHSdq/\nf7+eeOIJ3XHHHXMwZXKp7OGVV17R448/nvjAY+XKlXMx6qRS2cOqVas0MjIiSRoZGdHtt9+urKyk\n74CYe+ihh5STkzPp89O5rtMWnkz5zk8q+7jegQMHVF1dbTFaylL9s+js7Ex803y+fdUhlT309/fr\n/Pnz2rBhgyoqKnTo0CHrMZNKZQ/19fV69913lZeXp7KyMu3bt896zBmbznWdtrSm+zs/c+Vm5jl6\n9KheeOEFvfnmm7M40c1LZQ87duzQ3r175fF45K7dchtMlrpU9jA2NqYTJ07oyJEjunz5stavX68H\nHnhAgUDAYMKppbKH3bt3q7y8XJFIRO+//742btyovr4+LV261GDC9LnZ6zpt4Unnd37mUir7kKST\nJ0+qvr5e4XA46cvQuZDKHt5++23V1tZKkoaHh/XGG28oOzt7wtcl5koqe8jPz9fKlSu1ePFiLV68\nWA8//LD6+vrmTXhS2cOxY8f085//XJJUWFioe+65R2fOnFFFRYXprDMxres6XW9AjY2NudWrV7uz\nZ8+6Tz75ZMo3l48fPz7v3pR1LrV9nDt3zhUWFrrjx4/P0ZTJpbKH6333u9+dd59qpbKH06dPu698\n5SvuypUr7n//+58rKSlx77777hxNPFEqe/jRj37kmpqanHPOffjhh87r9br//Oc/czFuUmfPnk3p\nzeVUr+u0fpze1dXl7r33XldYWOh2797tnHPu2Wefdc8++2xizfe//31XWFjoSktL3dtvv53O06fN\nVPuoq6tzK1ascOXl5a68vNxVVlbO5bg3lMqfxafmY3icS20Pv/rVr9yaNWtcSUmJ27dv31yNOqmp\n9vDRRx+5LVu2uNLSUldSUuJefvnluRz3hmpra92qVatcdna28/l87sCBAzO+rj3OzbObewAZj3+B\nEIA5wgPAHOEBYI7wADBHeACYIzwAzP0fycBDEj0HfjwAAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0x10abf10d0>"
       ]
      }
     ],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}