Created
June 14, 2018 10:50
-
-
Save naomi-henderson/25d33e22a0f4bc9ac33f646239126517 to your computer and use it in GitHub Desktop.
CMIP5 multi-model stats notebook from my home linux box
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "cells": [ | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "## Re-create the IPCC_AR5 Figure12.5 from the CMIP5 ts data \n", | |
| "[IPCC-AR5 Figure 12.5](http://www.climatechange2013.org/images/figures/WGI_AR5_Fig12-5.jpg)\n", | |
| "\n", | |
| "### Method: load the monthly CMIP5 ts (surface temperature) files, do some data cleaning and plot the figure\n", | |
| "\n", | |
| "* this notebook is what I use for general multi-model statistics - not just global means. So the models are regridded to\n", | |
| "a common 2x2 degree grid and the global mean is only computed when making the figure\n", | |
| "* the Raw CMIP5 netcdf files were concatenated on our home machine using xarray.mfdataset for each model and scenerio (historical/rcp45/rcp85)\n", | |
| "* saved in zarr format, using to_zarr\n", | |
| "* uploaded to the Google Cloud Storage" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 1, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "import xarray as xr\n", | |
| "import numpy as np\n", | |
| "%matplotlib inline\n", | |
| "from glob import glob\n", | |
| "from os import system\n", | |
| "import pandas as pd\n", | |
| "from matplotlib import pyplot as plt\n", | |
| "from matplotlib import rcParams\n", | |
| "from pathlib import Path\n", | |
| "#import xesmf as xe\n", | |
| "xr.set_options(enable_cftimeindex=True)\n", | |
| "\n", | |
| "cloud = False" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "### In this next cell are the basic functions (all in one cell so you can collapse it on first read)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 2, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "if cloud:\n", | |
| " def listd(fs,path):\n", | |
| " return fs.ls(path)\n", | |
| " def openzarr(path):\n", | |
| " return xr.open_zarr(gcsfs.GCSMap(path))\n", | |
| "else:\n", | |
| " def listd(fs,path):\n", | |
| " return glob(path+'/*/')\n", | |
| " def openzarr(path):\n", | |
| " return xr.open_zarr(path)\n", | |
| "\n", | |
| "def find_models(data_fs,base_path,var,scenario):\n", | |
| " \"\"\"\n", | |
| " Search for all files in path matching the variable and scenario\n", | |
| " returns:\n", | |
| " models: list of model names for given scenario\n", | |
| " \"\"\"\n", | |
| " allmodels = listd(data_fs,base_path)\n", | |
| " fmodels = []\n", | |
| " for model in allmodels:\n", | |
| " run = listd(data_fs,model)\n", | |
| " for sce in run:\n", | |
| " if scenario in sce:\n", | |
| " model = sce.split(\"/\")[-3]\n", | |
| " fmodels += [model]\n", | |
| " \n", | |
| " umodels = sorted(fmodels)\n", | |
| "\n", | |
| " paths = []\n", | |
| " for model in umodels:\n", | |
| " path = base_path + model + '/' + scenario\n", | |
| " paths += [path]\n", | |
| " \n", | |
| " return umodels, paths\n", | |
| "\n", | |
| "def get_datasets(var, umodels, paths, toprint=True):\n", | |
| " \"\"\" \n", | |
| " Load all datasets\n", | |
| " returns: \n", | |
| " ds : list of all models for given scenario\n", | |
| " \"\"\" \n", | |
| " ds = []\n", | |
| " for idx, model in enumerate(umodels):\n", | |
| " path = paths[idx]\n", | |
| " dss = openzarr(path)\n", | |
| " start_date = dss.attrs['start_date']\n", | |
| " nt = dss.time.shape[0]\n", | |
| " dss['time'] = to_enso(start_date,nt)\n", | |
| " \n", | |
| " ds += [dss[var]]\n", | |
| " \n", | |
| " if toprint:\n", | |
| " fstr = '{:2g}: {:18} , {:12} , nt={:5g},{:5.0f}Mb'\n", | |
| " print(fstr.format(idx,model,start_date,nt,dss.nbytes/ 1e6))\n", | |
| " return ds\n", | |
| "\n", | |
| "def find_short(ds, century, umodels, toprint=True):\n", | |
| " \"\"\"\n", | |
| " Identify models which are useful, finding those which do not span the interval\n", | |
| " returns: \n", | |
| " bad_models: list of the bad models\n", | |
| " \"\"\" \n", | |
| "\n", | |
| " slist = century.split('-')\n", | |
| " [start_year, stop_year] = list(map(int, slist))\n", | |
| " \n", | |
| " bad_models =[]\n", | |
| " for idx, dss in enumerate(ds):\n", | |
| " model = umodels[idx]\n", | |
| " \n", | |
| " tfirst = enso2date(dss.time[0].values)\n", | |
| " tlast = enso2date(dss.time[-1].values)\n", | |
| " \n", | |
| " if (int(str(tfirst)[0:4]) > start_year):\n", | |
| " print('trouble with model',model,'since start date is past',start_year)\n", | |
| " bad_models += [model]\n", | |
| "\n", | |
| " if (int(str(tlast)[0:4]) < stop_year):\n", | |
| " print('trouble with model',model,'since stop date is before',stop_year)\n", | |
| " bad_models += [model]\n", | |
| " \n", | |
| " if toprint:\n", | |
| " fstr = '{:2g}: {:18} , {:12} to {:12}'\n", | |
| " print(fstr.format(idx,model,tfirst,tlast,))\n", | |
| " return bad_models\n", | |
| "\n", | |
| "def regrid_all(ds,umodels):\n", | |
| " \"\"\"\n", | |
| " Define common grid and use xESMF to regrid all datasets \n", | |
| " returns:\n", | |
| " data_2x2: a list of datasets on the common grid\n", | |
| " \"\"\" \n", | |
| " # regrid all lon,lat data to a common 2x2 grid\n", | |
| " import xesmf as xe\n", | |
| " ds_out = xr.Dataset({'lat': (['lat'], np.arange(-89,89, 2)),\n", | |
| " 'lon': (['lon'], np.arange(-179,179,2)),\n", | |
| " })\n", | |
| " data_2x2 =[]\n", | |
| " for model,dss in zip(umodels,ds):\n", | |
| " #print(model,'nt=',dss.time.shape[0])\n", | |
| " regridder = xe.Regridder(dss, ds_out, 'bilinear', periodic=True, reuse_weights=True )\n", | |
| " data_2x2 += [regridder(dss)]\n", | |
| " return data_2x2\n", | |
| "\n", | |
| "def concat_all(ds_2x2,umodels):\n", | |
| " \"\"\"\n", | |
| " Concatenates all of the good models into one DataArray\n", | |
| " \"\"\" \n", | |
| " dsall = xr.concat(ds_2x2,dim='model') #,coords=['time','lat','lon'])\n", | |
| " dsall['names'] = ('model',umodels)\n", | |
| " return dsall\n", | |
| "\n", | |
| "def compute_global_mean(ds):\n", | |
| " \"\"\"\n", | |
| " Weights each grid point by the cos(latitude), computes global mean, normalizing by global mean of the weights\n", | |
| " returns:\n", | |
| " list of DataArrays: global mean model by model\n", | |
| " \"\"\" \n", | |
| " coslat = np.cos(np.deg2rad(ds.lat))\n", | |
| " d_ones = xr.ones_like(ds)\n", | |
| " weight_mean = (d_ones*coslat).mean(['lat','lon'])\n", | |
| " ds_globalmean = ((ds * coslat).mean(['lat','lon'])/weight_mean).compute()\n", | |
| " return ds_globalmean\n", | |
| "\n", | |
| "# N.B. Once cftime is working properly the following functions could be replaced\n", | |
| "# (we need to be able to use resample ...)\n", | |
| "def monthly2yearly(century,ds):\n", | |
| " \"\"\"\n", | |
| " converts a DataArray on a monthly grid to one on a yearly grid, replacing the time grid\n", | |
| " returns:\n", | |
| " list of DataArrays: yearly mean model by model\n", | |
| " \"\"\" \n", | |
| " slist = century.split('-')\n", | |
| " [start_year, stop_year] = list(map(int, slist))\n", | |
| " start = to_enso(str(start_year)+'-01-16')[0]\n", | |
| " stop = to_enso(str(stop_year)+'-12-16')[0]\n", | |
| "\n", | |
| " ds_yearly=[]\n", | |
| " for idx, dss in enumerate(ds):\n", | |
| " print('year:',idx)\n", | |
| " dss = dss.sel(time=slice(start, stop))\n", | |
| " num_of_bins = dss.time.shape[0]/12\n", | |
| " dnew = dss.groupby_bins('time', num_of_bins).mean('time').compute()\n", | |
| " dyearly = dnew.rename({'time_bins':'time'}) \n", | |
| " dyearly['time'] = start_year + np.arange(dyearly.time.shape[0])\n", | |
| " ds_yearly += [dyearly]\n", | |
| " return ds_yearly\n", | |
| "\n", | |
| "def to_enso(start_time,nt=1):\n", | |
| " \"\"\"\n", | |
| " Parse the time grid of a Dataset and replace by an enso time grid (months since 1960).\n", | |
| " \"\"\"\n", | |
| " import numpy as np\n", | |
| " # get the reference year from start_time\n", | |
| " ryear,rmonth,rday = start_time[0:10].split('-')\n", | |
| " return (int(ryear)-1960)*12 + int(rmonth) - 0.5 + np.arange(0,nt)\n", | |
| "\n", | |
| "def enso2date(T0,ryear=1960,leap=True):\n", | |
| " \"\"\"\n", | |
| " Print the date corresponding to an enso-time (months since 1960). The reference year can be changed.\n", | |
| " \"\"\"\n", | |
| " norm = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]\n", | |
| " iy = ryear + int(T0/12)\n", | |
| " if T0 < 0:\n", | |
| " iy = iy - 1\n", | |
| " res = T0 - (iy - ryear)*12\n", | |
| " im = int(res) + 1\n", | |
| " if im == 13:\n", | |
| " im = 1\n", | |
| " iy = iy + 1\n", | |
| " if leap & (im == 2) & (iy % 4 == 0 ): \n", | |
| " id = 1 + int(29 * (res - int(res)))\n", | |
| " else:\n", | |
| " id = 1 + int(norm[im-1] * (res - int(res)))\n", | |
| " return str(iy)+'/'+str(im)+'/'+str(id) " | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "## Connect to Dask Distributed Cluster when done debugging" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 3, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "#from dask.distributed import Client, progress\n", | |
| "\n", | |
| "#from dask_kubernetes import KubeCluster\n", | |
| "#cluster = KubeCluster(n_workers=10)\n", | |
| "#cluster" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 4, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "#client = Client(cluster)\n", | |
| "#client" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "## Lets get started. We need to specify where our data lives in the google file system" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 5, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "if cloud:\n", | |
| " import gcsfs\n", | |
| " base_path = 'pangeo-data/CMIP5-ts/'\n", | |
| " data_fs = gcsfs.GCSFileSystem(project='pangeo-181919', token='anon', access='read_only')\n", | |
| "else: \n", | |
| " base_path = '/d1/nhn2/zarr/CMIP5-ts/'\n", | |
| " data_fs = ''" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "## For each scenario and time interval:\n", | |
| "* make a list of available models\n", | |
| "* load the datasets for all models\n", | |
| "* eliminate the models which do not have the full time interval\n", | |
| "* calculate the annual means\n", | |
| "* regrid to a global 2x2 degree grid (until we have xesmf, just calculate the global mean)\n", | |
| "* concatenate the models" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 6, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "SCENARIO= historical TIME RANGE 1861-2005\n", | |
| "total number of models available: 49 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CESM1-FASTCHEM', 'CESM1-WACCM', 'CMCC-CESM', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5-2', 'CSIRO-Mk3-6-0', 'CSIRO-Mk3L-1-2', 'CanCM4', 'CanESM2', 'FGOALS-g2', 'FGOALS-s2', 'FIO-ESM', 'GFDL-CM2p1', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadCM3', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC4h', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MPI-ESM-P', 'MRI-CGCM3', 'MRI-ESM1', 'NorESM1-M', 'NorESM1-ME', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "trouble with model CanCM4 since start date is past 1861\n", | |
| "trouble with model MIROC4h since start date is past 1861\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 47 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CESM1-FASTCHEM', 'CESM1-WACCM', 'CMCC-CESM', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5-2', 'CSIRO-Mk3-6-0', 'CSIRO-Mk3L-1-2', 'CanESM2', 'FGOALS-g2', 'FGOALS-s2', 'FIO-ESM', 'GFDL-CM2p1', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadCM3', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MPI-ESM-P', 'MRI-CGCM3', 'MRI-ESM1', 'NorESM1-M', 'NorESM1-ME', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "year: 7\n", | |
| "year: 8\n", | |
| "year: 9\n", | |
| "year: 10\n", | |
| "year: 11\n", | |
| "year: 12\n", | |
| "year: 13\n", | |
| "year: 14\n", | |
| "year: 15\n", | |
| "year: 16\n", | |
| "year: 17\n", | |
| "year: 18\n", | |
| "year: 19\n", | |
| "year: 20\n", | |
| "year: 21\n", | |
| "year: 22\n", | |
| "year: 23\n", | |
| "year: 24\n", | |
| "year: 25\n", | |
| "year: 26\n", | |
| "year: 27\n", | |
| "year: 28\n", | |
| "year: 29\n", | |
| "year: 30\n", | |
| "year: 31\n", | |
| "year: 32\n", | |
| "year: 33\n", | |
| "year: 34\n", | |
| "year: 35\n", | |
| "year: 36\n", | |
| "year: 37\n", | |
| "year: 38\n", | |
| "year: 39\n", | |
| "year: 40\n", | |
| "year: 41\n", | |
| "year: 42\n", | |
| "year: 43\n", | |
| "year: 44\n", | |
| "year: 45\n", | |
| "year: 46\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_48x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_240x480_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_56x64_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_60x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_108x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_73x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_143x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_120x180_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "SCENARIO= historical TIME RANGE 1850-1860\n", | |
| "\n", | |
| " same scenario, re-use ds \n", | |
| "\n", | |
| "trouble with model CSIRO-Mk3L-1-2 since start date is past 1850\n", | |
| "trouble with model CanCM4 since start date is past 1850\n", | |
| "trouble with model GFDL-CM2p1 since start date is past 1850\n", | |
| "trouble with model GFDL-CM3 since start date is past 1850\n", | |
| "trouble with model GFDL-ESM2G since start date is past 1850\n", | |
| "trouble with model GFDL-ESM2M since start date is past 1850\n", | |
| "trouble with model HadCM3 since start date is past 1850\n", | |
| "trouble with model HadGEM2-AO since start date is past 1850\n", | |
| "trouble with model HadGEM2-CC since start date is past 1850\n", | |
| "trouble with model HadGEM2-ES since start date is past 1850\n", | |
| "trouble with model MIROC4h since start date is past 1850\n", | |
| "trouble with model MRI-ESM1 since start date is past 1850\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 37 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CESM1-FASTCHEM', 'CESM1-WACCM', 'CMCC-CESM', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5-2', 'CSIRO-Mk3-6-0', 'CanESM2', 'FGOALS-g2', 'FGOALS-s2', 'FIO-ESM', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MPI-ESM-P', 'MRI-CGCM3', 'NorESM1-M', 'NorESM1-ME', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "year: 7\n", | |
| "year: 8\n", | |
| "year: 9\n", | |
| "year: 10\n", | |
| "year: 11\n", | |
| "year: 12\n", | |
| "year: 13\n", | |
| "year: 14\n", | |
| "year: 15\n", | |
| "year: 16\n", | |
| "year: 17\n", | |
| "year: 18\n", | |
| "year: 19\n", | |
| "year: 20\n", | |
| "year: 21\n", | |
| "year: 22\n", | |
| "year: 23\n", | |
| "year: 24\n", | |
| "year: 25\n", | |
| "year: 26\n", | |
| "year: 27\n", | |
| "year: 28\n", | |
| "year: 29\n", | |
| "year: 30\n", | |
| "year: 31\n", | |
| "year: 32\n", | |
| "year: 33\n", | |
| "year: 34\n", | |
| "year: 35\n", | |
| "year: 36\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_48x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_240x480_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_60x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_108x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_143x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_120x180_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "SCENARIO= rcp45 TIME RANGE 2100-2300\n", | |
| "total number of models available: 43 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5', 'CSIRO-Mk3-6-0', 'CSIRO-Mk3L-1-2', 'CanCM4', 'CanESM2', 'FGOALS-g2', 'FIO-ESM', 'GFDL-CM2p1', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadCM3', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC4h', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MRI-CGCM3', 'NorESM1-M', 'NorESM1-ME', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "trouble with model ACCESS1-0 since stop date is before 2300\n", | |
| "trouble with model ACCESS1-3 since stop date is before 2300\n", | |
| "trouble with model BNU-ESM since stop date is before 2300\n", | |
| "trouble with model CCSM4 since stop date is before 2300\n", | |
| "trouble with model CESM1-BGC since stop date is before 2300\n", | |
| "trouble with model CESM1-CAM5-1-FV2 since stop date is before 2300\n", | |
| "trouble with model CMCC-CM since stop date is before 2300\n", | |
| "trouble with model CMCC-CMS since stop date is before 2300\n", | |
| "trouble with model CSIRO-Mk3-6-0 since stop date is before 2300\n", | |
| "trouble with model CanCM4 since stop date is before 2300\n", | |
| "trouble with model FGOALS-g2 since stop date is before 2300\n", | |
| "trouble with model FIO-ESM since stop date is before 2300\n", | |
| "trouble with model GFDL-CM2p1 since stop date is before 2300\n", | |
| "trouble with model GFDL-CM3 since stop date is before 2300\n", | |
| "trouble with model GFDL-ESM2G since stop date is before 2300\n", | |
| "trouble with model GFDL-ESM2M since stop date is before 2300\n", | |
| "trouble with model GISS-E2-H-CC since stop date is before 2300\n", | |
| "trouble with model GISS-E2-R-CC since stop date is before 2300\n", | |
| "trouble with model HadCM3 since stop date is before 2300\n", | |
| "trouble with model HadGEM2-AO since stop date is before 2300\n", | |
| "trouble with model HadGEM2-CC since stop date is before 2300\n", | |
| "trouble with model HadGEM2-ES since stop date is before 2300\n", | |
| "trouble with model IPSL-CM5A-MR since stop date is before 2300\n", | |
| "trouble with model IPSL-CM5B-LR since stop date is before 2300\n", | |
| "trouble with model MIROC-ESM since stop date is before 2300\n", | |
| "trouble with model MIROC-ESM-CHEM since stop date is before 2300\n", | |
| "trouble with model MIROC4h since stop date is before 2300\n", | |
| "trouble with model MIROC5 since stop date is before 2300\n", | |
| "trouble with model MPI-ESM-MR since stop date is before 2300\n", | |
| "trouble with model MRI-CGCM3 since stop date is before 2300\n", | |
| "trouble with model NorESM1-ME since stop date is before 2300\n", | |
| "trouble with model bcc-csm1-1 since stop date is before 2300\n", | |
| "trouble with model bcc-csm1-1-m since stop date is before 2300\n", | |
| "trouble with model inmcm4 since stop date is before 2300\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 9 \n", | |
| " ['CESM1-CAM5', 'CNRM-CM5', 'CSIRO-Mk3L-1-2', 'CanESM2', 'GISS-E2-H', 'GISS-E2-R', 'IPSL-CM5A-LR', 'MPI-ESM-LR', 'NorESM1-M']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "year: 7\n", | |
| "year: 8\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_56x64_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "SCENARIO= rcp45 TIME RANGE 2006-2099\n", | |
| "\n", | |
| " same scenario, re-use ds \n", | |
| "\n", | |
| "trouble with model CanCM4 since stop date is before 2099\n", | |
| "trouble with model GFDL-CM2p1 since stop date is before 2099\n", | |
| "trouble with model HadCM3 since stop date is before 2099\n", | |
| "trouble with model MIROC4h since stop date is before 2099\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 39 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5', 'CSIRO-Mk3-6-0', 'CSIRO-Mk3L-1-2', 'CanESM2', 'FGOALS-g2', 'FIO-ESM', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MRI-CGCM3', 'NorESM1-M', 'NorESM1-ME', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "year: 7\n", | |
| "year: 8\n", | |
| "year: 9\n", | |
| "year: 10\n", | |
| "year: 11\n", | |
| "year: 12\n", | |
| "year: 13\n", | |
| "year: 14\n", | |
| "year: 15\n", | |
| "year: 16\n", | |
| "year: 17\n", | |
| "year: 18\n", | |
| "year: 19\n", | |
| "year: 20\n", | |
| "year: 21\n", | |
| "year: 22\n", | |
| "year: 23\n", | |
| "year: 24\n", | |
| "year: 25\n", | |
| "year: 26\n", | |
| "year: 27\n", | |
| "year: 28\n", | |
| "year: 29\n", | |
| "year: 30\n", | |
| "year: 31\n", | |
| "year: 32\n", | |
| "year: 33\n", | |
| "year: 34\n", | |
| "year: 35\n", | |
| "year: 36\n", | |
| "year: 37\n", | |
| "year: 38\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_240x480_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_56x64_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_60x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_143x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_120x180_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "SCENARIO= rcp85 TIME RANGE 2006-2099\n", | |
| "total number of models available: 40 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CMCC-CESM', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5', 'CSIRO-Mk3-6-0', 'CanESM2', 'FGOALS-g2', 'FGOALS-s2', 'FIO-ESM', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MRI-CGCM3', 'MRI-ESM1', 'NorESM1-M', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 40 \n", | |
| " ['ACCESS1-0', 'ACCESS1-3', 'BNU-ESM', 'CCSM4', 'CESM1-BGC', 'CESM1-CAM5', 'CESM1-CAM5-1-FV2', 'CMCC-CESM', 'CMCC-CM', 'CMCC-CMS', 'CNRM-CM5', 'CSIRO-Mk3-6-0', 'CanESM2', 'FGOALS-g2', 'FGOALS-s2', 'FIO-ESM', 'GFDL-CM3', 'GFDL-ESM2G', 'GFDL-ESM2M', 'GISS-E2-H', 'GISS-E2-H-CC', 'GISS-E2-R', 'GISS-E2-R-CC', 'HadGEM2-AO', 'HadGEM2-CC', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'IPSL-CM5A-MR', 'IPSL-CM5B-LR', 'MIROC-ESM', 'MIROC-ESM-CHEM', 'MIROC5', 'MPI-ESM-LR', 'MPI-ESM-MR', 'MRI-CGCM3', 'MRI-ESM1', 'NorESM1-M', 'bcc-csm1-1', 'bcc-csm1-1-m', 'inmcm4']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "year: 7\n", | |
| "year: 8\n", | |
| "year: 9\n", | |
| "year: 10\n", | |
| "year: 11\n", | |
| "year: 12\n", | |
| "year: 13\n", | |
| "year: 14\n", | |
| "year: 15\n", | |
| "year: 16\n", | |
| "year: 17\n", | |
| "year: 18\n", | |
| "year: 19\n", | |
| "year: 20\n", | |
| "year: 21\n", | |
| "year: 22\n", | |
| "year: 23\n", | |
| "year: 24\n", | |
| "year: 25\n", | |
| "year: 26\n", | |
| "year: 27\n", | |
| "year: 28\n", | |
| "year: 29\n", | |
| "year: 30\n", | |
| "year: 31\n", | |
| "year: 32\n", | |
| "year: 33\n", | |
| "year: 34\n", | |
| "year: 35\n", | |
| "year: 36\n", | |
| "year: 37\n", | |
| "year: 38\n", | |
| "year: 39\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_192x288_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_48x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_240x480_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_60x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_108x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_143x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_64x128_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_160x320_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_120x180_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "SCENARIO= rcp85 TIME RANGE 2100-2300\n", | |
| "\n", | |
| " same scenario, re-use ds \n", | |
| "\n", | |
| "trouble with model ACCESS1-0 since stop date is before 2300\n", | |
| "trouble with model ACCESS1-3 since stop date is before 2300\n", | |
| "trouble with model BNU-ESM since stop date is before 2300\n", | |
| "trouble with model CCSM4 since stop date is before 2300\n", | |
| "trouble with model CESM1-BGC since stop date is before 2300\n", | |
| "trouble with model CESM1-CAM5 since stop date is before 2300\n", | |
| "trouble with model CESM1-CAM5-1-FV2 since stop date is before 2300\n", | |
| "trouble with model CMCC-CESM since stop date is before 2300\n", | |
| "trouble with model CMCC-CM since stop date is before 2300\n", | |
| "trouble with model CMCC-CMS since stop date is before 2300\n", | |
| "trouble with model CanESM2 since stop date is before 2300\n", | |
| "trouble with model FGOALS-g2 since stop date is before 2300\n", | |
| "trouble with model FGOALS-s2 since stop date is before 2300\n", | |
| "trouble with model FIO-ESM since stop date is before 2300\n", | |
| "trouble with model GFDL-CM3 since stop date is before 2300\n", | |
| "trouble with model GFDL-ESM2G since stop date is before 2300\n", | |
| "trouble with model GFDL-ESM2M since stop date is before 2300\n", | |
| "trouble with model GISS-E2-H-CC since stop date is before 2300\n", | |
| "trouble with model GISS-E2-R-CC since stop date is before 2300\n", | |
| "trouble with model HadGEM2-AO since stop date is before 2300\n", | |
| "trouble with model HadGEM2-CC since stop date is before 2300\n", | |
| "trouble with model IPSL-CM5A-MR since stop date is before 2300\n", | |
| "trouble with model IPSL-CM5B-LR since stop date is before 2300\n", | |
| "trouble with model MIROC-ESM since stop date is before 2300\n", | |
| "trouble with model MIROC-ESM-CHEM since stop date is before 2300\n", | |
| "trouble with model MIROC5 since stop date is before 2300\n", | |
| "trouble with model MPI-ESM-MR since stop date is before 2300\n", | |
| "trouble with model MRI-CGCM3 since stop date is before 2300\n", | |
| "trouble with model MRI-ESM1 since stop date is before 2300\n", | |
| "trouble with model NorESM1-M since stop date is before 2300\n", | |
| "trouble with model bcc-csm1-1 since stop date is before 2300\n", | |
| "trouble with model bcc-csm1-1-m since stop date is before 2300\n", | |
| "trouble with model inmcm4 since stop date is before 2300\n", | |
| "\n", | |
| " number of good models with data in the specified time range: 7 \n", | |
| " ['CNRM-CM5', 'CSIRO-Mk3-6-0', 'GISS-E2-H', 'GISS-E2-R', 'HadGEM2-ES', 'IPSL-CM5A-LR', 'MPI-ESM-LR']\n", | |
| "\n", | |
| " calculating annual means\n", | |
| "year: 0\n", | |
| "year: 1\n", | |
| "year: 2\n", | |
| "year: 3\n", | |
| "year: 4\n", | |
| "year: 5\n", | |
| "year: 6\n", | |
| "\n", | |
| " regridding to 2x2 grid\n", | |
| "Reuse existing file: bilinear_128x256_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_90x144_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_145x192_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x96_89x179_peri.nc\n", | |
| "Reuse existing file: bilinear_96x192_89x179_peri.nc\n", | |
| "\n", | |
| " concatenating time series\n", | |
| "CPU times: user 11min 14s, sys: 24min 53s, total: 36min 8s\n", | |
| "Wall time: 2min 46s\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "%%time\n", | |
| "var = 'ts'\n", | |
| "recompute_all = True\n", | |
| "plot_global = True\n", | |
| "plot_nino34 = False\n", | |
| "save_netcdf = True\n", | |
| "\n", | |
| "if recompute_all:\n", | |
| " \n", | |
| " # Each of the following scenarios and time periods has a different subset of CMIP5 models available\n", | |
| " # These subsets are different for each variable chosen\n", | |
| " # YOU CAN PICK AND CHOOSE WHICH TO CALCULATE \n", | |
| " all = []\n", | |
| " all += [['historical','1861-2005']]\n", | |
| " all += [['historical','1850-1860']]\n", | |
| " all += [['rcp45','2100-2300']]\n", | |
| " all += [['rcp45','2006-2099']]\n", | |
| " all += [['rcp85','2006-2099']]\n", | |
| " all += [['rcp85','2100-2300']]\n", | |
| " \n", | |
| " ds_master = [] # list of datasets for each scenario,century\n", | |
| " scenario_last = ''\n", | |
| " for scenario,century in all:\n", | |
| " print('SCENARIO=',scenario,'TIME RANGE',century)\n", | |
| " if scenario in scenario_last:\n", | |
| " print('\\n same scenario, re-use ds \\n')\n", | |
| " ds = ds2; models = models2\n", | |
| " else:\n", | |
| " models, paths = find_models(data_fs,base_path,var,scenario)\n", | |
| " print('total number of models available:',len(models),'\\n',models)\n", | |
| " ds = get_datasets(var, models, paths, toprint=False)\n", | |
| " ds2 = ds.copy(); models2 = models.copy()\n", | |
| " scenario_last = scenario\n", | |
| " \n", | |
| " bad_models = find_short(ds, century, models, toprint=False)\n", | |
| "\n", | |
| " bad_models = sorted(list(set(bad_models))) \n", | |
| " for model in bad_models:\n", | |
| " idx = models.index(model)\n", | |
| " del models[idx],ds[idx]\n", | |
| "\n", | |
| " print('\\n number of good models with data in the specified time range:',len(models),'\\n',models)\n", | |
| " \n", | |
| " print('\\n calculating annual means')\n", | |
| " ds_yearly = monthly2yearly(century,ds)\n", | |
| "\n", | |
| " print('\\n regridding to 2x2 grid')\n", | |
| " ds_temp = regrid_all(ds_yearly,models)\n", | |
| " \n", | |
| " print('\\n concatenating time series')\n", | |
| " dsall = concat_all(ds_temp,models)\n", | |
| "\n", | |
| " sctype = scenario+':'+century\n", | |
| "\n", | |
| " if save_netcdf:\n", | |
| " dsall.to_netcdf('ts-'+sctype+'.nc',encoding={'time':{'dtype':'float32'},'lon':{'dtype':'float32'},'lat':{'dtype':'float32'}})\n", | |
| " \n", | |
| " dsall.attrs = [('sctype',sctype)] \n", | |
| " ds_master += [dsall.to_dataset(name=var)] " | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": {}, | |
| "source": [ | |
| "## Plot figure" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 8, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAm8AAAGBCAYAAAAwkFbnAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XlcVOX+wPHPw7AKguAGioqy6tWL5pKmJaVoaS6p5c0FQS/9WsztZaapuaS3NLuapd68aWpXXLCba1mZ+0JuYWqIioCYmCCLsjPM+f1xhrkM64Cs+rxfr/OCOec5z3nOmcPwnWc7QlEUJEmSJEmSpLrBrKYLIEmSJEmSJJlOBm+SJEmSJEl1iAzeJEmSJEmS6hAZvEmSJEmSJNUhMniTJEmSJEmqQ2TwJkmSJEmSVIfI4E2SJEmSJKkOkcGbJEmSJElSHSKDN0mSJEmSpDpEBm+SJEmSVAOEEE8KIcKEELuFEE41XR6p7hDy8ViSJEmSVP2EED8BY4AeQBdFUebUcJGkOkLWvEnlIoToIYTYLoS4LYTIEULcE0L8JIQYJ4TQ6NMECiEU/eJVTB5+Bbb3LbSPR4F0BfNRhBAPhBAXhBAThRDmJeRXcEmpxPOeL4So0Ded/H0LlrkmylEVhBDjhRDX9PdCpV3vhyWEGCqEmFbT5agt6vL1qMy/n1pIAAqQRwn/j4UQLwohdhb4zP1TCPFfIUSfai2pVKvI4E0ymRBiCnACcALeBfoC44GrwBrgxUK7PADGFpNVgH6bqV5G/WY6HDgNfAa8X0y6Sfp0+UvfchxDKichRDNgLXASeI7adb2HAnUyWKki8nrUTu8De4A3gOUFNwghzIUQXwO7gWxgCuAPzAQaAz8KIWyrt7hSbfEofpORqoAQ4hngn8DniqJMKrR5lxDin0DhD5L/AmOEEO8r+vZ5IYQNahD2DRBo4uHDFUW5rv/9R33t3BSKBnARiqKEmZinVEFCCCtFUbIBT0ADbFQU5XgNF+uxU+B9eKzLUJcpinISeLKEzauBV4GRiqKEFlh/BPhKCDFKUZT0qi6jVDvJmjfJVDOBJGBGcRsVRYlSFOW3Qqu/BloBvQqsewn1H/43D1GWM0B9IUSTh8jDQAjxqhDiihAiSwhxUQgxWAhxWAhx2IR9nxdCnBJCZAohUvXNG94lJG8rhDgkhMgQQsQLIRYKIcwK5OUhhPhaCBGtz++GEGKNEMKxAufkJYT4VghxV39eN4UQoflNT0KIDUKImGL2MzrvAk1W7YUQPwgh0oDtQogNQH66n/VpNpTnPIQQvvoy3tOnixRCzComzW4hRLI+zQkhxNNlnPsGYBzQvEATekyB7WXmWeC8ffTnna6/hkH67WP190ya/j11L2H/DqW95xUoT+H3ocxrbcL1eKh7oaLvU6HzL/U+0GsthNinv+axQoj3K/L3U+A8PEvLT5/WpM+Ghzn/Eq5JHyAYWFIocDNQFCWkovlLdZ8M3qQyCbUvmx/wo6IoWeXYNRY4inHTaQDwLZD2EEVqjdpHpHAem4UQefp/AiFCiJZlZSSE8Ac2A1dQawSXASuAIn31itn3eWCfvhwjUZs+2gPHhRDNi9llJ3AAtQkrBJiLce1hM+AWaq1if2Ah0Af4rqyyFGMv0Fxfpv6owXc2Ff+b34X6jX8wavPOB6jN1ABvoTZTf6B/XeZ5CCG6AacAd2AqMBC1Zte1QJonUJtknVD/kQ0H7gEHhBCdSynrB/pjJfC/JvSXKphnKOp7PBQ4B6wXQvwD9brOBIIAb9T3szilvucVKE/h98GUe6bE61FBRmV4iPfJpPuggG+Bg6jXciewADUozVfev59S8zP1s+Fhzr8Us4AMYGkF95cedYqiyEUupS5AU9ROtR+amD5Qn94DtU9cMmANuABa1H4bfvo0fQvvU0w+3qhN/I7A/6EGbjsLpOuE+sE6COiN+uF9F/gDaFJGWU8Cl9CPvNave0J/3MMF1s1X/1yM9j0LXAPMC6xrDeQC/yy8LzCz0P7/Ru3716CEspmj1loqQKeSylHMfo30+wwuJc0GIKaY9YeLO29gcjFp++q3+ZVRnuLO4ygQB9QrZb+fgQjAssA6jX7dzjKOuQG4VdE8C5x3QIF1jvr79x5gX2D9JH3aVuV9zytQniLvQ1nXurTrURn3wkO+T6bcB/nHDSq0/iLqF8ryXguT8sP0z4YKn38J5XZE/Yz7T3n3lcvjs8iaN6mqhQJWqIHVaOAO6oddeVxBDYiSUPuBbEYNCgFQFOVXRVGmK4qyR1GUI4qirACeRw06C/fPM9DXKHYBvlEUxTCCU1GU80B0aQUSakfhJ4BtiqJoC+wbjTqoo3cxu20v9HorYIdaW4cQwlII8Z6+mSZTf87H9GlLaootzj3gBvCRECJYCOFZjn1L8q2pCcs6DyFEPaAnsFlRlIwS8rBBvYahgE6onbfNUUfnHQCeKe8JVDDP7/N/URQlGfVLQZiiKPcLpLmi/9mimP1LfM8rWB6j96ES75nyMJThYd4nU+6DQvYVen0JMNSuV+BalJifqZ8NVXGfAn9FrSG/WIF9pceEHLAgmeIekInaf61cFEV5IITYidp06ob6Qa0TQpQnm5dQm0MeALGKCU23iqKcF0JcBbqWkqwRYIH6D7mwP8s4hCPqB3R8MdvuUPy1Kpxn/uv8JtYPgbdRm3tOop6vK+rAD+syymOgKIqib/KZr8+zoRAiGvhYUZQ1puZTSHHnWZKyzsMR9Z/TrVLycEKtvZirX4oQQpgpiqIrR7kqkmdyoSQ5JayD4t+j0t7z6AqUp/D7UCn3TDkVLMPDvE+m3AcFJRV6nY3xOZb3WpSWn6mfDVVxnzoUcxxJMiKDN6lMiqJo9R10/UXFRpdtQv2Wa4Y6eqq8Lin/G21aHvlzKJUkEfXbeXEDH5oCN0vZN1mft3Mx25xRA97i8rxR6DWozbsAfwM2KYqyKD+BEMKulDKUSFGUG0CAUKNkX2AisFoIEaMoyvdAFmBZzK4NSyh7eeaWK+s8kgEd/wtai5OiT7MK9f4pWqDy/UOsqjzLUtp7XpHyFH4fKuOeeZh74WGuqSn3QXlU2t8Ppn82VMU9lR+0FdfvT5IAOWBBMt1HqB/mHxe3UQjRWgjx1xL2/Qm1+ehfiqJcrqLyFS5PF9SOxb+UlEZRlDzUfmvDRYGqQH0n49al5a+oQ/TPAS/rm1jy920FPIXaobuwVwq9/hvqYIdL+tf1UP9hFBRUWjnKoqjC+d8cX+31P2OBpkKIRvlphTpisjKa2ko9D30T2XHUaWRsSih3OmqTly9wXlGUs4WXMsqQDRjlXQl5VkSJ73kllcfUe6bI9SigwvfCw5yDKfdBOVXa34+pnw1VdE+dR63dHCeEsCq8UQhRTwhRWouC9BiQNW+SSRRFOSrUGdr/KYRoi9rJ+SZq00cf4O/AKKDwdCH5H4QVqXEziRBiM2oT1HnUb8KdUEdr/YE6oW9p5gE/At8KIdaiNpfMR236LOsb81zUGsW9QojVqH2ZFgCpwCfFpA/WT0VwBnU03N+B+Yqi5D+ZYD/qB/ZF4DowDDUQLBd9EP0psE2fjwZ18IcWdXQdqH10PkAdoftP1POehVrj8LBMOY/pqAHuKSHEJ6hNZ22AjoqivK1PMw21Q/sPQoh1qP/QGqH2NdQoijKzlDL8DjgJId5A/SecpSjKxYfMsyLKes8ftjym3jMlXQ94+HvhYc7BlPvAVJXy91OAqZ8NlXpPKYqSK4R4E9gBnBZCrECtvXVA7SMYiHrdzlT4zKS6r6ZGSsilbi6oH4ahqB9Q+YMIfkR9Pp+ZPk0ghUaOFpOPH6aPNi0xH326WahBY6q+THGoM/+7mHhOo4BI1NqJy6h97H4Fvi2QZj7FjPJEHRhxCrVPYCrqNArehdLM159He+CQPu0d1H+YZgXSNULt0J6sXzaj9tlTgMDSylHoeE2AjahPvsjQv0dHgP6F0g1FrfXLBC4A/Sh5hKF5MccpdrSpKeehT9cJdXb5FH0ZrgDvFsqrrT6vu/r35xbqjPMDyrgGtsAW/te8HVOePEs6byCGQqMAKXQvl+c9r4TymHqtS7welXQvVOh9MuU+KOXcNxR6X029FiblZ+pnw8OefynXpTvq50kiar/KONRBEFMBh4rmK5dHY5EPppekQoQQrqjf3BcrivJBWeklqTAhxHzUmhsLpcBoZKluk58NUm1R7j5vQghXIcRnQp1VPkOoM1W7FUrTRQixVj9kO0OoM5NvFkKU2o9IkqqbEMJGqLOwDxdC9BbqDPo/odZYfVnDxZMkqYbIzwapNqtInzcP1E6451A7avYrJs3fgL8AK1Grmpuj9g86K4ToqChKXMWKK0mVLg91dOjnqAMy8jsgv6woSnmmx5Ak6dEiPxukWqvczaYF56wRQvwddcbw1oqixBRI01hRlIRC+7VC7VS+SFGUwg8UlyRJkiRJkkxQ7mZTxYQ5awoHbvp1sajP1qusOX0kSZIkSZIeO9U2z5t+eokmqM97kyRJkiRJkiqgWuZ50z/r7V+oNW/rSkn3GvAagK2tbWcfH5/qKJ4kSVKVi4tTu/q2aFHcI1ClR02F32+djrzMTLQpKZCXBxp1DnBhJufUr82UvDyuJyRkZObl2VbH8aprkt7PUecHG6ioD3culqIoa1Hn56JLly7K2bNVMeG5JElS9ZsyZQoAK1asqOGSSNWhIu+3oijcWrGCzKgozB0cEBoNuXfvYtmsGUKjKTsDqcZoU1N5asWKlLJTVo4qD96EEB+i1qaNUxTlx6o+niRJkiTVRRm//05WVBSWzs7kP5XLStbUSsWo0uBNCDEbmAlMUhTl66o8liRJkiTVVXlpadwNDUVjZ2cI3CSpJFUWvAkhJgGLgNmKopT1fElJkiRJeqykX7mCmZUVFo0aEb92LdqUFCwbN67pYkl1QIWCNyHECP2vnfU/XxBCJAAJiqIcEUL8DViB+qDgg0KI7gV2v68oyu8VLrEkSZIk1XG6nBz+3LCBvPR0zOrVQ5eTg0WjRjVdLKmOqGjNW2ih16v1P4+gPqT5eUDofz5fKG1+mgq7f/8+d+/eJTc392GykaQ6y9bWFldXV8zkCDRJqpMenDlDXkYG5g0bQl4e5vb2NV0kqQ6pUPCmKEqpDfKKogQCgRXJuyz379/nzz//pHnz5tjY2Mi+AdJjR6fT8ccff5CYmEiTJk1qujiSJJlIl5VFbmIiwsqKxF270NjbY2ZhARYWNV00qY6prqlCKs3du3dp3rw59erVq+miSFKNMDMzo2nTpsTGxsrgTZLqkOSffyZx9240NjZgbo65jU1NF0mqo+pc8Jabm4uNvOGlx5yFhQVarbamiyFJkonyMjJIPngQS2dnUBTMrKxqukhSHVYnO8zIplLpcSf/BiSp9tPev8/d7dvR5eaSevIkSm4uZpaWMnCTHlqdq3mTJEmSpNpMl5tLVnQ0WbGxJFtYoOh03A8Lw9zRsaaLJj0i6mTNmyRJkiTVVgmhoWTGxGBmbY2liwspBw8izM0xs7Ss6aJJjwgZvNUSbm5ufPnll9V6zMOHDyOEqLS+U0IIDhw4UGqapUuXMmHChEo5Xk0LDAxkzJgxlZ4+LS2NZs2aERsb+zDFkySpBmTfvs39sDA09eqBmRlmFhZYtWiBhax1kyqRDN6qmJ+fH3PmzCmyfsyYMQQGBhpenzlzhtGjR5eZ34EDByqtv9NTTz1FfHw85ubV03r+4MEDlixZwsyZMw3rpkyZgru7OzY2Njg7OxMQEMDdu3cN2+fPn48QosgyYMCAailzTbCzsyM4OJiFCxfWdFEkSSoHXVYWCTt2IMzNocDntJDzMUqVTN5RtUTjxo2rdRRtdnY2lpaWODs7V9sxQ0JC8PLywtPT07DO19eXr776ioiICHbt2kV0dLRREDt9+nTi4+MNy40bN6hXrx7Dhg2rtnLXhICAALZs2UJqampNF0WSJBMoWi23Vq4kMyoKcyenmi6O9IiTwVstUbDZNCsri+DgYJo0aYKNjQ0+Pj7s3LmTmJgY/P39AQw1UBs2bADg2rVr9OvXDxsbG5o0acI777xj1Bzq5ubGkiVLGDZsGPXq1WP16tXFNptu3bqVDh06YGVlRYsWLVi8eDGgBnsBAQG0aNECW1tbnnjiiTKbSAvbsWMHAwcONFoXFBTEM888g5ubG08++STvvPMOx48fN2y3s7PD2dnZsJw8eRJFUXjllVdKPE5gYCCjR49m1qxZODk50bRpU9avX09KSgovv/wydnZ2/PWvfyU8PNxov48//pgWLVpgZWVF9+7dOX36tNH2Tz/9lKZNm+Lg4MC0adPQ6XRG2zMyMnjrrbdo3LgxDRo04MUXXyQmJqbEcq5YsYLWrVtjZWWFq6sr8+fPN2xzd3enRYsW7N27t8T9JUmqPTKuXCH71i0smzSRo8GlKvfIjDb18/MrM82LL77I9OnTDekDAwMJDAwkMTGRESNGlLE3hvRVbeXKlZw7d47vv/8eJycnrly5Ygimtm/fziuvvEJ8fDwADg4O5OXlMWTIENzd3Tl9+jRxcXGMGzcOR0dH3nvvPUO+S5cu5aOPPmLZsmVYWlpy/fp1o+P++OOPBAQE8NFHHzFo0CASExP5/Xf1MbRarRYvLy+mTZuGnZ0d//nPfxg6dChRUVE0bdq0zHNSFIWwsDCmTZtWYpqUlBRCQkLo1atXiWk2bNjASy+9hH0Zj5LZs2cPb731FmFhYYSGhvL6668TGhrKqFGjWLx4MTNnzmT8+PGcP38eUGsF58+fz9q1a3niiSf45JNPGDBgADdu3MDe3p5Dhw4xffp0PvvsM/z8/Fi1ahXffvstQ4YMMRzz9ddf586dO3z//ffUr1+fJUuW8OKLL3LhwgU0Go1R+U6fPs28efPYvn07bdu25datW9y4ccMoTefOnTlx4oRJzemSJNWslCNH5BQgUrWRNW/VYOnSpdjZ2Rkt27ZtKzF9XFwcnTp1onPnzrRu3ZoXXniB5557Do1Gg6O+02t+TZSNjQ0//fQT0dHRbNy4kQ4dOjBgwAAWLFjA8uXLjfIdPHgwwcHBtGnTBldX1yLHXbx4Ma+99hrTpk3D09OTHj16GAYX2NraMmfOHDp27IiHhwfz58/H1dWVH374waRrkJycTFpaGi4uLkW2rV69Gjs7OxwdHYmNjS3x2ty6dYuDBw+aFEC3aNGCDz/8EC8vL2bOnImlpSVt2rRh7NixeHl5MWPGDH799VcePHgAqAHzW2+9xejRo2nbti1r1qzB2tqar7/+GoA1a9YwfPhwXn/9dXx8fFixYgVOBZpGYmJi2LZtG9u3b6dLly54e3vzxRdfEB0dzalTp4qULy4uDmdnZ/r06UPLli156qmnigxmcHFx4ebNm2WeqyRJNUebkkLyoUNkXLmCxsGhposjPSYemZq3w4cPVzh9o0aNyr1/eQQHBzN16lSjdfk1gMUZO3Ys/v7+hIeH079/f4YPH07nzp1LTB8ZGYmnp6dRMNGzZ08SExNJSkoyrO/UqVOp5bx06RJvv/12iduXL1/Ohg0biIuLIycnh8zMTOLi4krNM192djYAVsV8Mx09ejT+/v7cvHmTBQsWEBwczDfffFMk3aZNm2jWrBl9+vQp83jt27c3/K7RaGjYsKHRuvzawoSEBOrXr09kZCQzZswwbLewsKBr165ERkYCcPXqVQICAozyLPieXL58mdzc3CJBcWZmJjdu3ChSm+jv78/s2bNxd3fnhRdeYNCgQQwYMMCoucXa2prMzMwyz1WSpJqTuHs3qSdPYuHkJJtLpWrzyARvtZmjoyMeHh5G6+zs7EpM361bN6Kjo9m3bx/79++nZ8+eLFq0qMSAT1EUkz40HuZ5sJs3b2bevHl89tlndOzYEVtbW1566SVyc3NN2t9J/8GWkpJSZJuDgwMODg54enri4+ODq6srFy9epEOHDkbpNm7cyNixYzEzYeSWRaEHPQshjNblX6+C/dZKu4ZlXeO0tDRsbGyK9KMDin3+qL29Pb/99hsHDhxg//79jB8/nieffJLdu3cb0iQnJ9OoUaMSjylJUs3Ky8wk7fx5rJo1QxTqGiFJVUk2m9ZSTk5OjB07ls2bN7Nw4ULWr18P/C8oycvLM6T18fHh6tWrJCUlGdaFhYXRuHFjo9q4srRv354jR44Uu+3UqVM899xzjBs3Dl9fX5ydncvVpGdlZYW3tzcRERGlpssPpgoHaKdOneLq1auMGzfO5GOWh7e3N2FhYYbXWq2Wc+fO4ePjY9hecACDTqcz9JcDddRsRkYGGRkZeHh4GC0l9c+ztLRkwIABrFy5kj179rBnzx6jaVIiIiLw9fWt7FOVJKmSpF++rH6xk4GbVM1k8FYLLV++nNDQUK5du8bFixf54Ycf8Pb2BqBVq1YAfPfddyQmJpKdnU2/fv1o3bo148eP59KlS+zfv5958+YxZcqUch139uzZrF27lk8//ZTr169z+vRpvvrqK0Ad/Xjy5EmOHTvG5cuXCQwMLDLasiz+/v6cOHHC8DoxMZGFCxdy9uxZYmNjOXr0KGPGjKFTp06GoCnfhg0b6N69u+E6VLbJkyezatUqtm7dypUrV5g4cSIZGRmGfmhvvPEGO3bsYO3atURGRjJ16lTu3btn2N/Hx4dhw4bxt7/9jR9++IHo6GiOHj3K22+/bZQu3969e1m1ahUXL17kxo0bbN26lUaNGtGwYUNAbWY+d+6cYXSxJEk1L+3CBdL1X0AVrZbkH3/EzNq6hkslPY5k8FYL2dra8sEHH+Dr64ufnx8NGzZkzZo1gDrlx7vvvktQUBCNGzdmy5YtmJmZsWvXLtLT0+natSsBAQGMHTvWqA+XKfr168f69ev517/+Rbt27Rg+fDh37twB1JGUffr0YcCAAfj7+9OrV69y1woFBQWxc+dOQ1OrhYUFZ8+eZeDAgXh5eREQEIC3tzf79u0zGp2ZlZXF9u3bq3Sk76uvvsq8efOYPn06vr6+hIeH89133xlqzZ599lmWLVvGnDlz6Nq1KwAvvfSSUR6bN2/mhRdeYPz48fj4+BAYGEhubm6xzdUNGjRg27ZtPP300/z1r3/lzJkz7N2713De+/fvp2XLloZjSZJUsx6cO8cfq1aRevQoAEk//EB2fDzmZYx8l6SqIBRFqekyFKtLly7K2bNni6yPiIigbdu2NVAiqTL079+fUaNGVVnz56Oib9++ZT5OS/4t1C35NeErVqyo4ZJIFRG7eDFa/aTZbu+/z405c7BwdFSfplCMBfv2ATCv0NyW0qNJm5rKUytW3P4zJ6d5dRxP1rxJ1WrFihVyRFYZ0tLS6Nu3L6NGjarpokiShDowIefPP9HY26PLzib15ElQlBIDN0mqavLOk6pV27ZtZW1RGezs7Iye/ypJUs1QtFpy7txBl5WFMDNTv3gqCqnHjyMKjWiXpOokgzdJkiRJKiDnzh3MrK1JPXGC5IMHcejZk/wuRkKjIffuXSybNavhUkqPMxm8SZIkSZKeoijEr19PbmIiSm4uik7H/VOnDKNKNfb2KHl5cnoQqUbJPm+SJEmSpJcTH0/OnTto7OzQODhg2bQp2pQUNLa2AJhZWmJpwvOcJakqyZo3SZIkSdJ7cO4cmJkZPWTe2s2t5gokScWQNW+SJEmSBGTfvk3K4cNy7jap1pPBmyRJkvTY0z54wB+ffYbQaIxq3SSpNpLBmyRJkvRYKjhJfdq5c+Slp2Pu4FCDJZIk08jgrZZwc3Pjyy+/rNZjHj58GCEEWq22UvITQnDgwIFS0yxdupQJEyZUyvFqmp+fH3PmzKn09LGxsTRv3pwHDx48TPEkSSqFotMRt2QJ16dNI/ngQZIPHUJTv35NF0uSTCKDtypW0j/sMWPGGD2r88yZM4wePbrM/A4cOFBpTyh46qmniI+Px7yaZgl/8OABS5YsKXYCWq1WS9euXRFCcP36dcP6wMBAhBBFljfffLNaylwTWrVqRd++feVjlCSpCmXFxpIdH4+mfn0SvvkGbVISmmKeQyxJtZEM3mqJxo0bY2NjU23Hy87OxtLSEmdn52o7ZkhICF5eXnh6ehbZ9sEHH9CwYcMi6z/99FPi4+MNy+nTpwEYNmxYlZe3Jo0dO5a1a9ei0+lquiiS9Eh6cOYMCKFO/dGkCeaOjjVdJEkymQzeaomCzaZZWVkEBwfTpEkTbGxs8PHxYefOncTExODv7w9gqIHasGEDANeuXaNfv37Y2NjQpEkT3nnnHaPmUDc3N5YsWcKwYcOoV68eq1evLrbZdOvWrXTo0AErKytatGjB4sWLATXYCwgIoEWLFtja2vLEE0+U2URa2I4dOxhYzEOaz5w5Q0hICB9//HGRbQ4ODjg7OxuWH3/8kRYtWvDcc8+VeBw/Pz/eeecdgoODqV+/Pm5ubuzbt49bt27Rt29fbG1t6dmzJ7GxsYZ9tFotM2bMMFxzf39/rl27ZtiuKAqzZs3C0dGRxo0bs2TJkiLHTUhIYPTo0TRo0IBGjRoxevRo7t27V2wZ8/Nr3rw51tbWtGnThi+++MKw/dlnnyUpKYlTp06VeJ6SJJVfzt273A0N5f6pU5g3aACAMDeXgxSkOuWRmOdtypQphIeHV8uxOnbsWOXNWStXruTcuXN8//33ODk5ceXKFUMwtX37dl555RXi4+MBNbjJy8tjyJAhuLu7c/r0aeLi4hg3bhyOjo689957hnyXLl3KRx99xLJly7C0tDRqngT48ccfCQgI4KOPPmLQoEEkJiby+++/A2pw4+XlxbRp07Czs+M///kPQ4cOJSoqiqYmTFipKAphYWFMmzbNaH1mZiYBAQF88cUX1Dehv8mmTZsYO3YsZmalf+9Yu3YtCxcu5Pz586xYsYKAgAC6dOnC5MmT+fzzzxk/fjxTp07lv//9r+HabNy4ka+++opWrVoxc+ZMBg8ezKVLl9BoNGzYsIHPP/+c9evX85e//IV58+Zx7tw5evXqZTjmiBEjcHV15dixYwghePfddxkzZgzff/99kfKFhoYSEhLC9u3bcXV15caNG6Smphq2azQafH19OXHiBD179izzukiSVDYlL4+yaZQWAAAgAElEQVQ7mzaRFR2Npn59zOTzSaU66pEI3mq7pUuXFgn4srOzS+zjFhcXR6dOnejcuTMArVu3Nmxz1FftF2zu3L9/P9HR0Rw/fhwnJyc6dOjAggULmDdvnlHwNnjwYIKDgw2vCwdvixcv5rXXXjMEWJ6envTo0QMAW1tbo7578+fPZ+vWrfzwww8EBASUeQ2Sk5NJS0vDxcXFaP2MGTPo06cPzz33HDExMaXmcfLkSa5evcq4cePKPF6vXr2YPHkyAHPnzmX16tX07duXQYMGATBp0iTeeustQ/qVK1cyb948Q83gxo0badGiBfv372fgwIGsWbOGt956i5dffhmAL7/8kmYFnm149OhRIiMjOXjwIBr9Y3P+/e9/07x5c27duoWrq6tR+eLi4vD09OSpp55CCEGrVq2KnIOLiws3b94s81wlSTLN/bAwsmNjsXRxqbS+w5JUEx6J4K22d+wODg5m6tSpRuumT59eYvqxY8fi7+9PeHg4/fv3Z/jw4YZArjiRkZF4enri5ORkWNezZ08SExNJSkoyrO/UqVOp5bx06RJvv/12iduXL1/Ohg0biIuLIycnh8zMTOLi4krNM192djYAVgWaJg4ePMj+/ftNrjXdsGEDPXr0wMvLq8y07du3N/yeXzNYeF1SUhJ5eXmkpaXx559/0r17d8P2hg0b4u3tTWRkJAMHDuTq1atGAy0cHBzw9vY2vL548SIJCQk4FDPNwI0bN4oEb8OHD+eTTz6hbdu2vPDCCwwdOpTevXsbpbG2tiYzM7PMc5UkqWSKopBz5w4WTk7c27sX8wYNZOAm1Xnl7vMmhHAVQnwmhDglhMgQQihCCLdi0lkLIT4WQsQLITL16Z+pjELXNY6Ojnh4eBgtdnZ2Jabv1q0b0dHRTJkyhdjYWHr27MmyZctKTK8oikkfRvUeYiTV5s2bmTdvHtOmTePQoUOEh4fTrl07cnNzTdrfyckJIQQpKSmGdUePHiUqKgoHBwfMzc3x8PAAwMfHh3nz5hntn5WVxfbt241G6JbGokBzSP61KW5dwXmeyrqGpW1PS0vDw8OD8PBwo+XatWt07dq1SHo3NzeuXbvGokWLSEtLY9CgQUUC5+TkZBo1alRqmSRJKl1uQgJxS5dyZ+NG8tLTDQ+Yl6S6rCIDFjyAV4Bk4Fgp6dYBwcD7wItAPPCDEKJjBY752HFycmLs2LFs3ryZhQsXsn79euB/AUheXp4hrY+PD1evXiUpKcmwLiwsjMaNGxvVxpWlffv2HDlypNhtp06d4rnnnmPcuHH4+vri7OxcriY9KysrvL29iYiIMKx78803+e233wyBznfffQfAnj17jJo0Ab799luys7MZOXKkycc0lYODA02bNiUsLMywLjk5matXr+Lj4wOAl5eXYaQrwP3797l69arhta+vLzdv3sTe3r5IoF7SKGJbW1tGjBjBv//9b7788kvWrVtntD0iIgJfX9/KPFVJeuxk3riBLjubtPBwOQGv9MioSLPpUUVRmgIIIf4O9CucQAjhC4wCxiuK8pV+3RHgMrAQGFzhEj8Gli9fjqurKx07diQrK4sffvjB0ESX3zfqu+++o0ePHtSvX59+/frRunVrxo8fz6JFi7h16xbz5s1jypQp5Tru7NmzGTRoEB4eHgwcOJCkpCQuX75MUFAQ7u7ubN++nWPHjuHk5MS8efPKPY2Fv78/J06cYPz48QA0adKEJk2aGLbn10Z6enoarQe1yXTo0KHFNktWhsmTJ7NgwQLc3Nxo2bIls2fPpmXLlvTv3x+AN954g8mTJ9O5c2f+8pe/MH/+fKOauH79+tGhQweGDRvGhx9+SPPmzYmKiiI0NJS1a9cWOd7GjRtRFIUnn3wSjUbDzp07jZphb9++TVxcHM8++2yVnK8kPS7SLlxAY2trGFkqSY+Ccte8KYpiyn/swUAusK3AflpgK9BfCCHHZJfC1taWDz74AF9fX/z8/GjYsCFr1qwB1Oa2d999l6CgIBo3bsyWLVswMzNj165dpKen07VrVwICAhg7diwzZswo13H79evH+vXr+de//kW7du0YPnw4d+7cAeD111+nT58+DBgwAH9/f3r16lXuWqGgoCB27txpclNrvj/++IMDBw6Y3GRaEe+88w4BAQGMGzeOLl26kJaWxu7duw2DDwIDA3nrrbf4+9//ztNPP42vr69RP0QzMzP279+Pt7c3w4YN4y9/+Qtvv/02DUr4h+Hg4MDq1avp1q0b3bp1Iykpia1btxq2h4aG8vzzzxcZ4CFJUtmyb91C0WpRtFoyIyPRlNJNRZLqIlGwz0+5d1Zr3v4NtFYUJabA+q1AJ0VRvAulfwU1oGuvKMrl0vLu0qWLcvbs2SLrIyIiaNu2bYXLLNWs/v37M2rUKJNGjD6udDod7dq148svvzSaiqQw+bdQt+TXhNf2AVZ1XWZUFLc+/ZSmY8di3qABf3z+OZaFavKrw4J9+wCYV8zcltKjR5uaylMrVtz+MyeneXUcr6om6XVC7RNXWFKB7UUIIV4TQpwVQpxNSEiooqJJNWnFihVypFcZbt++zdtvv11q4CZJUlG67Gziv/wSFIXkAwdI+v57hJzLTXoEVdVUIQIorkqv1P/aiqKsBdaCWvNWBeWSaljbtm1lbVEZXF1diwzYkCSpbJnXrpGXno5F06bkxMeDomBhwiTiklTXVFXwlgS0LGa9Y4HtkiRJklRp7p8+jbCwQAiBWb16hscIStKjpqqaTS8DrYUQhScWawfkANeL7iJJkiRJFaPLyiL9t98M04GY168vBypIj6yqCt52AxbAy/krhBDmwEjgR0VRsqvouJIkSdJjRNHp0N6/zx399DtCP0Jckh5lFWo2FUKM0P+aP1fCC0KIBCBBUZQjiqKECyG2ASuEEBZANPAG0Boo/oGekiRJklRO8evWkX7xIgiBRePGNV0cSaoWFe3zFlro9Wr9zyOAn/73IGAxsAhoAFwAnlcU5XwFjylJkiRJBrqsLNIvXcKicWOEWVU1JElS7VOh4E1RlDJ7gCqKkglM0y+SJEmSVKmyYmIAZOAmPXbkHS+V29KlS5kwYUJNF6NSzJ8/v1zzqZmaXqfT0bZtW6PnpUqSVLnSIyJAjiaVHkMyeKtifn5+huHqtra2dOzYkW3bthVJl5CQwKRJk3Bzc8PKyoo2bdoQGBjIlStXAPXZnvn5mJub07p1a2bPnl3so6Z27tyJEIIxY8aUWrb853MWXIYOHVrqPg8ePGDJkiXMnDnTsG7KlCm4u7tjY2ODs7MzAQEB3L1712i/r7/+mvbt22NjY0Pbtm355ptvSj1OXWdmZsaMGTOYPXt2TRdFkh4piqKoi05H2q+/Yl6/fk0XSZKqnQzeqsGUKVOIj4/n0qVLjBo1ilGjRhEeHm7YnpCQwJNPPsm5c+dYu3Yt165dY8eOHbi4uLBgwQJDOhcXF+Lj44mJieGjjz5i+fLlLFmyxOhYd+/eZerUqSbXJnXr1o34+HjDsmHDhlLTh4SE4OXlhaenp2Gdr68vX331FREREezatYvo6GhGj/7fuJS9e/cSHBzMzJkzuXz5MtOnT2f06NH88ssvJpWxrnr55Zc5deoUkZGRNV0USXpkpBw6xJ2NG7l/+jTapCSElXxUtvT4kcFbNbC1tcXZ2ZnWrVszY8YMHB0dOXLkiGH7e++9h1ar5eeff6Zfv360bNmSJ554gg8//NDwQHpQa3OcnZ1xdXVl5MiRjBkzhr179xodKzg4mGnTpuHu7m5S2SwsLHB2djYsJT1IPd+OHTsYWOhZfUFBQTzzzDO4ubnx5JNP8s4773D8+HHD9pCQEF599VXGjBlDmzZtmDBhAoMGDSr1GY/5zZOffvopLi4uNGjQgMWLF5Odnc3rr7+Ovb09Hh4e/PTTT0b7bdq0CQ8PD6ysrOjQoQPff/+90fZt27bRqlUrbG1tCQgIIDMz02h7Xl4ec+fOxdXVlfr16+Pn58dvv/1WYjm3bNmCj48P1tbWODs789prrxm22dnZ8fTTT7Njx46SL6gkSSbTZWWR9P33PDh9mrubN2Pu6Cgn4ZUeS1X1hIVqdTIO7mVUz7Ea1oOnWlRsX51Ox86dO0lKSsLS0tKwLjQ0lFmzZmFtbV1kn9KCKRsbG3Jycgyv169fz/3795k4cSJBQUEmlenChQs4Oztjb29Pv379WLx4MQ76SS4LUxSFsLAwpk0reQxKSkoKISEhRjV/2dnZRc7D2tqao0ePllq23377DU9PTw4dOkRYWBhBQUEcOXKEgQMHcu7cOVasWEFAQACxsbFYWlpy8uRJxo8fzz//+U/69+/P119/zdChQ4mMjMTNzY1r164xevRo5s2bxyuvvMKOHTv46KOP8PX1NRxzwYIFfPfdd2zZsgUXFxe++uor/P39uXbtGvb29kblu337NkFBQWzcuJHu3buTkJDAuXPnjNJ07tyZEydOlHqekiSZ5n5YGLrsbCydndFlZmJWzGemJD0OZM1bNVi6dCl2dnZYWVkxfPhwmjdvzogR6lR5CQkJpKam4u3tXa48f/31V0JCQujduzcAMTExzJkzh/Xr15v8TbR79+5s2rSJn376iU8++YTDhw8zZMgQFKX4x8omJyeTlpaGi4tLkW2rV6/Gzs4OR0dHYmNjjfr19e3bly1bthAWFoaiKBw5coRvv/2WO3fulFo+c3Nz/vWvf+Hj40NgYCDe3t5YWFgwefJkPD09mTt3Lnfu3DE0S65cuZJhw4YxadIkvL29WbRoER07dmTVqlUA/Pvf/6ZLly7MnTsXb29vZs+ebRS4ZWVlsWzZMjZu3MjTTz+Nh4eHIZjdvXt3kfLFx8djZWXFwIEDadWqFV26dOH//u//jNK4uLhw8+bNUs9TkqSy5SYmkrhnD+YNGiA0Gvn0BOmx9kjUvFW0Jqy6BAcHM3XqVG7fvs20adNYsGABjSswmeTt27exs7MjLy+PnJwchg0bZugTN27cON5//31at25tcn7PP/+84fcOHTrQrl07PDw8OHfuHF26dCmSPjtbfTCGVTF9TEaPHo2/vz83b95kwYIFBAcHGwYl/N///R8RERE888wz6HQ6WrVqxdixY8vsX+fp6Wl0rKZNm9K+fXuj16AGwACRkZGMHTvWKI+nnnrKENxdvXqVbt26GW3v1q0bp0+fBiAqKorMzEy6d+9ulCYzM5MbN24UKZ+vry++vr60adOGAQMGMGDAAIYOHWqoVQW1hrFw06wkSeWjKAp3vv4aADPZx02SHo3grbZzdHTEw8MDDw8PQkJC6NmzJxcvXsTZ2ZnGjRvj4OBgUqf2pk2bcuzYMczNzWnWrJlRkHD06FFOnDjBxIkTAbU5FmDr1q1kZ2ejMeGRMe7u7jRo0IDo6OhigzcnJyeEEKSkpBTZ5uDggIODA56envj4+ODq6srFixfp0KEDZmZmrFy5kk8++YS7d+/i4uLCe++9V2agaWFhYfRaCGG0Lr+GMf9cFUUptdaxrO1paWkAHD58uEgzr5OTU5H05ubmHDp0iKNHj7J//35mzJjB0qVLOXnypOG9SU5OplGjRqWdpiRJZciKiSErOhqLJk1quiiSVCvIZtNq5uXlhZ+fH4sWLQLUQQjDhw9n9erVZGVlFUmfmppq+F2j0eDh4YGbm5tR4AZw8eJFwsPDDcvgwYMZOHAg4eHhJgVuADdv3iQlJQU3N7dit1tZWeHt7U1ERESp+eQHU2aFJs60sLCgefPmKIrCt99+y6BBg0wql6l8fHyKzKsWFhaGj48PAN7e3oZatnxnzpwx/N62bVssLS25ffu2IdjOX4oL3kB9T5599lmWLFnC6dOnOXfunNFI4oiICKOmWUmSyi/12DGEubkcnCDVKF12NpnXr5MeEUH6pUtkXLuGtpjKjOogg7caMHHiRNatW0d8fDwA//jHPxBC0KNHD3bv3k1UVBTh4eHMnTuX119/3aQ827dvb7Q0aNAABwcHo2bGPn368Pnnnxtez5gxg+PHjxMTE8OhQ4cYNmwYPXr0oHPnzsUdAgB/f3+jDviJiYksXLiQs2fPEhsby9GjRxkzZgydOnUyBE137txh3bp1XL16lbCwMIYMGUJmZiYzZswo13Ury6RJk/jvf//LqlWruHr1KvPnz+f8+fO8+eabALz22mucOXOGxYsXc/XqVT788EOjQMve3p6JEyfy5ptv8s033xAdHc2pU6d47733uHz5cpHj/fLLLyxZsoTz588TGxvLpk2bsLKyolWrVoY0J0+exN/fv1LPU5IeJ7lJSTw4dw7zMkbCS1Jl0WVlkXH1Knn6Li+67GzSf/+dpB9+IP3iRTKvXCEzKoqMS5dIOXSIlKNHUbTaai2jDN5qQO/evfHy8mLZsmWA2hx6+vRpevbsydtvv027du146aWXiIuLM5rn7WFFRUWRmJhoeB0bG8vLL7+Ml5cXQUFBdO7cmV27dhWpMSsoKCiInTt3GiYHtrCw4OzZswwcOBAvLy8CAgLw9vZm3759hho/RVH44osv6NSpE/3798fa2ppjx47RsGHDSjs3UPu3rVu3juXLl9O+fXu++eYbdu7caahJ9PDw4D//+Y+hLJcuXeKNN94wyuPjjz/mzTffZPr06Xh7e/PKK68QFxdXbFnt7e0N07u0bduWLVu28N///tfQF+/ChQskJiYyePDgSj1PSXrUJR88SM6dO+iysohftw6h0SBMbEGQpIpQFAWdvl93XloaGZcvo0tPByDn9m0yIyOxdHamQZ8+NBw6lEZDhuD0wgvUa9cO7b175CYlVWt5RUkjC2taly5dlLNnzxZZHxERQdu2bWugRFK+/v37M2rUKMaNG1fTRanV/v73v+Pm5sacOXOqJH/5t1C3TJkyBaDU+Q0lyMvM5MaMGWBmZnhmqXnDhnWuyXTBvn0AzCs0L6ZUO2XdvEn6hQs49O6Nub09efqpaIQQ5KWno+TlYV5ouigARacjJz4eYWlJr9Wrb/+Zk9O8OsorByxI5bZixQqjvmJSUTqdDnd3dyZNmlTTRZGkOiXn9m3QaLBwcgIh5EPnpXJRdDoyrlwh988/cejd2+j+yUtPJys2Fo2dHUpODhoHB/Lu3yczKkp9Xb8+Gv3j1jQ2Nob9NLa2JR5PmJlh1bw52gL906uDDN6kcmvbtq2s8SmDmZkZs2bNquliSFKdkxUXBzqdbCaVykXRatGmpJB+6RLa5GQsGjVCmJmhy8kh9+5dLJ2dEZaWZN+6ZWgOzWfu6AiWltj5+laohjcvLY3suLjKOhWTyOBNkiRJqjUyIyPlkxOkEhWe8iknIYGMy5fRJicDIKysqN+1K1auroA6ufODs2dxeuEFzKyscOzbl7z79xFWVuTcuoWSl4eNt/dDNcvnJieTee3aw51YOcngTZIkSapxiqKgS08nMyoKs3r1aro4Ui2h6HRk/P47ls7OmNnYcP+XX6j/xBMIKysenD1LXkoKwsoKGx8fzO3tsWjSBLOC84FqNDTo3Ruhn15LmJkZRi7beHpWShktnZ2x79kTSnkOdmWTwZskSZJU4+6fPMndLVvUjuElPF9ZerQpWi0ZV69i4eSEpbMziqLw4Nw5cm7dIvvWLex79FATmpmh5OSoNWgWFjg8/bRRH7WCLPWj/6uSmYVFtT/5QwZvkiRJUo1SFIXkn39G4+CARta6PVYUnQ5dZiZmNjag0ZAVFUVmZCTWrVujsbcn59YtrN3dsXZzw9zengbPPmto4nTs2xd4PB+ZJoM3SZIkqUZlRUeTm5goH3/1CFMUhcxr18hNSEBTvz71fHwws7QkKyqK9EuXcPDzw8LREacXXiD999/JiooCwKJJE2w7dDAEbAX7pj2OQVs+GbxJkiRJNep+WJg6EW8dm8tNKl1eejo58fGY2dqSdf06uYmJaOrXJzchAW1KCg5PP421hwfC3ByNnR0Awtwc2w4dQKcj+9atCo8AfdTJ4E2SJEmqVpnXr/MgPJwmI0agaLWk/forGtnPrc7TZWcjLC0RQpBx9SoZERGgf9a1sLDArlMnrFq1QnvvHum//66uFwLr1q2N8hFCYNexo1rjJqeMKZac/VAqt6VLlzJhwoSaLoZJhBAcOHCg0tOHhITwwgsvPEzRJOmxlXbhAik//4w2JYWsmzfR5eQYjRCU6h5FqyXl8GEyr1wB1Bo0K1dXGvTpg/1TT+HYrx/Wbm4IIbBo1IgGzzxTZo2aDNxKJoO3Kubn54cQAiEEtra2dOzYkW3bthVJl5CQwKRJk3Bzc8PKyoo2bdoQGBjIFf0fwoYNGwz5mJub07p1a2bPnm14xmhBO3fuRAjBmDFjSi3b/PnzDXnmL0OHDi11nwcPHrBkyRJmzpxpWHf79m1effVVmjRpgr29PWPGjCElJcWUy1NnjRw5ksjISI4dO1bTRZGkOif90iWUvDxST50i7eLFmi6OVAG67GzSfvuNjKtX1bnXzM2xdHbG0tkZAJs2bajfuTPm9vZYNm2KmX6qDqlyyGbTajBlyhTeffddMjMzCQ0NZdSoUXh7e9OxY0dADdyefPJJXFxcWLt2LT4+PiQmJhIaGsqCBQvYsmULAC4uLpw/fx6tVsuJEycICgrCxsbG6NmZd+/eZerUqfTq1cuksnXr1o1du3YZXluXMTlmSEgIXl5eeOrnx9HpdAwdOhQnJyd+/vlntFotU6dOZcyYMezdu7dc16ku0Wg0jBo1ijVr1vD000/XdHEkqdZTdDqEmRna1FRy793D0tmZe3v2oOh0WMqBCrWeLieH9N9+Q5uaipmlJdrUVBStFszMsGrWDI2dHXa+vjVdzMeGrHmrBra2tjg7O9O6dWtmzJiBo6MjR44cMWx/77330Gq1/Pzzz/Tr14+WLVvyxBNP8OGHH7JmzRpDOjMzM5ydnXF1dWXkyJHFBkjBwcFMmzYNd3d3k8pmYWGBs7OzYWmgn7ywJDt27GBggQctX79+nTNnzvDFF1/QoUMHOnXqxKpVq9i3bx+RkZEl5iOEYP369Tz33HPY2NjQrVs3oqKiOHz4MO3bt8fe3p6xY8eSlZVl2OfPP/9kxIgR2NnZ4ejoyIQJE0gv8JiTlJQUhg0bho2NDV5eXuzfv7/Icc+fP4+fnx82Nja4ubkxb948tFptsWW8d+8eL7/8Mk5OTtja2uLr68upU6cM2wcOHMjOnTuLrf2UJOl/sm/dIvYf/yAvLY2smzdBCMysrLBq3hzrFi0e61GDtZUuKwtF319Nl5tLyuHDZN+6hZmNDYqiYNGkCQ2efZaGL76IWSnP/pSqxiNT87an5DjBoKUD+Dr/L71XQ/BuBFla+Cmq7P3z01eUTqdj586dJCUlYamvQtbpdISGhjJr1qxia71KC6ZsbGzIyckxvF6/fj33799n4sSJBAUFmVSmCxcu4OzsjL29Pf369WPx4sU4lNBxWFEUwsLCmDZtmmFddna2oSz58s/j5MmTeHt7l3jsRYsW8emnn7Jq1SomTJjAqFGjsLW15auvviInJ4chQ4awdu1aw8Pdx44dS2pqKkeOHCEjI4OgoCCmTp3K2rVrAZg8eTK///47Bw8eRAjB5MmTjY537949/P39effdd/nyyy+5desWwcHB2NjYGDUD55s7dy4PHjzg6NGj2NjYcOHCBcP7BtCpUyeys7P59ddf6datW6nXWZIeV4pWy53//IesGze4t38/2XFxsi9TLaXodOTEx2PRpAmpx45h36sXGhsb7p86hS4jA4devbBo9BD/BKVK88gEb7XZ0qVLWbFiBdnZ2Wi1WlxdXRkxYgSgNpmmpqaWGuQU59dffyUkJMTQry0mJoY5c+Zw4sQJk4dVd+/enU2bNuHh4UFMTAyzZs1iyJAhHDp0qNg8kpOTSUtLw8XFxbDO29sbV1dXZs2axcqVK8nLyzM04965c6fU47/xxhsMGjQIgEmTJvHqq69y+vRpunbtCsCIESM4cuQIkyZN4sqVK/z0009cvnyZdu3aAbBy5UqGDBnCxx9/jBCCkJAQ9uzZQw/9LNyLFy/G39/fcLxVq1bx7LPPMmPGDAA8PDxYsGAB77//frHBW1xcHD179qR9+/YARWozra2tadCgATdv3pTBmySVIP3yZXL++APrVq1IPnAAMwsLLKph1nupbIpOR1ZMDOh0aOzt1XnY7t6lfrduCEtL8lJSQKtFe+8eNj4+MnCrRR6Z4G1Q+WIfo/TW5uXfvzyCg4OZOnUqt2/fZtq0aSxYsIDGjRuXO5/bt29jZ2dHXl4eOTk5DBs2jAULFgAwbtw43n//fVoXGnJdmueff97we4cOHWjXrh0eHh6cO3eOLl26FEmfX8tmVaCJw9LSktDQUAIDA7G3t0ej0fDmm2/StGlTzMxKb5XPD4oAmuo/zAuv+10/nDwyMpL69esbAjeAnj17otVqiYqKwszMDK1WaxREFQ6oLl68yO7du7HTzycEkJeXR25uLjqdrkh5g4ODGTlyJD/++CP+/v6MHDmySJBtbW1NZmZmqecpSY+ztAsXEBYWaof2pk3V3+W8XTVO0Wq5HxZGbkLC/1YKQb327bFq3hzLZs0Mqx38/AzPA5Vqh0cmeKvNHB0d8fDwwMPDg5CQEHr27MnFixdxdnamcePGODg4lNo/LF/Tpk05duwY5ubmNGvWzKgJ7+jRo5w4cYKJEycCanMswNatW8nOzkZjQjOFu7s7DRo0IDo6utjgzcnJCSFEkZGk3bt358qVKyQmJmJhYYGFhQWfffZZmYGkRcGHB+s/zAuvyz8PRVGKfOAXfK0oSpF1hdOnpaXxt7/9jffff79IWYoLNAcPHsyNGzfYs2cP3333HYsXL2bTpk2MHDnSkCYlJYVG8tuoJBWhvX8fTb16pF+6hMbeHkCOOKwl0i9eJDc5GW1SEnadOmHRqBF5mZmY16+Pmb7bS8HPTwtHxx+qQqkAACAASURBVJoqqlQCOWChmnl5eeHn58eiRYsANWgYPnw4q1evNuqcny81NdXwu0ajwcPDAzc3N6PADdRapfDwcMMyePBgBg4cSHh4uEmBG8DNmzdJSUnBzc2t2O1WVlZ4e3sTERFR7PZGjRrh4ODAN998g6WlJX31z52rDD4+Pty/f9/o2KdOncLc3Bx3d3fc3d0xNzfn9OnThu0Ffwfw9fXl8uXLuLu7G4Lp/KUkLi4uvPbaa+zcuZMJEyawceNGw7bY2FgyMjLwlSOsJMmIotNx69NPiV+3Dl12tpzDrYbkZWSQm5hotE6Xm0vOnTto793DrmNHrN3c0NjZYdm4sSFwk2o/GbzVgIkTJ7Ju3Tri4+MB+Mc//oEQgh49erB7926ioqIIDw9n7ty5vP766ybl2b59e6OlQYMGODg4GDVD9unTh88//9zwesaMGRw/fpyYmBgOHTrEsGHD6NGjB507dy7xOP7+/pw4ccJo3bZt2zh+/DjXr19n3bp1vPHGGyxcuBAnJ6fyXJZS+fj40K9fPyZMmMD58+c5ceIEkyZNIigoCAcHB+zt7Rk1ahRTp07ll19+4ZdffmH27NlGebz11lvcuHGD1157jQsXLhAZGcn27dsNgXRh8+bNY+/evdy4cYOzZ89y/Phxo2bTkydP0q5dO5oVaF6QJAmyYmPJvXuXtN9+q+miPLbyMjJIPXKErOhoAHLu3CH94kXMLCxo0KcPjn37Yl3CF3Wp9pPBWw3o3bs3Xl5eLFu2DFCbQ0+fPk3Pnj15++23adeuHS+99BJxcXGGPm2VISoqisQC38JiY2N5+eWX8fLyIigoiM6dO7Nr165S+6oFBQUVmR7j5s2bjBw5knbt2vHxxx+zbNky3nnnnUord75Nmzbh4uLCM888w8CBA+nVqxfLly83bF+xYgXe3t707t2bUaNGFWkebdGiBUePHjUMROjatSvLli2jZcuWxR7P3Nyc6dOn065dOwYOHMiTTz5pFOht377d5FG9kvQ4eXDmDGg0WLq4YFGB/r1S+SlaLdqUFHS5uaRfvkzKwYMoeXnYeHqiKArZcXHq3Gz6+fY09evXdJGlhyDy+wr9P3v3HV9lffd//HWdfZKTkz1ZiUQ2BQTZKg5QW1dRbu9SqVLFH22tRe+7WtQ+KM66btFHax2l1WppLW0Vd60WEVCGyN4rzJBB1jnJ2df398eVnCQkYeZkfp6Px3mQM67r+70OGe/znR3NqFGj1Ndff93k8e3btzNw4MB2qJGoc+WVVzJ9+nRuvfXW9q5KuykoKGDs2LHs3LmzxaVVYk1+FjqXOXPmAMaHjK5MhcPsmzsXU3x8t+4unf/BBwDMa7AuZiwopUApKleuJFxaCiYT6Dq2nBziBgzAUvv7KeLzYXI4ZLJIjIQrKxm/YMHRomCwR1uUJxMWxBlbsGABa9eube9qtKvDhw+zcOHCdgtuQnRUvn370INBmZ0YQ6HyckIlJTjy8lB+PxXLlqFCIRx9+6LCYRx5eU0mGZgbrMUpOr+YhjdN0yYA84DhgAPYA/xGKfWHWJYrYmvgwIHdvsXndLcfE6K7qd60CaR1p9UoXSfi9WKpnbELEC4tpWbrVpznnYeuFNb0dCxJScSd4XqhovOKWXjTNO1bwKfAKmAWUAPcBCzUNM2ulPrdyY4XQgjROeh+P5rdDkrhWbeuUdAQ56Zm2zZ8u3fjGj4cW3Y2JocDR34+tp490SwWLG437jFj2ruaoo3FsuXtvwEzcK1Sylv72L81TRsG/ACQ8CaEEJ2cHghw8OmnSb70UmO9sOpqbLKDQqtQ4bCxA4LJhHfDBrStW3FPmIA1OVm6Qbu5WIY3GxACTlx+vgKQFf+EEKILKF+6lOCRI5S8/TYmux2zbFLeagKHD6NCIdzjxxOprsZkt0urpgBiu1TIa7X/vqBpWo6maUmaps0CLgeea/kwIYQQnYHu91P+r39hzcoy7gcCmBtsPyfOTriqisCRI1Rv3ow5KQlrRgbO887D3qMH2mkuui66tpi1vCmltmiaNgl4G/hx7cMhYLZS6q/NHaNp2p3AnUCLa28JIYToGHx796J0HZPViqkVF+Xu7gIHDuDbswezy4V7zBhZ3kM0EcsJC+cD/wC2ArMxuk+vB17SNM2vlPrziccopV4BXgFjnbdY1U0IIcTZq96xA93rpWbXLrSTLOotWqaUQq+uxhQfj6ZpeL75BmtqKo4+fXDk52PNyMCSktKt18oTLYvlmLfHMVrarlFK1S3H/5mmaanA85qm/UUppcewfCGEEK3Is349zrw8it54g4jXi8lmi246L06fHgrhWbuWcFkZKbWL+KpgkHBZGap3b8xOp0xIECcVy49MQ4GNDYJbnTVAKpARw7JFDD311FPcfvvt7V2NVpGbm8vvf//7Vn/9l19+yZAhQ9B1+XwiuoZgURGFv/89B554gojHg8lmMwbR22ztXbVOIVRaSsWyZVStXk35v/9NqLgYZ79+0S5R99ixuEaMkC5ScVpiGd6OAcM1TTvxJ3sM4AfKYlh2hzFp0iQ0TUPTNOLj4xk+fDhvvfVWk9eVlJRw9913k5ubi91u57zzzuO2225jx44dALz22mvR81gsFvLy8njwwQcb7TG6cOFCBg4ciNPpJC8vj0ceeYSTbX/2q1/9KnrOutsNN9xw0uvxeDw8+eST/OIXv2jyXDgc5sILL0TTNPbs2dPouQ8//JBBgwbhcDgYNWoUa9asOWk5nd348eNJTU3lz39uMjpAiE7J8/XX0cHy1tRULElJ2Hu0yU5AnZ7/4EEqly9H9/kIl5VhdrlIvPhi4vr1a++qiU4qluHtN0Ae8J6maddrmjZF07TfAN8DfqeUCsaw7A5lzpw5FBYWsmXLFqZPn8706dPZsGFD9PmSkhLGjBnDunXreOWVV9i9ezd///vfyc7ObrQxfXZ2NoWFhRQUFPDrX/+a5557jieffBKAL774gtmzZ3P//fezfft2nn/+eZ5++mkWLlx40rqNHj2awsLC6O2111476esXLVpEv379OP/885s898gjj5Camtrk8V27djF16lSmT5/ON998w/jx47n66qspK+va+X3GjBn87neynKHo/FQ4TOWKFViSkrC43WgW2VnxdKhAgPKlS/F+8w3W9HSSr7iClKuvJunii7HKBA9xDmIW3pRSfwe+DdiB32NMXpgI/AT4eazK7Yji4+PJysoiLy+P++67j+TkZJYtWxZ9/oEHHiAcDvPZZ58xZcoUevfuzQUXXMATTzzR6I+/yWQiKyuLnj17cvPNN3PLLbfw/vvvA7B69WoGDx7MbbfdRm5uLtdddx2TJ0/m66+/PmndrFYrWVlZ0VvSKfYj/Pvf/853mtloee3atSxatIinn366yXMvv/wyI0eO5KGHHmLQoEE8//zzuFwu3njjjRbLyc3N5emnn2bq1KnExcUxaNAg1qxZw+bNmxkzZgwul4vvfOc7jQJgdXU1d9xxB8nJybhcLm688UaKioqizwcCAe644w5cLhe9evXi9ddfb1Luvn37uPbaa3G5XOTk5HDXXXdRU1PTbB39fj+zZs0iIyMDp9PJgAEDeOedd6LPX3PNNXz11VccPny4xesUoqNTSlG6ZIkxxs1ub+/qdHihsjJqdu8GQLPZMDudWJKSSBg9WkKvaDUx/U5SSn0EfBTLMgAOPfssNTt3xroYAOL696fX//zPWR2r6zrvvPMOZWVl2GrHiei6zuLFi5k7dy4Oh6PJMScLU06nk2DQaMAcO3Ys8+bNY8WKFUycOJFt27bx5Zdf8tvf/vakddq4cSNZWVm43W6mTJnCY4891uJm60opVq1axb333tvocZ/Pxw9+8ANefvllEhISmhy3evVqLrvssuh9TdO47LLLWL169Unr9swzz/Dss8/y61//mvvvv58ZM2aQk5PDk08+SVJSEtOmTePRRx/l//7v/wC45557WLZsGUuWLMHlcjF79mxmzJjBJ598AsDjjz/OBx98wD//+U8yMzOZM2cOxcXF0fKCwSBXXnkl1113Hc888wxVVVXcdddd3Hvvvbz00ktN6vfCCy+wbt06PvroI1JSUtixYwf2Bn/csrKyyM7O5ssvv+S//uu/TnqtQnRE4aoqSpcsoWrVKtk14TT5CwoIlZSAUqBpuMeObe8qiS5IPga0gaeeeooFCxYQCAQIh8P07NmTm266CTC6TCsrK+l/hhsKr1+/nkWLFnHLLbcAcNFFF/HSSy9xxRVXEIlEiEQiPPbYY0ydOrXFc4wdO5Y//elP5OfnU1BQwNy5c7n++utZunRps4Nmy8vL8Xq9ZGdnN3r8vvvu4/LLL+eyyy6joKCgyXHFxcVkZDSen5Kens66detOeo033XRT9Pruv/9+xo0bx6OPPsqkSZMAuP3221m8eDFgjMX74x//yJIlS7j44osBY5zg4MGD2bp1K4MHD+Z3v/sdDz/8MFOmTAHg1VdfbdT9+9e//pWkpCSeffbZ6GPPPfccl156Kb/97W8xn7A45qFDhxgxYgQjR44EIC8vr8k1ZGdnc/DgwZNepxAdVfGiRVRv2YItK0uWBGmGUoqI14teU4NmtWJNSSF+yBDQNKj90ChELHSJ8Ha2LWFtZdasWdxzzz0cPXqUe++9l/nz55Oenn7G5zl69Cgul4tIJEIwGGTq1KnRMXFbt27lwQcfZMGCBUyYMIHNmzdz9913k5+fz7Rp05o931VXXRX9eujQoQwaNIj8/HzWrVvHqFGjmrw+EAgANGpd+s9//sPHH3/caAzfiU42aeJkhgwZEv06s/ZT/4mPlZSUAEZ3ZzgcZmyDT7mDBg0iKSmJnTt30qNHD0pKShg9enT0+fz8fJKT63dq27x5Mxs2bMDVYIV4pRTBYJAjR440WTh6xowZTJ48mQ0bNnDllVdy4403RoNcHYfDgc934g5xQnR8EZ+P6h07sGZmSnBrhtJ1qjdvxr9vHwDWzEwSx4+X2beiTXSJ8NbRJScnk5+fT35+PosWLYqGq6ysLNLT00lMTGTnaXT7ZmZmsnz5ciwWCzk5OdGuV4AnnniCyZMnM3v2bMAIY3v27OHpp59uMbydqG/fviQlJbF///5mw1tKSgqaplFRURF97IsvvmDv3r1NuloHDBjAgw8+yPz588nMzGzUPQnNt8adyNpgccq6lsATH6tbiqMuIJ5qmv3Jnvd6vVx00UW88sorTZ47sbURjMke+/fv54MPPuDjjz9mwoQJPProo/zv//5v9DXl5eWkpaWdtE5CdBT+Awfw7dtH8qWX4t+7F5SS4NYMFQ5TtWYNoaIiHH37YsvKwpIsW3aLtiM/lW2sX79+TJo0iUcffRQwJiHceOONvPjii/j9/iavr6ysjH5tNpvJz88nNze3UXADY9zZid16pjP8pXvw4EEqKirIzc1t9nm73U7//v3Zvn179LEf//jHbNq0iQ0bNrBhwwY+/PBDAN577z1+8pOfAEbI+eyzz6LHKKVYunQpY8aMOaP6nUzfvn2xWCysWrUq+tjOnTupqKhgwIABJCUlkZGR0WiJkr1791JeXh69P2zYMHbu3EnPnj2jYbvuZm1hlfOUlBRmzJjBn//8Zx5++GH+8Ic/RJ8LhULs27ePYcOGtdp1ChFL3k2bKF60CH9BAZ5162QfzVrhykpUOAwYv7+qVq8mVFRE/LBhuL71LWwZGbITgmhTEt7awV133cXChQspLCwEjIH0mqYxbtw43n33Xfbu3cuGDRv45S9/GW1JO5Wrr76aN998kzfffJP9+/fz3nvv8cILL3DNNddEX3P55Zfzm9/8Jnr/vvvuY8WKFRQUFLB06VKmTp3KuHHjmnT9NTR58mRWrlwZvZ+RkcGQIUOit3616xadf/750Za1O++8k3Xr1vHYY4+xfft25syZg8fjYcaMGaf/pp1CQkICP/zhD5kzZw4rVqxg/fr1zJw5k8mTJzNo0CAAZs+ezfz58/n000/ZuHEjd955J84Gq5h///vfx2azcfPNN7N27Vr27NnDe++916glraHnnnuOxYsXs3v3bjZv3sy//vWvRmMXv/76a5xOZ7OtmEJ0RDXbt6PZbBx8+mmq1qzB0sLkpe5EKUXVqlV41q4lcPRodIHd+G99C+d557V39UQ3JeGtHVxyySX069ePZ555BjC6Q9esWcOECRP46U9/yqBBg/jud7/LoUOHGq3zdjJ33HEHjzzyCA8//DCDBg3i7rvvZtasWcydOzf6mr1791JaWhq9f+DAAaZNm0a/fv2YOXMmI0eOZMmSJSdtsZs5cybvvPNOo8WBT2XAgAH84x//4M0332T48OEsX76cDz/8sNk14c7Fs88+y8SJE7nmmmu46KKLyM7ObrQcyQMPPMDVV1/N9ddfz7e//W1uueWWRl23CQkJfP7559hsNiZPnsywYcN46KGHmu0yBWMJmEceeYRhw4YxadIkUlNTGy3t8re//S0aCIXo6PRgkMCRI9iysrBlZWHPyemWS1uEKyqo+uorarZvx79/PyocJmHkSJz9+6NZLOjV1ZjdbhzNTFASoq1oZzuYPNZGjRqlmlujbPv27QwcOLAdaiTqXHnllUyfPp1bb721vavSYXk8Hs4//3y++uqrZmehtgb5Wehc5syZA8CCBQvauSbN8xcUcGjBAmxnMZmqK6lYvpxwWRnUjqdNueaaRl2i/gMHsKSkYGlmWaSG5n/wAQDzmlkXU3Q94cpKxi9YcLQoGGyTbUe638cqcc4WLFjA2rVr27saHdqBAwd49tlnYxbchGht/kOHooGlO9JDIULHjhEuLSV+6FAcubnooVCTsWyOPn3aqYZC1JPwJs7YwIEDpcXnFOrGAArR0emhEFVffUXpP/+JucEyOd2J7vdTsXQput+PyeHAkZeHZjZj7obdxqJzkO9MIYTopsIeD4eeeYZwWRmW5ORuuf2VUorKFSvQg0ESRo3Cmp4us2xFhycTFoQQopuq/OILwuXl2LKyulVwU5EI3k2biFRXo0IhIl4v8d/6FvZevTA1s02hEB2NtLwJIUQ3o5QifPw45f/5T7dcXDZcUYF//35smZlYMzJIveaabjmzVnRe8t0qhBBdnFIKIpFoQCl95x0q/vMfNJutW23nFCwsxJKcjDU1lZSrrkKz2YxdVyS4iU5GvmOFEKKLq/j8cypXrqTn3XejmUxUfvEF1rS0btHaVLNrFyoYROk6/r17sefmkjBiRLfqJhZdT9f/yRVCiG4scPgwpe+8A0px+PnnievfH9WgFa4r8xcUULN1a/S+IzeX+KFD27FGQrQOmbAgzthTTz3F7bff3t7VOKWCggI0TWPPnj2t/vrHH3+cH/3oR+daRSFiSuk6RYsWYbLZsGVmEvF4qPj8cyxJSe1dtZiLeL14N27Emp5O8lVXkTx5Mq4RI7pFaBVdn4S3GJs0aRKapqFpGvHx8QwfPpy33nqryetKSkq4++67yc3NxW63c95553HbbbexY8cOAF577bXoeSwWC3l5eTz44IONtqlauHAhAwcOxOl0kpeXxyOPPMLJdtD41a9+FT1n3e2GG2446fV4PB6efPJJfvGLX0QfO/EcdbeGm8B3NXfddRdvvfUW+/fvb++qCNGsYFERlV9+SeDQIcy1e5RaEhOx9+jR5bsMlVJ4N2xAM5lIGDUKs9PZbdewE12TfARpA3PmzOH+++/H5/OxePFipk+fTv/+/Rk+fDhgBLcxY8aQnZ3NK6+8woABAygtLWXx4sXMnz+fv/zlLwBkZ2fzzTffEA6HWblyJTNnzsTpdPLQQw/xxRdfMHv2bF599VUmTZrEpk2buOWWW8jOzuaOO+5osW6jR49myZIl0fuOU0yTX7RoEf369eP888+PPlZYWNjoNc8++yxvvfVWl96Q3e128+1vf5tXX32Vxx9/vL2rI0Qj4YoKDj71FCoUwpyYaAzK7+IiNTWgFOb4eDxr1xIqKSF+2DBZ+kN0SdLy1gbi4+PJysoiLy+P++67j+TkZJYtWxZ9/oEHHiAcDvPZZ58xZcoUevfuzQUXXMATTzzRaKNzk8lEVlYWPXv25Oabb+aWW27h/fffB2D16tUMHjyY2267jdzcXK677jomT55Mc/vDNmS1WsnKyorekk7RnfL3v/+d75ywV1/D47Oysnj33XeZMWNGixvc13VP/vOf/2TUqFE4nU6uuOKKaGDt27cvycnJ3HPPPY1aDnfv3s2UKVNwOp1kZGTw85//nHA4HH3+4MGDXHbZZTgcDoYPH97stX/66afRMvv168dvf/vbFq913759XHXVVbjdbtxuN2PGjGnUpfqd73yHxYsXn/T9EqI9+PbsQYVC2DIzMTud7V2dNlG9cSOVy5ejdB1bZibxw4bJ5vGiy+oyLW8777zzlK9JvOgismbMiL4+9dprSbv2WsIVFey9775THl/3+rOl6zrvvPMOZWVl2Gqn5+u6zuLFi5k7d26zrV4nC1NOp5NgMAjA2LFjmTdvHitWrGDixIls27aNL7/88qThBGDjxo1kZWXhdruZMmUKjz32GIm1XSwnUkqxatUq7r333hbPt3LlSnbt2nVam9Y//PDDPP/88yQmJjJt2jSmTZuGy+ViyZIlHDhwgKlTp3LppZdy3XXXEYlEuP766+nbty9r1qzh0KFD3HrrrSQnJ/PAAw8AMGPGDMLhMKtXr6a4uJi77rqrUXk7d+5k6tSpPPfcc1x66aVs27aNH/7wh6SmpvLf//3fTep31113kZmZydq1a9E0jbVr1zYKpCNHjmTPnj0UFRWRmZl5yusVoq14vvmmy3eNAqhwmFB5Obb0dOK/9S0iHg+ayST7j4our8uEt47sqaeeYsGCBQQCAcLhMD179uSmm24CjC7TyspK+vfvf0bnXL9+PYsWLeKWW24B4KKLLuKll17iiiuuIBKJEIlEeOyxx5g6dWqL5xg7dix/+tOfyM/Pp6CggLlz53L99dezdOnSZrtZysvL8Xq9ZGdnt3jO119/nfHjx9OvX79TXsODDz7IJZdcAsDtt9/OAw88QGFhIZmZmQwZMoRLL72UZcuWcd111/Hvf/+b/fv3s2LFClJSUhg6dCjz589n3rx5PPDAA2zbto0vvviC7du3M2DAAAB+/vOfM2vWrGh5Tz75JHfeeWd0ssV5553HnDlz+P3vf99seDt06BDf+973ov83J15T3ftw6NAhCW+iw4j4fNRs397lF9/VQyGqVq4kUl1N6ne+gzk+HnN8fHtXS4g20WXCW/9XXjnr11uSks74+DMxa9Ys7rnnHo4ePcq9997L/PnzSU9PP+PzHD16FJfLRSQSIRgMMnXqVObPnw/A1q1befDBB1mwYAETJkxg8+bN3H333eTn5zNt2rRmz3fVVVdFvx46dCiDBg0iPz+fdevWNTteLRAIAGBv4RO9z+fjb3/7G0899dRpXU/DjdszMzNJT09vFIIyMzMpKSkBjFaz888/n5SUlOjzEyZMoLS0lLKyMnbt2kVCQkI0uIExnq+hzZs3s3nzZl566aXoY+FwmJycnGbr9+Mf/5g77riDN998k8mTJ3PzzTfTq1ev6PN1LaU+n++0rleIWCtftozSf/wDTKYuuT+n0nUChw4RrqggePQoeiBA/LBh7V0tIdpclwlvHVlycjL5+fnk5+ezaNGiaLjKysoiPT2dxMREdu7cecrzZGZmsnz5ciwWCzk5OdGuV4AnnniCyZMnM3v2bMAIY3v27OHpp59uMbydqG/fviQlJbF///5mw1tKSgqaplFRUdHs8W+//TbBYJCbb775tMqzWq3RrzVNa3S/7rFIJAIYXbYnG3R9qucBvF4v99xzT5NlTiwtLB3wox/9iCuvvJL33nuP9957j3nz5vHxxx9z0UUXAUZLJEBaWtpJyxWiLShdp/yTT7AkJnbJQfp6KIRnzRpCxcVgNmNNT8eVn4/tLD4IC9HZyYSFNtavXz8mTZrEo48+ChiTEG688UZefPFF/H5/k9dXVlZGvzabzeTn55Obm9souIHR+mM+4ZN2SxMGWnLw4EEqKirIzc1t9nm73U7//v3Zvn17s8+//vrr3HDDDS2OmTsXAwYMYNeuXZSVlUUfW7VqFenp6aSkpNC/f3+qqqoaheC1a9c2OsewYcPYuXNnNEjX3Vq6XjC6Vn/2s5/x6aefcskll0Rn/gJs376duLi4RjNvhWgvvj17iHg8XTO4+f1UrlhBqKQE14gRpF57LYnjxklwE92WhLd2cNddd7Fw4cLoEhuPP/44mqYxbtw43n33Xfbu3cuGDRv45S9/GW1JO5Wrr76aN998kzfffJP9+/fz3nvv8cILL3DNNddEX3P55Zfzm9/8Jnr/vvvuY8WKFRQUFLB06VKmTp3KuHHjGDlyZIvlTJ48mZUrVzZ5/MiRI3z66aenNVHhbEyZMoW8vDx++MMfsmXLFj7++GPmzZvHnDlzABg0aBAXX3wxs2bNYuPGjXz22Wc8/fTTjc7x85//nA8++ICHHnqIbdu2sXXrVl577TVefPHFZsu85557+PTTTykoKGD58uVs2rSp0djEL7/8kkmTJrXYcidEW6r84osuuQCtHgpRuWIFEY8H99ixOHJzu8XSJ0KcjIS3dnDJJZfQr18/nnnmGcDoDl2zZg0TJkzgpz/9KYMGDeK73/0uhw4dio5pO5U77riDRx55hIcffphBgwZx9913M2vWLObOnRt9zd69eyktLY3eP3DgANOmTaNfv37MnDmTkSNHsmTJkpO22M2cOZN33nmn0eLAAG+88QZZWVlMnjz5TN6K02YymViyZAnV1dVceOGF/OAHP2DGjBnc12CW8BtvvIHZbGb06NHMmTOHRx55pNE5Ro4cyb///W+WLVvGyJEjmThxIn/84x9bbHkLhULceeedDBgwgO9973tMnz690QzWv/3tb8ycIH/JLgAAIABJREFUOTMm1yvEmfDt2YN3w4auOUlBKUxxcbjHjcOWldXetRGiQ9BOtgJ/exo1apRqbp2u7du3M3DgwHaokahz5ZVXMn369Ji1snUGK1eu5Pbbb2fr1q1NuqvbivwsdC51rcQLFixo1fNGamo49OyzRLxeLG53q567PQUOH8bsdmNxu09rTGtHM/+DDwCYd8K6mKJrCldWMn7BgqNFwWCPtihPWt7EGVuwYEGn+0Xa2srLy1m4cGG7BTchwBgLdvSllwiVlnaJ4NawMcG/fz/e9evRQ6Fu//tGiBN1vQESIuYGDhzY7Vt8Go4lFKK9HP/4Y3z792Pr5OsMKqUIFRURqakhcPgw7tGjcY8bh+73d8lxfEKcK/mpEEKITsa3dy/+Aweo+OwzbOnpnb5lKnD4MN7aYTKm+Hg0ux1N02QzeSFaIOFNCCE6maJFiwgWFmJJSur0i/HqwSDVmzdjSUrCkpqKLTOz04dRIWKtU4Y3XdfPeA0zIbqSjjrRSMReuKKCUEkJtpycLhFyfHv3ogIBXOPHYznJXs5CiHqdLgHFx8dz5MgRgsGg/AET3ZJSiuPHj0e35xLdi2//ftC0LhHcVCSCf98+bFlZEtyEOAMxb3nTNO3bwC+ACwAd2AXcp5T6z9mcr2fPnpSWlnLgwAHC4XD0cY/HA0BCQsI511mI9naq72eHw0HPnj3bskqig6jevLlTd5Xqfn90TJv3m29QwSCO/Pz2rpYQnUpMw5umaf8P+E3t7RGMlr7hQNzZntNkMpGRkUFGRkajx2O1hpIQ7UG+n0VDodJSqrdvR/f78Xz9daddjFfpOpUrV2Jxu0m48EJ0vx9bTg5W2R9YiDMSs/CmaVousAD4uVKq4V+gf8WqTCGE6GpUOMzRhQsJFBSgWa3GxvMn7G3ckSldJ1hYiDU1Fc1ux9m3b3T/1cSLLmrn2gnROcWy5e2HGN2kL8WwDCGE6NLKP/uM4KFD2Hr06HTj3CI1NVQuX45eU4NrxAgcubk4WtiOTghx+mI5YWEisAP4b03T9mqaFtY0bY+maT+JYZlCCNElqHCY8v/8h9J338WSltapgluorAzvhg1UrVyJCoVwjx2LvU+f9q6WEF1GLFvecmpvTwMPAHuBacBvNE2zKKWeP/EATdPuBO4E6N27dwyrJoQQHVfo+HEOv/AC4bIyrOnpmKzW9q7SaQtXVBihTdfRLBYSRo/GdsIYZSHEuYlleDMBCcBtSql/1j72n9qxcHM1TXtBnbDWh1LqFeAVMDamj2HdhBCiQ1K6TvFf/kK4qgpbVlZ7V+eMqEiEqjVr0KxWki65BLPT2d5VEqJLimW36fHaf/99wuOfAJlAdgzLFkKITqd6xw4OP/881Tt2YE1Nbe/qnDHf7t3o1dW4Ro6U4CZEDMUyvG1t4fG6gRt6DMsWQohORQ8GKXr9dYJHj3a6LaKUUiilCBw5gq1HD2zp6e1dJSG6tFiGt7dr/73yhMevBA4rpY7FsGwhhOhUPGvXEqmuNvYr7UTb/1Vv2UKopARN03CPHYtr+PD2rpIQXV4sx7x9CCwFXtY0LQ3YB9wETAFmxrBcIYToVJSuc/yjjzC73e1dlTOiBwIEDh3C7HIBYI6Pb+caCdE9xCy8KaWUpmk3AE8A84FkjKVDvq+UWhSrcoUQorPx79tHpKoKW2Zme1fllMJVVdRs24YtOxtHnz4kXXopWidaNFiIriCm22MppaqAn9TehBBCNKNq9eoOv1+pHgoROHyYmq3GcOa6mbB1uyUIIdpOzDemF0II0bLofqVJSe1dlRYppahauZJweTnmxETcY8ZIF6kQ7UjCmxBCtCPvli2oSATN0nF/Hfv37ydcXo5rxAjsffp0qpmwQnRFnWdKkxBCdBERr5e6Ncorly3D1IHXRAsWF1O9ZQvW9HQJbkJ0EB33o54QQnRBgcOHOfzCC2Tecgu2jAz8Bw5g7YDbR0Wqq/Fu3EiouBhzQgIJo0ZJcBOig5DwJoQQbUSFwxx58UUiNTWUvv02ZpcLzWrtkKHI5HRijo/Hkp+Ps18/TDKjVIgOQ8KbEEK0kXBlJRGzGXt2NsFjxwiVlmLtYMuDeL75BkduLtaUFFzDhrV3dYQQzZDwJoQQbSRUUoKWlgYQDW0dpdUtVF6OOS6OUGkpmsWCNSWlvaskhGiBhDchhIgxPRBAhcOESkowjx0LdJzQBhCpqaFy2TLiBgwgefLkDlU3IURTEt6EECKGqnfsoOj11/EfPIjSdUxWa3tXqZHQ8eN41q4FTcPeq5cENyE6AQlvQggRI+GKCgpfeYVwVRUqEmnv6jSilMK/dy/VW7Zgiosj6eKLZeFdIToJCW9CCBEjJUuWoHQdcwdbx00pRc3Wrfh278aWk4Prggs6XIugEKJlskivEELEgHfrVjxr12JNTW3vqjSilMK7YQO+3btx5OWRMHq0BDchOhkJb0II0YqUruPbu5fiN97A4najmdr316weCFC5fDm+/fsBqNm6lUBBAc5+/YgfNkzGuAnRCUm3qRBCtKLSd9+l4rPP0Ox2rMnJ7V0dNLMZpRRmhwOl64TKynDk5RE/eHB7V00IcZYkvAkhRCsJlZZSsXQp1rS0DrPRvGaxkHTxxdH7iRMntmNthBCtoWP8dhFCiE5O6TqlS5agaVq7Bzfvxo2oUAg9GCR+4EAsDVoA27sbVwhx7iS8CSFEKyj78EM833yDLSurXesROHoU/7590ftqwIB2rI0QIhYkvAkhxDnybNjA8Y8+wpaR0W4tW3owSPWWLQSPHsXsdpMwciSR6mrZ5kqILkjCmxBCnIOaXbso+tOfsCQnt1t3qVIK7/r1BAsLseXkEDdwIJaEBCxJSe1SHyFEbEl4E0KIs1T+xReULF6MxeVq14V4fTt3Ejx6lLjBg4nr16/d6iGEaBsS3oQQ4iwEDh+m9B//wJqSgslma9e6mOLisOfm4jz//HathxCibUh4E0KIMxQ8dowjL72EyW5vt+AW8fkIl5Vh79EDR+/eOHr3bpd6CCHanoQ3IYQ4A7rfz+EXXkCFQu06pqxm2zaChYVYMzJkeyshuhkJb0IIcQaqVq0i4vViy8xs87KVUkQ8HixuN/FDh+Ls10+CmxDdkKzWKIQQpyl0/DjHP/qoXVrclK7j37ePis8+I1Raislmw5KQ0Ob1EEK0P2l5E0KIU1BKUfH558YOCmYzJre7zcoNHj2Kb/duwlVVEIlgSU7GkpraJuULITomCW9CCHESKhKh9P33Kf/kE6zp6W3aTVmzYwe+HTswJyTgyM3FmpaGLTMTTdParA5CiI5HwpsQQrRAKcWxP/4Rz4YNRmgym9us7GBREb4dO7D36oXrggtkT1IhRJSENyGEaIF340Y8Gzdiy8pq09auiM+H5+uvMbvduIYPl+AmhGhEfiMIIUQz9GCQ0n/+E0tCQpt3UwaPHEFFIiSMHt1uW24JITou+a0ghBDNqFq1ilB5OfasrDYpTw+FCBUXY0lOxpmfjy0nB3NcXJuULYToXKTlTQghThDx+Tj+/vtYk5PbrszKSjxr1hCurASQ4CbEOVC63t5ViClpeRNCiAaUUpS+8w56IIAlMbHNyrWmpZF02WWYJLSJbkwpBbre7OQgFQ4TOn4cTCY0sxlrSkr9c5EIaBqayYQeDBIqLkazWrEkJmJyONryEtqEhDchhKillKL800+pXLECWxt1l/oPHkSFwzjy8to0LArRUShdj07KCZWWooJBzC4XlsRE9GCQiNdrhDNdJ/Xb3yZ+yBAOv/ACofJyVDBonETTIBLB7HYTrqwk43vfA7OZkr/+FUtyMuGKClAKzW5HhcNQGxKtLSy9o3TdCIMddFmeNgtvmqZ9DFwJPKaUeqityhVCiNPhP3CAsk8+wbthA9b09DaZ4Rk4ehTv+vVYU1Jw5OXFvDwhTkfDFq66kIPZDJqGNTW1SaCJeL1EqqsxxcVhSUggXFVFpKoKW04OmslExOtF9/uNxaV13Wg50zT0YJBweblRhqaBpmHLyCB92jSOvfYaofJydL8f94UXYoqLw33hhThycwHI+N73KPnHP0j57ncx2+04+vbFs3YtlcuXk37TTSRedBEAnjVrqNmxg6zbbsPRpw9HX3oJa3o6mbfcwrE//IFgUVH0Q5NSinBZmREU61oA7Xbjd0Ftax+6DhYLKhgkUlWFNSuLSFUVEY+nTf+P2iS8aZr2PWBYW5QlhBBnyrNuHcdefx3MZmNZkDYIbqGyMjxr1mBJSjJmlXbQT/ii86vbncOcmNjke1uFw4RKS7GmpaFZLMauHiUlJF1yCfEDB1Kzezf23r2xpaZS9sknVG/eDCYTFrcbk8NB2ONB0zRSr72WqtWrCRYXY05IIH7wYHz79mFyONADARznnUfNjh3R2dN19Ui99loSx48ncOgQStdx5OZijoujx113cej//o/EiRPJuPnmJj8fCSNGGMvoNHg8ZcoUUqZMafS6zOnT8R88iPvCCwHoPXcuJqsVzWIh+YorOPrKK+iBgHHtoRBx/fuTPHkyjj59qN62jZrt29GDQfSaGiOg2mzGuFS7nYRRoyj/9FPMCQnogUCr/7+dTMzDm6ZpScBzwD3AoliXJ4QQZyJw+DDH3ngDS1ISJrs9JmUopRrd14NBPF9/jcnhwD1hgmwuL1pU16XYXDee0nV0n8/Yss3hQEUihMrLQdcx1y5xE6qowJaZiS0tDe/mzViSk1GhEOb4eDSzmVBZGc78fHz79qFpGkopXMOGkX7jjWgmE/FDhkTLy541C/++fQQKCylZvBjN60WFw/T82c9w5ueTfPnl6DU1aFYrwaIiDj75JChFj5/8BEdeHv6DB7H36EG4rIxIdTXW1NToPsFxAwY0um57Tg595s6NXkdzTucDjy0zE1tmZvS+2emMfh0/eDCOPn2IHzLECKWaRvbtt0fHyCWMGEHCiBHN/78ohQqF8B84gL1XL8o+/PCUdWlNbdHy9hSwVSn1F03TJLwJIToMFQ5z7I03MFmtMQlu4cpKfLt3Ezx2DD0+HqXr1OzciW/PHlQoJMFNnJSKRAgcPmyECU1DhcNYEhKIVFcb3YyAvVcvgkeORMeGJY4bh8nlombrViI1NaRdey1JF1+MyeGgau1aShYvxpaTg3//fpSuY7Lbyb79dqMLs7ISs8tldI020/qsaRrOvn1x9u2LvUcPPGvWYE1Lw5mfbzxvMmF2uaL1Sr7iChJGjIh2dTprhwac7nhSa4z38NUsFnrfdx9gtNqpcPi0JzdomoZms9Hrf/4H//79lH30USyr2kRMw5umaROBHyBdpkKIDkYPBil9912CR4/GZHJCsLgYz+rVoGmY3W6Uz0eopISabduwZmYSP2hQtNVBdHx6IIBmtRqzGQMBYxxXKITF7T7lsi5KKUIlJcaYsbS00+4iD1dWknDhhaRdey3mhATKPv2UymXLyL7jDlAKW48e2NLT8WzYwLGFC8n8/vdJHD/eOPiGG5qcz33hhSSMGmW0yJWWUr1lC5akJMwuF2ZoNHvzVJx5edEw1hxN00j/7ndP+3ztTbNYzmpB7LrA2tYfwmIW3jRNswIvA88opXae5jF3AncC9O7dO1ZVE0J0cyoc5sgLL+A/eBBrWlqrnz/i8eBZvdoYZD1+PJHKSkKffQa6TuLFF8e8RUGcm3BFBSoUAosFS2IioaIiTHFxRMrKjD/WCQkkXHghtqwsjr/7LhGPxwhA8fFNzqUiEYLFxSRccAEqEomOGTOeVNHB+ioSwRwXhwqHja9dLlQoROK4cdEPF2nXXUfyZZdhSUhoVEbC8OE45s07re/luuBoTUsjadKkc3ujBIDx/97GO6HEsrT7ASfw2OkeoJR6BXgFYNSoUeoULxdCiLNS/vnn+A8ciEmLW6S6mqrVq8Fkwj1+PGank0hlJfbaD6QS3DqecGUleiCACofRzGbMCQlk3XYbRW+8QeDIEZIuuYSMadM4/tFHRKqqSL/ppmg3e9z55+PdtInyTz5Bs1jQfT5MteOqVDhMsKiI5MsuI+273wVNI1xeTsTrRTOZ8G7ahNntJn7gQMLl5RQvXoy9Vy+sSUlUfvUVJrsdZ79+0XpqmtYkuNWJxYcQcXpMcXFdo+VN07TewIPAHYBd07SGg0nstZMYPEqpSCzKF0KIloRKSyl7/31j2YLWPndZGVUrVwKQMHZso8HRmtkcswkR3ZFS6pTdj6Hjx9EDAWyZmU0WfVW6bnRlApakJHL+3/8DpTj25ptk3Xorzrw8Mr7/fUrffpu0665Ds1hIu/baJmXYe/bE3rMnZpeLojfewNm3L8GSEmMigdVK1syZJIwcWd/ilZIS7Z609+wZPY81NZU+c+dG7yddeinBoiJMNtvZvUGizdS1xralWLW8nQc4gDebee5/a28jgA0xKl8IIZpV+v77AK3+STlSU0PVqlVoNhuJEyc224UmWkeorAy9piY6Xqs5kepqNKuVxNGjqVy+HM1iweRwGH9kIxFCpaUkTpiAa/hwbDk5WNxuAHLnzYsGLdfQocQPHHhaY6ESJ0zA5HTiGjYMFQyS4PFgcjpxjxp1Vtdojo/Hed55Z3WsaHvm2u+fthKr8LYBuLSZx5diBLqFwJ4YlS2EEE2ocJjKlSvxfP11o6UDzvm8ShEsLMS7YQNEIrgluJ0VpevGCvgntpAphV5dbSyL4fdHX5P1wx9SvGhRNJTpgQDh8nJMdjuW5GQi1dVk3XYbrmHDcA0dCkpR9NZbhEpKMNlsxA8dStrUqU1atk5szTvdQeya2RxdSwybTfam7WbaeneUmIQ3pVQF8PmJj9f+UBxQSjV5TgghYiV47BiFr71G8MgRrCkprbYIb6SmhqqVK4l4vZhdLhImToy24IjTpwcCxsr2gMXtRrNaCVdURBdOtWZmYklMxDV0KOGqKhy9euEaPhwViVD05pvYsrIIl5WRdPnleNauJVJTA0Bc//7GWmWDBwOQm5+PHgphdrlkUWTRqqzJyW1anuxtKoTo0iI1NRx+4QX0YLBVJyhEfD4qly9HhULGzMPs7GY30xYtU7oOuk6orIzM6dOxJCdT/Je/EPF6SZk8mZrdu8Fkosfs2c2OF3SPHo13/Xpqtm4lbsAA0q6/HktiIkVvvIFr+PAmrV8mh6NLblIu2p/lDJZZaZXy2rIwpZR81BFCtKmyjz4i4vW2bldpJIJn9Wr0QIDEiy5q80/dnVWkuhqT3Y5msUT3v8RkIq12iySAvIcfjr4+pW5/yRZCsWYykTN7NuHjx42V+E0mEi64gNK33ybhLMeaCXE2LImJKGizVTKk5U0I0WV5N22ifOlSbBkZrXZOPRCgavVqwuXlJIwZI8HtBErXjZXqbTbCFRVgNmNJSEBFIoSrqjBZrdGNv3vecw/muDhsOTnNnkvTNGND9JPQahe+rWNJSiJ92jTihw5t1esS4mTM8fEo0NuqPAlvQoguKXD4MIV//CPW5ORW687UAwEqly8nUl2Na9Qo7C2Eju4o7PGg+3wAmOx2QqWl2NLTjRa26moiPh+JY8eSev31xqxcTSOuwRpmrSnp4otjcl4hWuLIzcUbiXjaqjwJb0KILkfpOkV//Wt0w+7WEK6sjA6GT5wwodsviqrCYSI1NVjcbmMmaE0N2bNmET9wIAA1O3Zg79OHYGEhR19+GRUIkDhxIpaEBFImT27n2gvRuswuF35d97VVeRLehBBdjuebb4wdFFppnFvw2DGqVq9Gs1pxjxvXrYObHgyi+/1EvF5MNhsRn4+Ix4N79GhjSY5a8UOGAGBJSKDP3Ll41q/H3qdPe1VbiC5FwpsQoktRSlH28cdYWmk5iFBpKVWrV2N2u0kcP75b75KgwmFCpaXEDx5sbHButXL05Zdx9OpF6jXXtHicNS1NWtuEaEUS3oQQXYp//35CxcVYW2GSQsTrpWrVKsxxcRLclCJUWkrylCmkX3999LEeP/oRcQMGdOv3Roi2JuFNCNEl1OzaxbHXXjM2F7dYzrnVTff7qfzyS9A03OPGdatwoodChIqLsaSkYHY6ifh8hCsqiB80iNSrroq+TtM0XMOGtWNNheieJLwJITo975YtFL76Kua4OGMdsXPc5UAPhaj86it0v9/Yp7SF/TO7EqXrRDweIj4f6DpJF19M1apVRDweNLOZ7DvuwDVsmOxMIE4qGAFvEFKc7V2Trk3CmxCiU/MfOMCxhQsxJyRgdp77XwwVDuNZvZpIZSXusWOxtvHK6e0h4vMRLi8nbsAAXMOHY01JIW7QINzjxxMqLsbes2er7k4huo6wDocroU8SaBrsLYP95TCxD7jbsLE6EAaTBtY23OSkJgRVAUiLg4gOVruzza5YwpsQotMKlpRw5MUX0ez2cw5uStcJFRdTvX07kYoKXCNHdunAUreYrgqF0GtqyJk9u9FsUQBH7944evdupxqKlihlBKVY05URiFp6TFewvhCOeSElDuxmOFAJJhNsOAYX9T55PZWCQAQc55hEvEFYfRhcNriwR339VO1+By3V4Vzexwq/UWZIhx4JxnXEJ6a02YrdEt6EEJ2SUopjr7+OCofPuXVMD4Wo+uorwsePR5cD6crBDYxZtJrJBEDOT35C3Pnnt3ONuidVu6fSiSGpJTUhIzCN6QkW05mXF4rA4Sro4QZbg1aqiA7m2vMFwkb4qgrAuF5GKKore+VBSHTAgDQ4VAkl1cZxxdVGmFMK4m1QHTRCVcJJ2qIqA/D1EeNarGbjek68pmDEOK/DYpx7/THITTK6ZSM6bC+BQ1VGHUprjHrHW6GnG7YUQ3aC8fWOEij0GuUMTDPe7+2lMPYk7+MxL+w+bpQdZ4WB6cZxYd0oR9MgwQZHPUY92yZSGyS8CSE6Jd/u3QQOHMB6jmu5Raqr8axdS7iiAteIEdh79epyG8wrpVCBAErXMTmdqGAQzWIhd948NJOpW4zp6wiUgohqHBYKKoygM+Q0vo11BRuPGSGl0g+pcWdWfk0IVh02yiupgVE5xmObimqDWk8jmO0vh+M1RhTZdRwuyDbqvrXYCC7lPiPEKQVxNiNEHaw0Qp/TWl/ecV/T8Naw5a64GrwhWHPECD92M5yfaoQ/m9k4/7pCsGgwvrcRpg5Xgj8MY3oYAaqo2ghWJs04pi5E7i036hWMGO/3/gqjboFw7TlNxrXvLYP+aUZZh6qg0AOD0o06bC8xWtb8YeO86fGQEQ9Hqoxj60Kty2a0vLUlCW9CiE4lXFGBZ+NGqr78EpPDcdYD6JVSBA4epHrTJtA0EkaP7rLbXYWKi7EkJmJyOAgVFwOQftNNWM5xYoc4MwUVRmga3cO4H4wY4cikweCM+nabsA77yo2QkJcMfRKN50qqodxvhJTjvvrwFtGNEBPRjZYmi8kIMIcqjX97JRqv+/qo0fLmthsBcPkBo6ywDhpGS9S3MqGg0ghlGnDMA8cTjdBTXG0EFU0zytI0o+6ayQh/VnN9MLOYjCCUm2QEo8qAEQqPeSHbBQPSjRarBJsRaJ1Wox6bi2vP0aBFskbBN4VGoIy3QYUPlh0wQlW8tf59M2lGkKt7T0ya0fq3q9R4z+pa9upa89x2433OdBnXubfceH7lQeiZWHv+2oCmgH1lxvi2feVG0KzTHnN4JLwJITqN6h07KHz1VVQggGazYTnL7lIVieBdv57AoUNY09JwjRyJOe4MmzE6KKXrRsuazYZmMhEqK8Peqxc9f/pTNLudYGEhZpdLglsrCetGELCdorE2GDG64MLKCBQumxHmdGXcPMH6Af77y41Q57DUt3blpxivt5iMUFLkhX6pRqhac8QIT7oywtl5KUZZu8uM4GQ1GS1D3gC4astw2eq6+ozgpJQRjj4vALNWH8JsFiM4qdrX1QUVc4PWQ00zWtgaZhib2RgXVhUwumkLKoxzOq1G9+Vxn1F+nBWsWv0xzb2PujKCq91S32qpq/qWr+Y0rF91qPFrG5ZhMRnvXzhiBDVTbTA9WGFcex17bWDecAx84ZOX3RYkvAkhOoWI10vRa69hstsxn8MYN6UUnjVrCB47RtzAgTj79+/Uy19EqqsJV1UZ11A7hs2amkrw2DEsSUmocJisW2+N7vHaVVsXY6XSb4SQXolNx6XVhGDtESMcnJdsjANrqG7AvFJGGIsoI+AcqjS66g5UGAHNH4YynxHedGUM+o+z1rcU7Sw1QlmZr75lyRs0ws/+ciP4JdiNkLenHNCMViSXzXhsb7nR/eewNq6f7YTWo7pA0vDHwWYGf8g458nG2J343tS1yq08aFxTXTAC4+tA2AhYp/OjZ9LqW8DgzGaUxlmN8lsqx177/jus9fUzm+pDbsPribMarY/2DjCqQsKbEKJTKH3/fSI+H7Zz3DnBv3cvwWPHiB86FGd+fivVLnZUJBINZSoQIOL1osJhMJmiEw56zJ5NuLISU1wczr59sbjdVH39NUd/9zuSJ08+5/esLdV1sbntpz+I/3SFa7vSTve8lX746rDREpPogCSH0b12sMoYL/b1UeMPf5zVaFnKTzECjlJG0ApGjK83FRmtT/E24/6ByvpuwrrWpGMeo3u0tMYIanVhpa4rsK7lqmEI2VpsdJc6a/+SW0xGKNpRYpzfpBmtbtUhIzTGn6K1qKWAc2LoO111LXrNndfeRulD04yWxJM53dmuZhPEncUkkViQ8CaE6PBCZWVUffUV1tTUczqP/8ABqjdvxpaVhaNv31aqXewoXSd47BiaxQKahiUx0Vh7Lj0da0oKpvh4zC4X9uzsJscmjBxJ5i23kDBiRDvU/PQoBRuLINFePy5re4kRbnq5YWjmmQU4X8gINMOzm7YSlfuM7r9eiUZ34+koqjbqaNKMFheTZnSbhZUxk7Gu+7PuWsp8xoD27aVG8APj+o5W1Y8VQwNTxLjOaNekGcr8cMRjjKsyn1B3s6lp8IqzGmPKNK3iIgCNAAAgAElEQVTx6+OsjcOSphktRa0dhE9XJ27U7tAkvAkhOryKpUsBznoWqNJ1qrdswb93L9aMDBJGj+4UXaWhkhISL76Y5EmT0AMB7L17n3a9NU0j+dJLY1zDM1MVMEJMXZCo8BvBprC2m8+s1Y8nOlxltHT1STq9c+uqfvZhhgt6J9Y/F9GNVjJdGYPN+yQaLT9KwZ4yYxJBWpzRclbuM8o1acaA+rquxaMeY5yZyQQOzWgJi2vQImXSjAkGGlBQbnTtKWXMwowGt1oOq9FiV9f9VhewNh5r2kXYEq2F1zX37dGWC9eKtiHhTQjRIdXs2MHxjz7CHBeHd9Oms+760/1+qlavJlxWhuO884gfMqRTLAWiBwKYnE7Sb7ghOl6tM4voxqKm2a76mZUFFUarUV0XYqTBIHSHxQhWzY01g6bdq0eqjJmYCTbjuJ7u+uOOeY3z160/VlBhjDmrmxhgry3rcO0SEKlOIzT6QvUBrToEqPqxUKYTxoDZLUYdCr3GQHcN4xoT7c0HqhO76qzm0+via6gTfP4QMSLhTQjR4SilOP7++/gPHkSzWrFlZUXHd50JPRCgcsUKIjU1JFx4IfaePWNQ29gIV1aScvXV7RLcdGUEnizX2Xe3lfvq1+sCo0sxpBtdogl2o3Wr0Gss9QBNuzmtZqNbcuOx+sH5I7KN48p9RtgqrTEWa+2fBjuPG61XVrMR0I55ISfBCHn7yutbn5xWo5UPjH/jrEaAtNYuIZFgM1oINxU1HmNmPWEB2RPra9LA3cx/1ZkErLNZdFd0TxLehBAdTuDwYfy1C/Cebfdm2OPBs2oVkZoaEidMwJqWduqDWpnSdVQohNL15p8Ph9GDQcxxcSiloteqIsYaDoljx7ZZXevUdT8eroLxvYzuxDNV101pNRtrmsVZjfFcJs1oodpaYoQh6ylmGzosRgizmY2lKtYcqV/KwVQ7O7IqYCw823BAvtVsLJWR5TJmdnoCjScA2MxGi1u8rX68mKbVD6J3NjNAv60G2AtxOuTbUQjRoahIhOPvvw8Wy1kHt0h1NZXLloGmtVtwi3i9RLxetMQkzC6XMengBKHSUjSbjXB5udGVq2lgNqNCIVImT8aSdJoDvs6AUrCj1Agx+SlNB8cXVxszH621i6ymxRldiVuLIdlpdEfWdfmp2oVge7gbn+e4r379s3VHjYVfj3nql2NwWIznTzXLz2ICS4NxXXVdig2Pq5vR2JCtttVu7RGjLs4TBvHbzKcOjkJ0ZBLehBAdyvEPPqB669az3ltU6TqetWsBSJo0CXN8fGtW77TogQCRmhqyfnYva8159EgAy//cRuDgwehrwlVVWDMy6PHjH1O9ZQuu4cOJeL2UffQRSilSv/OdJuf1h+u38zlVF1vdYPychMYD23cfN7oRNc0YTN8n0diSqNBrdDvuKavfZ/Jo7VZB+8qMUFdaY5Q/KsdY3b+kBjYUGTMlh2XWh6FDlfULsnqC8OUhI9w1XH3/bLoITwyadZoLYU6L0SrntDR/nAQ30ZlJeBNCdBjezZsp/+QTbBkZZ9XqppTCu2ED4fJyEkaPbp/gFgoRKisje+ZMDifm4Ss3xladuPeh7vOR9YMfYE1NJW7iJewoheSERLJvv6PZcWbBCHx1yAhDdXssnkyR11iOosxndF1qmtF9uLe88Uryu8tqt/4pN/6F+i2H9LDRxXm0qv6YYMToVh3fy1g81lE7UN9qMoJehd+Y8Vk3li3e2nQ/z7ZgNrUc9oTo7CS8CSE6hHBFBUV/+hPmpKQzng2qlCJcXo5v1y6ChYU4+/fH3qNHjGraWMTvJ1hWji0pkYjXC0DGf/0XleePZN8xI/TUhIygBEZ3nlWFMNlsOPv1I6LXbm/kNwbzV/iNdcgq/EaXZV2GLffVLy9xpKo+vClVvx5Zenz9Qq1bSoyyS2uM9c1Sncb4sYaL1NbN9NzVYLC/or5Mp9Uoq2Grmc1szML8vMC4H2c1ji2oMIKiL2y8pu4cmmZsLC6EaD0S3oQQHcLxjz9GDwaxneE4Lz0QoGrVKsJlZWgWS3TLqzNVt8+jN2iM1YqzGrdAxFhwNc5qDJCPKKMLUQEOLUK4qAzvBZf8f/beOzqu8z7Xfb7Z04AZYAAMegd7E6vYRImUbMoqlkTLkigXNcuxHDtHLic5yUnuXcn1vfec/JGbkxX7OEt2LMtKrOIiW0psFauYKpYpiZQodlEs6L1Pn93uH9/sGQAESLCA9XvWwgIxs8u3N8CZd37l/VHaspvgmqux1l7HAV8pg905l3sA3G5M02ZwOImWiFJ51WpcHg8DcRkRC2ZGIzUPy+PHdFhcJgd7CyFFmFNU3xOTNWNul3x8Z6dcT0NI2nB82COfdzop++PymELk3Pgd3C7Z/TlZtG8qz7HJCvqDXnmvNKF8xRSKmUaJN4VCcV4Z3b6d2IEDRN9/H09Z2Snta6XT0gokGiWwdCm++npcnunN8tEz4kcg/cJGkmOiRchomUeT2wkgaUrhltQzY5awSfb1kl51Ldx0J+3mVgQyAuZOH2/M6goESc5ajDudAEOnf/ZKZiFFVXabTJ2Ybsl04/5M2nNJuawv82ac8m0bDg/A7BLZSOBzy7Rl26gUdrqZ8ydzjnkizob7vhAnH86uUCjODkq8KRSK88bojh30/PSnCL8fLRQ6pXSpbVlE3nkHMxKhcP16vBUV097XsHIWFJqAuAFuTYo0Bw0ZCfNk6qYsWxbAA7htAzHQi7l0DeKm29E0QeAkS9ct0IvK0P/r/4DeLoZL6xlJynWMtaEYW8wf9Eoxlu4abxjreJW1jMiIW/Zxt/x5Og79CoXi4kWJN4VCcV4whofpffJJ3CUluHy+U9rX0nUi776L3t9PcNWqUxJuIOvJIJcq9UxR2D42IuUS4IqNIlJJsG2MjTdjXnMDaNN7Gc3aWXh9UNuI25C1bk56czKEkBG4oWTOsd9ZizPkfGxRvirQVyjOHaYl//+eDw9AJd4UCsV5YeD557FN85SEm23bJD76iMThw9iGQXDFCvz19Sfcx7RlCtQG/JoUbE5K9FQQo8PYbg/6F/8MK1wOBWfmweZzy8YC7SR+Y46Am2BlJoWlagRQXOKY1vEfSmz7xP9nbDtT2iCm/4FGN8c35ti2tM/J90w+4ixlSCNr3YKy/HMv4JR4UygU5xRjdJSen/6U2N69eKuqTmnfxKFDxA8cwFtZSd78+XhKSibdLpqWX6YtvwRjxM/p2FbEoyAE+v3fwg6f3ozVyZjuC74QSqcpzi0pI/OB5wR/o7qZKS3I1GJadm4s2skMmMfiNN9MZDAha0IrAtIgGuSHr5ZhqCrINBBlol8JQzbmgFyHlfl/Xx+SZQaWLUsPAl752pDQ5Vxcn1uuuSdTexr0yskcbpfsoh5JSvEG8lxCSPPqSFqWXPjdsh41ML1S27OGEm8KheKcMvTKK9KEt7p62l5utmUR37+fxMcf46urI7hq1ZT7Os0HWqYZYWJK9GSf2o/DNBGxKOl7v3FWhZtCcarYmSjymdQ0pgxpOVPgyzXm6JlpF2MF10BCipw5JePLB9KmfDxpjKkBdUFTkew27osfv764Lms7BXIWbaFP/j/1ZcygW0agKiijXMeGcx6EeiZ61hOT2wUzx83L2NOALIHoy4i2PHdOSPo0eQ09MdmFHddlpKyxSIq8pCHLEUI++T3gkTWwI6mcTU9Jnpyl63gkdkbk8dOmfK4kT15TJH3um3WUeFMoFOcMIxJh5M03T8mE1zHeTbW04G9qInDFFeP2TWe6QL2a/AQ9WfPBWE5JuFkWYrAPc9112I1zT2FHhWI8hgWxtBQuE/8GJ36gcCJHaTMXlcpzy0hUX1wKkImRLaem8mR/34MJKVB6YjJypI8Zu1tXKAWUEDIVGBpjIRNLS1ETSWfOgxQvPk2uUXNBfkbEORHl7qgUYZF0bhxZV1Reg2FBgVcew+OS12Pa8niWDR43FLkh5JeTQDoicr1zSuSINoegVwoqTeTW7qC5pMC0Mq8PVUF5nsYieW9bR6RwK/DKSJ5LjJ/lK8T4kWyFPnkt1QXy3w5FfikGzyVKvCkUihkn2dJC9+OPY4yMyAHsk8z5nHLfI0dItbSQN38+gUWLxj9nyDeCsW9+UzUfTAsrkxdJxhExabhrLVyOsfHmMzio4lLDzogrbUwn8onsVkaSUrQUeKUYsWwYiEuhoFtynFiRP9fRPFZQOTQWydShR5MCJ5KSgmQ4Kf8fpDIpTCeVl+eRxsxCyGibYcmIWGUwZxytZyJIAa80Y24blUKxPpSZ/6rJa+2IyGibJuT2xf6Mzc4k1+wIN9uWx4/rcvvSfHmPRlJy7UV++bjmgoYx5aM1hccfsyEko2WTndPnnrr8oCATqXNq38ZGx7yavKfmhIaDqWrkPJoUeKdW6DFzzJh4E0LcCXweuBIoB1qBXwH/07btyEydV6FQXFgYkQjt3/0uwu3GXVgIp2AHku7pIbZnD96qKvwLF5I0pHmt06E5nMx0gboykwHOaKE6YqAXXBp2cSnGjTdhzlkEReEzOariHHGiyJNhSSHhdx//vO50HE/4szQtmc4TyDmuBd5cdGUkJUVJfUiKmqEE1IWkuLAzEbOEIT9IRFJSuOW5c1GdpCFTekFvzhrGSf0FPFLgOWbHbpdcixNpK/RJgdeZeRcVSKFW6JP7GJY8fn9cCiS3kNGluA5NmchUoW985Ajk+oeTxwtR3ZLirzQ/JwangxDymBOjikV++XUqaC6oCJ7aPs4aJl7nWE53xu6FwExG3v4CKdj+BmgHVgD/F3CdEOIq27Yn+WyhUCguNaK7d2Pr+pTNBZNhJZOMvvsuxsAAroJCIvOvZDgi3wGcJsuhjJu/80ZzWsIt698BYmgA4xO3Ya69DtweqQgVM0pcl4XjBb6MgIlJ4aNlisVLMkXqY6NcDoaVSZlbUtA79VfFfjkmzMG2ZWQrlTEurgzmIjCjKZnay3NLoRHXpRArD0jxlTbltp0RuZ6x82mLHYHlypkn98dlRM1G7p9ARq18mkz1OdfgdslokmOe3FA0tYicjEKfPGc0nfv3RJyCfZCCcOyotclwu8anDB28mjSDPl1OY0SxYhrMpHi71bbtvjE/vy6EGAQeB64FXpvBcysUigsA27YZef31UxoQbxsGo3/8I0YkgnvBYiJlDbjcbjwT3gTOuD7Y0BGDfeCSR7LqZ2Ou+wR4lMPtmTJVFCyuS0EU8kkh1D4qRUZfXAohp+YpZeaO0R+XIqWhKJeuHEzIKJEjvQW5EV+OZ97YNZTkyWjUSEqOH8v3SNGXMKRwGxvVSRi5iE2eW4qa7qgUemX5Ulz53TnBFPDmCvTzMnVRfre8xqQho1Xh/PHRLK/GcX/Ap2pdM7HBYCJjz1dwajaKiouAGRNvE4Sbw3uZ7+dmYrRCoTivpDs6SHd34ymfXpdmVrgND8PKdYwWV+EWM/Dp3TQQg30Y19+OeeUmsEzwnWIu5xLFtnMp6NO575YtrRRiY2qdnONE01IEBb0yyuYSskg+pudqoxyh4aTbvFquM7I/IYUbSHFUkEkV+tzjxYrj0aUJKcxCfggha8Z6M+PDXEIKq7Hry/fA7OIx9ZMZQVVVII9zsjFiY4WclvHnK5v+5xaFYtqc64aFTZnvB87xeRUKxXlg4PnnEW73tDpL9cFBoh98gDk6ili6ilhx1Zk1H0yGZUEsgojHMK69GXP95swT59ikaQYwM7YKJ7rVCV3WPyV0GQ0q8sv0o23nLBBaR+RjfrcUVk6ReVyXjxVnrBPGRrbsjJ/eQFx24zk1iQMJGUUr8kvBFs6TX5pLpi9tpDibbPaqcx2hMZraKUB3olsnIuCRonBsDZdXG9+pOBlT3b+zMf9VoThbnDPxJoSoAf5v4BXbtndMsc1DwEMA9SdxTVcoFBcmZjRK/3PPgctFbM8ePNMYXaUPDTHy1lu4vF7y1q1nMFB59oWbaSIGerGa5mEuWI656uqzfIKZxbJlCtG0ZE3W2Boww5LpwCJ/rm4pocsUoCvjUWXaMo3npBa1jOjqGJXbNxTJ+q2gVz43kJC1Xo7dgibk93yPFE5HhzKdh3lS2LVljhPI1JUFvVLM9cVlwXugePyaTzVNCFLkTSb0JiJEJtqmgqmKS5RzIt6EEEHgOcAAvjTVdrZt/xD4IcCVV145cRqMQqG4wEl3d9PxL/+CMSrfybVA4IRRNyuZJLZvH+muLoTXR3r1JpI+P5ztdiZdRwz1YWz4FOZ1t563ZgTbnnzcTzQto15TRbUMS9aHOQXz0bQUTjayJiuU6R50TEx7Y7n0IkhR5tEmr7+qCEpx5ZyvMlP/pbnkcUCeqyxfnt+JeOV5ctfhEnLtIf/4iFg4X/7sm6TLU6FQnD4zLt6EEH7gP4BZwCbbtttn+pwKheLck2pvp/173wPLwltWdtLtbdtm9L33MAYHMcOV6LMXIbx+7ClG5ZwWhoEYHZKD5G/cirl64zlREaYlI10jSWm1kO+RAimmw+FBGSEryZPX2RuTqcy8TEoyrstoWG2hFEhDSbmNQD7mdkFPdIK1hE9G4xx8mvy5wJtrApjqsqcqfC/Jk8d1BKKzrUN1Qe7fJ4qInck0AIVCMTkzKt6EEB7gGWANsNm27T0zeT6FQnF+SBw7Rsf//t/Sy63o5APbkwZE9h1A9PeTWrgSUduQa747S9pKjA6BYWKuvBpzzSbskpMLysmYaMg6nJRu81om4qW5cu7xsbRMbSbGuK0LZNoybeW6IoeS8svvlvfCsbcQmbE+eZ6cYPJqUqCF83OdlPUhKQrtjHP8RGE2Nl14JtV80x3qrVAozi0zadLrAp4APgl82rbt7TN1LoVCcf5ItbfT8f3v4/J60YIndtJMGDAUM6DlCN6jH2FW1eOqObv1rWJkCNIp7LwA+lf/8oTzSG1brmkkKdOWNrK2K+jNmQHrmfFbznDqkVTOOHU4KR9zZfzmDEsKuXBezlJCE9JqIphRUdUFMgXp2F1M7HjUXMeP/wlOiF45XZgKheLyZCYjb98H7gL+BxATQqwb81y7Sp8qFBc/qc5O2r/3PYSmnVS4pU3oH0rif/c1RDqFVVaJtWjF2UljWhYkYohkAru8GmPNtdi1jVMKN9OSBflDCSnYXJmRQq7MOp25j343eL3yeUeoFfmhIiCjcc6cx6Qhf/a75fMTOxOLJzjTezRZa1YeULVgCoXi1JlJ8XZT5vv/kfkay3eQ0xYUCsVFSrqvj/bvfhdsG3codNzzNlLQQMaiImHjPbIP9DTGyg0yjXk2lItpIgZ6sKoasGtD6J/+PARzoSvdlILM7ZLfk4ZMWVp2rsg/6M0JLtuWhfmOa75DRaamLBshE6c+5mciSrgpFIrTYSZNehtn6tgKheL8YgwP0/kv/4JtGHiKi497Xrdk7ZdhAYaB5+AH+HvaEYDZOO+EqczpL8IAlwsx0IOx4QbM624BlwvTgv6MI77IpDIn4pinTlaoL8TUjysUCsWFwLk26VUoFBc5yZYWuh59FDMaxRMeP7TdQoq2pCFTkG4stA//gBgexKptAo8Xq2nemS3AthFD/dgIUsFizHU34LruVoQQxNPSWd+wZKelEDLiFvBk6tEydWjKcFWhUFzMKPGmUChOij40RNejj2JFo+gDA7jy8o4TbgCjSSncHINd1+GDuIYHMZZciV1Vd1bWIkYGGV62ie51t5GyZdW+GJRCzbLluccO/VYoFIpTxbJPf0TcuUCJN4VCcUKMSISO734XfWQEdyCAp7wcMcHkNm1Jm4xICtyZLkjR1Ybr2EdY1fWnLdxsIFI1l3hJNVgWwtCJltUTq56LR0BlQEbRxjYMFPpUZE2huJRJZqx4TjYibTImE2VO17lt5+x4Ymm53Vgz6gsJJd4UCsWU2JZF75NPog8O4s0Ml7fsnCUGSMuLqC5f6NzpOO6jBxB93aCnsYvCmAuWn9a5U4EiWtfeTqKketzj3lSUsnxpIuu8ABf6TvcKFQrFhcDYqSJj0TNTRZxxaropX2uc16HJDL0tW04h8WrHC7yEnpup69jtJAwoyryGDGa6ymsLZUPSvl4p4Gyk0bYQ0lDbWWc8fXqj3s4UJd4UCsVx2LbN4EsvMfrWW+jDw3grK7PPDSTkFIB8j/yK6hkhl0rifu91KdoqarDzAlj1c0Cb3itbOj+EsEw8ySjRsgZa1t8BQN3rTxGK9qDf9RB2Xj4ucWJLEoXiYkU3ZdT4Qoz0nClJQ3Z7BzKRLGcKiVeTwi2asd1xDKqdrm8nUuYIrpQBK6qkcNvdIzvFHSGVMmSzlGVDYxF0RaRY87tzZRVCwMpK2Nk1XvitqJLPfdQvm53ml+b229sru8sXlUkboYG4fNyJzKXMjMi0HQk68yjxplAoxmHbNoO//S0DL7yAJxzGW1mZnU8azww792ryxTiuZ14AR4dw7/sADB1jzSYoOPmUBQdLuIhUzqF17WewXS7yhntIFFXgiwzS+NbT+FsPkf7Sf0Xk55+t4QuKS4ykId9IfWfwjmbbp1bf5Bg4+06jAcbMCAy3K3fOhJ4TNWlzfM2mYeWMoKezLsHZEYBm5rzTuS+6Of56HEHmrNmjwexiODYsBVwiUxsbzaQnS/PlXN1DA7lZu4U+WFwu1/FOxhk2z5Obv9sVlZF/J9VpWLCgVEbGKoNSwO3pkdZAzoi68oD0WGwIQcuIPHddKPe3c0XF+OuqD0FNQe4+1Iekv2PLsDz/yip5bdE0RIb6Bk73Xp8qSrwpFIostm0z9NJLDLzwAt6KCkQmamYhX+SGk/ITKGREm23jOnoQ19GD4PFiLl0zLeGm+4OkCsKM1C5koGk5uDT8Q90E+ltJFZZR+vF7VHzwO9yDvZgrN2DXNs3YNV8q6KZ8wz6f9X5Tpb7OFCcN5nePnyxh21K4eTW5TdKQ29i2jOpMV/DE9dza8zwn38e0pBAoyZPiYez8VkfwTCSp5+6LaUuBEctEhUxLXsP6Oinc3mrNiUnDktdlZ2o6PZrc3sgIQI+WiyA598O25ZrGehfC+N+LbuaiVAJ5HM8Y8eXUkdp2xrB6jJh0asRA7ut8mBNCCjMhZDSqJE+KJIScJOLVZATtyKAUa1dWQ+uI/FpUJkfH1YVkvZlvwu+6oUgKu1XVuTVeUQ5/bJfbezK+jPWh3P0IemFdLRwZgo8H5D2ryczknRuW1zAYlyLvRIwVwkJIUXlFhVyz85zfDaaR1k98pLOHEm8KhQLIRNyef56B55/HW16O0DRSprT+sDOGtoadq3VjdBj3/g8QkWGs6nrM+UvBfeIWz3ReIb0LNjDUtAzbpYFlUdL8IfkD7YTa9qONDiJMAwwDu7wKY/0nMT+5ZeYv/iLHtI5/85wOti3fwMe+aZ8JCUMOu3AJKYLOloiL61BbAH1xKZq8mhQHhiX/LldWSaHxVqt8LKHL88czM2gtO9fQ4tHkNnFd3ivdkvNkV9dA+wgcHc6l6UDeG68mtwe5v25KMTEvDG+2yDX53PLxhJH7gAO5e+yk7UxLirQiP7SPwp5eud3SCrk+vxtK82TtlS8jipZWgFvA3j4pVNwued1+N3RGwNJyAiqcJ+fgHhqQ1wdSJAqkoLPtnGisKZCPWRZ0ROQ1uoTcxuWCDXVy3X9oldfo0eRxYroUX7OLYX+fFLALy6AvJsfH+TL3eE6J9FMcy8JSOR4u4M1E40rkl4NLSBE3kblhub6qMVUTeR7Y1CDTo3t6oa7w+Bo4IeQ6KgLQG5P3BuT9urJa/k5O9wPP+UxvK/GmUFzGGJEIvU8+SaqzEzMSwdZ1Kdzc8qVhJJn7ZD40JupGKon7/T+AcGFcsRq7ouak79SDjcvoWHEjCEHJsV2E2g/gjQ3hjY3IIfJpHatxDrbbg7HpZuzqhpm9+IsIZ1KF8ybjpNY0lxQXCV2+uQ0npdjOc49/YzGsXCrK5869wTlp72happXMjLg4mc2KE9nSMqkkR6C4BKypk3VDTqpq4hujU/Q92WxWJ6KVlxFVTnTJiXRE0vDHNnntJXnyzbq6MHeOuSXyTXxeWEZFuqPQNiK31QTs75fbJnSYVSxTXwEPLK+SUZoFZTKl1heX98UtZNRmNCVFUVMx7OqWa6oukNssq4TtbfJ+pwxYXiH30TPm0LGMILqyWn4fTcpxaSCL4ntjcpvyMSJnQRns7JRCcHGZFFlCSOERSY8fwVZbCB92S1GnuWTasMCXqffKROGK/fL3FU0BQgqZlVXj/8vOLpHC69iQ3GdeaW6m7ooqGbkaSeWOt6pKnm91jRSKs4plBOvooKwLKw9IgTcRIY4XdNPBmxHYE9FcUBuS93Sy5gWHAt/kovBi7UxX4k2huEwx43Ha/+mf0AcHcRcW4i4qyqZJQb75pMzcC6LLBoGN6O/FdXg/mCbG2mvGjaI67hweH22rbiEersHIKyDYc5TaHb/FmxjNbiMGerHqZ2Nsvh27+uwOqb/Yse1cd5yNfDO1Mo9VFUhREdPlm+HsYrlNyzB8NJATTk5UbkFYRlOcmiIbKVzW10FrZh8hwK/l0o8T15I2c1Eor5ZpVslEogJeWFQqhdLqGllM3huT6/VqudSgmYlCRY1c0bcT5bJsCPnlzNmAV6bcgl4plDSXFC0LSqVQWVR+/Buvk/5yvlcGc/VRti2jS6MpGSWaVSzvmVcbL2KK83LiCqAkHz7ogqWVci1JA9pGIZQRAiV5MlXXPALLKqSQKAtI8eICbpg9vm6sZMyxhZBRNScy51Dog02NOSHr4NHG7w9S0F3XJP8OxorlxeWyTsyrSeHrEvJ3EklJYTvxs5bIjHtbUcVxlGeuZzABB/qk2B2bLlw6pk5sQZn8OteMTV1fDijxplBcpgy9+ip6fz/eitwrr9MOr1u57i8HERtF27cT1+gwtseLuWTVpMJttGouw7ULMX35JIoqML15hNr24yo/MJsAACAASURBVIsOUn7wbYQ9Zl5VLIIdKkH//NfAc3m9+tp2LmUzWdDSSW+V5cs3x7czEae0KYXHgjK5jTM5wmF2iTzux4MyAhfX5f51mfGzNQXQE81FQLwazCqRxdfFfqgvktGtZCYFOrZ9zhFuLgHXNMhI3+4eGX1yokMgBb+Tyjw8KKNAQa8UcwtL5Xn74zICZGXqpwYTUgg0hKB5GGoKJ/fxajp+GlsWzTU+BTcWITL1UuSK06fT4BD0wtX1uWtrLJJrG/s7W1AmH3NEn8+dK6I/WWrNo01uNeGknqeDELkomUNJnhSAY20t1tbI/9uTRT2nQ0kebFCfry4IlHhTKC5D0j09DL3yClpJmL64TCn5MgXUjgGmJsakIUaHcL//NgiBsWgFdlW9DONk0H0BkkUV9M1dQ7RyNu5EFE8iQqCvldLD7xEYaB+/AENHREbAstBvf+CyEG5juwp9bhk98mfE1Vi7A5CiLK7LN8sVVfL3MD8M73fLeqg54dwxJ3vjn10iBVHbiIwu1Y7R2B5NRofG4hKwtlb+zjWXjMZ90AXFQZmazPfk6sZeb5ECy6vJaMwnm6bOmLuETGHOC8tr6ovLa3JnomgOIXIRMmf9M8HpmLrC+OubLOXrdo2P1l0o5E8Qf1OlqxUXH0q8KRSXEbF9+3AXFdHz5JPgdtOfdpPOmGDGnILssZGCZAJX2xFcLYfB68e48mrIz73Lmm4fXUs/weCslQBoqThVH75C6eH3xkfYxiCGB8ClYS65EnP1JuzK2pm41AsOj0sWWxf5ZQSqsQiWlMsaqvZRKQAcqwiQImdlVU5AVxbANV6ZUjtZI4Ark46bU5LzqjoZY9/UnbTdZKytGR/lmW5TghDj67oUCsXpo8SbQnGZkDhyhM5HHgEhsG2bdHEFqeTU7uCipwNt7w6wLOyqOsx5V2B7/QzXLSQdKMadijEwZzXJglJKP36Hgu4j5Pe3o5mTdMvbthRt6TR2QYj0n/wlFE7fC+5iwmkOmOjEDjLitqZGFvMX+eVzi8rk76A0X3buFXplTdHEKJFLyHqwU2Fi5OVsEJ6kCF2hUJxblHhTKC4DjOFhun78Y1yBAOQHSFuCkYSYsh5HdLai7duJXVSCuXgV5AexhIu2NbcxUrc4u51LT9H01lMU9DYffxDbhlQSkYxDOoU1exHWvCVYDfMuCuHmeGaNjVw5ru+OEarN+M7OdMZpfV5YpkXbRidPXY0tOve5ZQQOZEfj+fZqUygUFz5KvCkUlzj64CCdjzyClUhghMIMxKToEEwQCekUpBK4+ntwHd6PXVKGuXwdtubG9PhpW/sZIpWzqdzzGqUfv0s6UISWTuJJxcafMJVEREfBtrFDJZhLVmOXV2GuvBrc5/clx5mH6KQIY+ncvXAsL0AKLmcEWFwf/1zIL6NlAY80Fz08CKaRMytdXSOjaLYNhX5pm2FOnkE+jvMxI1GhUFx8KPGmUFyiWKkUQ6+9xtBLL2G7NKL5YaIxGdkZayKKoeNq/hhX8yFERqFY4XKGN99N7xXXkSiqwPL4wLao2flbwsd2AeCPTDIJJh5D6CmM9ZuxFq2Ylv/b6ZLORMB82uQdfU4UzKvJfwuRs6twHNm9GmxskB5YTgPBkSFZg6Zlivgh57Pmdo2vDZtdIrsMnRE9hb6cIBZC1rVVF0h7ht/NQApToVBcnijxplBcYti2TWTHDvqffRYzEoGiEgZ0D7o+PrIjhvpxdbYgersQho5VVYdZVoVVUETnVZ9lYM6VeKNDFLfsxhMfJTDQfnzXqINlQTKOiEbQ7/kvWLMWnNbaTWt6ruV6ZmB1TYFMTTqRszy3jKQlMj5kTcXSpqKpKDeupzwA+/rk4ysqM1YN5GrMllbI45rW9LoT/W4p0KbCq8k6sUtx2LhCoTg/KPGmUFwiDP7ud8T27MFTVsboO+9gFhSRDlUQSWVmF46dzzfUj7bzD6Bp2GWVmHWzSNTNx/T46Fm8kWh5E6WHtlO55/e4pugaBSAygkgmQAisqnrMDZ86beGW0GUkTRszWsmpOzMzkS9vZv5iypAWGlUF0hG/PyatMXozGdy6QukQ79WkqetEVlZBNJwbHzQRVZSvUCguZJR4UyguAdLd3Qz89rcIt5tkWxtGaSUDSRciPclw7lgUbdd2yA8weMuXiNbMJR6uJV5aB8gmhNodv6Gk+cMTnzQ6Cl4fxrW3YM1eiF0yfVt1Z7ySYeXqwTQXXFUHB/ulC75jCFsZlHVmAlljFs9MF3B8wRz3d8g5/59oTI7DRFNThUKhuFhQ4k2huMjRBwfpevRRLM0NhSWYNgwkZARrnGgzTcRwP9pHu7E1N4e++X1iDYsQpo43Nkz1rt/hiQ2TN9yNNxGZ/GS2DdFRGW0D0lu/gl0/55TW60xwyHPLwv7SPPC65c8hvxzl82aLjKQtLJMpT6dsrqlYpkunEmczYY2hUCgUFxpKvCkUFymWadGx7W1Gf/MfWIbBSH4YOy6fGyfcDB3Xx/twdbUhTAPTm8fhh/8Zo7ichrd/QWHX4SkNdSciBvqwqmoxF6/CDhROS7gZFsTTUnBZSIG1okraYkzWyxD0SgHndUF5cPxzx0URFQqF4jJEiTeF4iLl4LadDP/7k1hFYez80PGRNoBkAveONyERIzlvOR2b7yNZN5eyIzso2vs7XJY5/RMmE9iBIPp93wLP+BCXaUmR5s54lCWN3MSAhC6HiPfHpWBrLDp58X7t1LPuFQqF4rJHiTeF4iLBSiaJHzoECPo9xQw++yxacTFu/xS2+5FhtL07sQ2Djgf+lt71txHsOcrsN5/ErSend9J4DJGIYQcLESNDGLd8YZxws2xpw2FacrbjaBLSlhRpMV0Kt6Yi+TXrBAPFFQqFQjF9lHhTKC4CUu3ttH/veySiSRKGtMpwC6C0fPyGto0YGcJ17CNc/d2YXj9H/ux/ka6ZTdlHf6Ry77Zpp0ixLERsFHPZOlwf78NcsR7ziiuBTAeoJaNsxXkwt0R+T5swEJcNBM4gc596lVEoFIqzinpZVSguUExL9gdYsSjtj/yQoQQkguW4mCTtaFno0Rje/TtwR4axvH46b/kKfdduperAGxT/7pWTnzAWRSTjsmtUCKyBfoxFV8Jt92BasoFAAKSlxcaSajmfc2yq1qvJTlAANSxAoVAoZgYl3hSKCwzDMIhEIkTb+3n/lfcp3PEqQ/0RUqHScV5t2DZWZBRXRwtaTzsePUWsbgGtW/+C0UXr8CdGmPXHX+Af7Tv5SZMJRCqJNXshrkN7sGyBCJeR3nwHhi6F5PxSCPmkqe3YSQIKhUKhOLco8aZQXEDous6//du/ceDAAepSFqmf/it9Xh9W8XjhZlsm1sF9+DuOYOQXMrjmBtI1symJdFDduZfqzj1oRvrEJ0slEaPDct6oZZG44yukmhaRd3Qv1qG9XPHZT+EKF9A8DAEvzC6esUlXCoVCoTgFlHhTKM4ztm3zxhtvsHjxYrZt28b27dsxDANwQbgMl8+HC7DTSUxLoBeGoeUIBZ1H6bvhHoy62djePEoPvzftejZraBDbpWFt/DTu3g6Sa64jWTWLAi8MNS1l+fqlhENy2+K8Gbt0hUKhUJwGSrwpFOeZl19+mZ///OeEw2Gi0SgNDQ18+OH46QZ2ewvuw/vw6ikcLRW58hMUWTFo2X3yk1gW6Gnw+bGGBiFQQPhr36SF4qzp7ZWVUqg5DQcKhUKhuDBR4k2hOI98/PHHPPfcc8yZM4d4PE5BQQGpVAqXYYDtwdLc6JEYwQPvM7x0I1ZVHe7oMD5bxx+YpsKKRRDxGHZ+EGugHS1cyuK/fJhgRTFzLBhOyho2b6bDoCJ44sMpFAqF4vyixJtCcZ4YGhriscceo7S0FI/HQygUIhKJ8KMf/YhlwWJI9JMciVK063Ui81bimjufwGg/5HuByQdzOvYc2fFRloUVjRG555voZXWU73udK65dTl5FGJDblaoh7AqFQnFRocSbQnGOOXz4MEePHuXll1+mv7+fDz/8kP7+fmpra/noo4+IRCKMVjZQAIT2vEVi/nK81XW4RvtPeNyscLNNDNOFZ3QAM63jWbmWDVfPoT8GdYs/ddLpBgqFQqG4sJlR8SaEqAP+CbgeaRH1CvAt27ZbZ/K8CsWFSDqd5sknn+S9997D7XZjWRZPP/00tm1TWVnJjh07qKmp4Y477mDXrl24CgvR12/G4/NNfVDDwIyM4nK7Mf1Bio0R3HqS4ZQgMXsp7qs2sXz9LDxeOTNUoVAoFBc/MybehBD5wGtACrgfsIH/F/i9EGKpbduxmTq3QnE+SKfTeL2TKyTTNHnyySfZuXMnDQ0NfPTRR/zqV7/C4/Hw8MMPU1JSkt02lUqxa9euSY9j2dKuQ1gWDA1g4sJ7xTJSfQMUdB/B59eo++9/he3Pp8ddQnlQTBxDqlAoFIqLnJmMvH0FmAXMt237MIAQYjfwMfBV4H/N4LkVinPK6Ogo3/3ud9m4cSMbN2487vlt27bx7rvvEgwGefTRRzl06BClpaU88MAD44TbiTDSOu6RQSwhsE0T+4rVlN/xWebUFWCZJv2/+AUuvw9/XR0AjWfzAhUKhUJxwTCT4u02YLsj3ABs2z4mhPgDsAUl3hSXAMlkku3bt/Pmm2/S19fHr371KxYsWEB5eW7maFtbG8899xwDAwP86Ec/wu12c9ttt3HVVVfhdp/8v6BhAakEnugI9bfdALEY/vnzCS1fhisz5kBza1R8/nMzdZkKhUKhuICYSfG2GHhuksf3AXfN4HkVijMmEonwzjvvUFBQwMqVK/FMknuMRCL84Ac/4NixY5SUlNDQ0EBnZycvvvgi9913H6Ojo3R2dvLjH/+Y/v5+fvnLX9LY2MgXv/hFQiHpgDs0NMTw8DAul+wiEELg9Xrx+XyQsrGHBwkE/fiMJLXf+CqFS684p/dBoVAoFBceMyneSoChSR4fBIon20EI8RDwEEB9ff3MrUyhOAGdnZ388Ic/pL+/H13XaW1t5a67cp83enp6KC8v5ze/+Q2tra3MmjWL5uZm3n33XTZt2sTOnTu56qqr+MEPfkA6ncblcvHcc89RVlbGgw8+iN/vxzRNuru78Xg8fPOb38QwDEKhEC+++CItLS0sXryYaAJqttyCb7iHwIIFSrgpFAqFAph5qxB7ksemnI5o2/YPgR8CXHnllZPtq7hEGB0dxefzcfToUVpbW8nPz2ft2rW43e5sFGoqTNPk8OHDNDQ04Pf7z9qaTNPkrbfe4plnniEYDNLU1IRpmrz22ms0NDSwaNEinnjiCXbt2sXmzZt5++23KSkp4fnnn2fbtm1YlkVXVxcbN27k8ccfB+SHkMcff5x0Os19992H3+8nlUrR1tbGNddcw0033URxce6zzP333086neav//qv8fuhccuNZ+36FAqFQnFpMJPibQgZfZtIMZNH5BQXGb29vZSVlfHmm2/y+uuv09jYiKZpbNmyhUQiQWtrK7FYLFsPtnbtWjZu3MixY8f4+c9/jmXJOZyapmGaJs888wymaVJUVMSNN97Ihg0biMViPPHEE7S2tlJTU0N9fT179uyhpaWFsrIyvva1r1FVVXVG19HS0sKbb77JkSNH6Ovro7q6WqYtM2urqanhqaeeor6+npaWFhoaGnj55ZdJpVL85Cc/IRqNsnLlSsrLy3nxxRdZtGgRAE1NTWzfvp09e/Zw8803U1FRQSQSob+/n3vuuYcNGzYctxZN08jLU8NEFQqFQjE1Myne9iHr3iayCNg/g+dVnIR4PE5vby91dXVompyJNDQ0xMGDB3G73RQUFNDU1ITP5yORSNDc3MzChQtJp9PZyFh3dzd///d/z8KFC9m3bx8lJSXs27ePWCzG8PAwHR0djIyMoGkaxcXFVFZW8s477/DHP/4RwzCoqanJCiQH0zRxuVwkEgl++tOfsmPHDkZHRxkcHKSsrIy2tjY+/vhjAoEAc+bMoa+vj0cffZS/+Iu/mFYErr+/n5aWFvbs2cPBgwe59tpriUajbNu2Db/fT2FhIU1NTePuUzKZpKSkBL/fT2trK3WZTs7y8nK++93vAvDtb3+bmpoaTNNk27ZtNDc3s3r1ag4dOsSvfvUr5s+fz6ZNm4hGo4yMjPCNb3yDefPmna1fp0KhUCguM2ZSvP0H8P8JIWbZtn0UQAjRCGwA/vt0DjAyMpIt7D5bJBIJ2tvbGR0dZf/+/fT29vLJT34Sr9dLfn4+ZWVlBCbMjEwkEui6TmFhIQDNzc0UFhZO2+LhfJNOp/nwww9Zvnw5R44cyUaL/H4/+fn5zJ07l7179xKNRrNizu/309TUREtLC8PDw9x6663s2LGDefPm8ZnPfIann34ar9fLwYMHKSkpyd6b0tJSDhw4QH5+/jghBFBbW3vCdTrndvbt7+8ft99EsVdWVkZrayv/9E//xE033cTy5csnPa6u6/zsZz/jnXfeQQhBfn4+wWCQF154AZCpTefcDocPH+aJJ55A13X+6q/+itLS0uxzkUiExx57jGQyyZ/92Z9RXV2dXf+cOXM4dOgQx44d47HHHqOiooJ77rkHkLVyDz/8sBJuCoVCoTgjZlK8/SvwX4DnhBD/J7L+7f8B2oAfnGznVCrFd77zHb761a8SCATYtWsXXq+X66+/HiHEuO0OH5ZuJI4VQzKZpKCggIqKCnp6elizZg2hUIjnn3+ePXv2ZPf1er14vV5+/OMfI4RACIFlWdTU1DBnzhzWrl3Lrl27ePXVV6moqOC//bf/RiKR4Pvf/z61tbU8/PDDJ63Pmozm5mZ6e3uZO3cuxcXFjI6O8tJLL5FKpVi3bh2zZs067ri6rrNt2zbWrl2bFUoA0WiU3/zmN9TU1LBhw4bsfrqu88tf/pIjR47g8/k4ePAgq1evZv/+/YTDYcrKyjAMA13X2bNnD/n5+ePsLVKpFJ2dnQQCAYqKivjNb35DKBTizTffZMeOHRiGMWlTiRCCxsbGU74nE9E0jaKiopNuV19fz8jICI899hh/+7d/S0lJCbqu4/F4EEKQTqf5yU9+wu7du2loaBh3Xydbv23b7Nq1i6effpri4mLi8TgvvPACW7duZXBwkF//+tccOHAAIQT3339/VrglEgn6+vqoq6tj7969PProoxQVFfHQQw+Rl5dHS0sL119/PQsXLjzje6NQKBSKy5sZE2+2bceEEJ9Ajsf6d2SjwqvI8VjRk+0/ODiIbdv88z//M5qmoWkaqVSKeDzOZz7zGSDnWv+HP/whm5LatWsXxcXFtLe388EHH+DxeHjnnXdwu934fD5qamqOi7KMjaBZlkUsFuOtt95i27ZtCCGorq6mra2Nt956i/b2dtLpNIcOHWLnzp2sXr36hNcxPDyMZVmUlJRg2zbvvvsuP/3pTwFwuVysWbOGnTt3Yts2mqaxfft2qqurueeee7LiQtd1nnzySX7/+99z7Ngxtm7dyq5du3jppZeIRCJompbtivzCF76AZVn85Cc/4cMPPyQcDhOLxVi4cCEHDhyguLiYYDAIgNvtxu12T1pj5fP5xkW65syZk71Xuq6f1UaBMyUUChGNRnn22Wfp7u6mu7ub0tJS1qxZw0cffcSRI0dobGwcJ/on46OPPuKXv/wlQ0NDNDU18eCDD/LKK6/wxhtvsHjxYn75y1+SSqXYvHkzS5cupbq6Gtu26ejowOv18qlPfQqAF154Acuy+NKXvkRBQQGRSITCwkJuvFE1HygUCoXizJnRbtPMDNM7TmdfXdepra3N1kEJITBNk5deeonCwkL6+/t59913SaVSNDY28swzz2RTnhMjNpZlZcXRyXC5XBQUFFBQUDDu8crKSp566qmsAEwkEvz2t79l1apVU0bfmpubeeSRRxBCcPfdd/OHP/yB/fv3U11djd/vR9d1du7cSVlZ2bixSgMDA/zDP/wD69evZ/bs2bz++uu0tbWxePFi9uzZw969ewGoqKggHA5nI4Zvv/02fr8fIQS7du2iqalpnGA5G/YrjpC+0HBmgxYXF9PY2EgkEuGll17C7/dPGgm0LAuXy5XtEP3www/Ztm0bFRUV3HXXXVlvt+uvv57du3fz2GOP4XK5+OY3v0lNTU32OO3t7cybN497772XgoIC5s+fz+9//3uuvPJKysvLMU2T3t5eHn74YdWIoFAoFIqzwkxbhZw2jugYKxQ0TaOuro6f/exn5OfnU1lZOc6h3rYndxc5ndTmRJzaMIdgMEhzczOHDx/O1jAlk0nS6TSFhYV88MEHPPbYYxQVFWEYBo8++mi2lsu5No/HM04IOITDYYqKiti5cyfbt28nEAhkBchUKUmXy0V9fT1vvPEGtm1TX19/0kjT+cK2bQ4cOMDHH3+MbduUl5dTXV1NTU3NcWa4yWSSjo4OIpEI6XSaw4cPk06nCQaDFBUVMWfOnGyXqxMdBCYV4MPDw/znf/4nhw4dIpFI4PP5MAwD0zQRQrBkyRLuvvvucVFFv9/P1q1b+eEPf8inPvWp7O/LNE06OzsJh8M88MAD5OfnA1JEXnvttTQ0NADQ0dHBxo0bVbpUoVAoFGeNC1a8TYXP57tgCr4DgQAvv/wyc+bMIZVK8cgjj9Dd3c3atWt59dVXqaqqykZbwuHwKR1b07RsPdV0cbvd58Xc2LIshBAYhpGtHTx69Gg2fX348GGSySRLliyhpKSE5557jpaWFrxeL0IIUqkUQDZKZts2lmWh6zrt7e0YhpE9VyAQyArnWCzGiy++yJYtW7jmmmtOuMY9e/bws5/9DNM0WbVqFcFgkFQqhcfjobKykjlz5oyrJXTo6OigqamJv/u7v8s2sti2TXNzMxs3buSGG27ICjcgKyqTySQg6ypvu+22M77HCoVCoVA4XHTi7UIiHA6zb98+nn32WQ4dOkRXVxehUIjXXnuN+vr6ac2tvJAwDIPDhw/T29tLYWEhS5YsyV6DbdskEomsUOnq6uLgwYP09PSwa9cuLMvCsix8Ph95eXkMDw9nj+v3+3G73ezatSv789atW7MpZ8daZO/evXR2dmZTs5qmsXbtWhYuXEgoFMLtdhMOh7OR1EQiwRNPPMHzzz/PwoULx3WEOiSTSV599VV+//vfU1dXxz333DNtId3X14emaXR1dWUtQpxrX7FiBXffffek0c3Zs2dz4MAB4vE4t99++zhxp1AoFArFmXJxqYsLDJfLRUNDA6+88gqhUCgb9Trb9iZnimmagPQt27lzJ7quY9s2/f39DA8PEwwG6evro7e3N7stQHFxMcuWLcPlcrF//356enrYsmULgUCAn//85+i6jtfrZdWqVQQCAbxeL4ODgwwNDfHpT3+aRCJBWVkZs2fPBqC1tZW+vj7mzp07ri6xuLiY4uJilixZckrXlZeXx1133cU//MM/8MILL3DvvfcC0N3dzX/+53/S1tZGPB4HYP369WzZsmXaglrXdWKxGA8++CA//vGPszVynZ2d+Hw+7rjjjinT0vPmzeO1116jrq6OdevWndI1KRQKhUJxMpR4O0PcbndWnFxI2LbNkSNH+OCDD9i9eze6riOEQNf17DYFBQWEw2E6OjooLS1lwYIFNDU1UV9fT3t7O6+++ipvvvkmtm1TXV3N7NmzefbZZwGorq7mwQcfpLCwcNo1hY2NjWfFRmQsoVCItWvX8uabbzI4OMj777/Pyy+/jM/nY+nSpYTDYWpra8fVK56I0dFRotEo8Xicz33uc6xYsYJVq1axd+9ePB4P5eXlPPTQQ+NGWk2kqqqKuro6vva1r6kmBYVCoVCcdS5b8WYYBu3t7Qgh6O7upq+vj4qKCoQQNDU1ZVNrtm3T2dlJQUEBHo+HQ4cO0dnZmTX2PRuYpsmhQ4cYHByksrIyW4A/HRKJBIlEYpzdSX9/P7/4xS84cuQIXq+XJUuWZGu8Nm3alO1Qdb4mY8GCBSxYsCAbiXNGWO3evZv8/HxmzZp1XHPB+eKqq67ijTfe4B//8R9JpVIsW7aMz3zmM8c1LEzEEbRONM4Rbhs3bqS8vJz169cD8OlPf5oPPviAeDx+UuEGUqT+5V/+5aQ1dAqFQqFQnCmXtHhzuk+FENmaqv3799PV1ZW1HnFwbCMciouLKSgowOVy0dzcfNyxOzs7eeCBB07ZNiOVSvHOO+/Q1taWnSAxPDzMsWPHstuUlpayZcuWKTsUU6kU7777Lu+99x6dnZ0ANDQ0cMcdd2Q7Wy3L4vbbb2f16tVnJDIndvuuWLHitI81U4TDYZYtW8aRI0f43Oc+xxVXXDHltqZp0tHRgWmaeL1edF3H7XZnrWS+/vWvHxdJLSsr49Zbb6W5ufm4qRFToYSbQqFQKGaKC1a86brO3r17CYfDdHV1oWkajY2NU9aTDQ0NcejQId577z1M00TXdXRdR9M0SkpK6Onpybr/b9y4ESFEdixSYWEhVVVVDA0NYRgG+/bto7e3l6GhIUZGRrjlllsAKQZramro7e3l2WefZdeuXaxatWra19Ta2srPf/5zuru7KSoqori4mEOHDmEYBlu3bmXevHkcO3aMl19+mccff5zvfOc7WaPc3t5ejhw5wtDQEH/84x9JJBLU1dVx0003oWkab7zxBt/73vcwTZPCwkIeeuihcRMTLhUsyyKZTJKfn49lWXR0dOByufj85z+PEOKEKVzTNGlububqq69m7ty5zJs3j2g0ys6dOxFCcPXVV08ZVdu8efM4sa9QKBQKxfnighVvQ0ND/OQnPxn3mBCC6667jptvvvm47TVNIxgMUlZWRkFBAV6vF4/HQyqVoquri7Vr12YL66fC6VasrKw84drmzJnDiy++SHNz87TEWzwe5xe/+AV79uwhEAjwla98hfnz5wNSUNi2nU3drVixgmAwyA9+8AM+/vhjlixZwvbt23n22WezlhlLlizhuuuuy3qJOfv97Gc/DmVW7wAAFgNJREFUIxwOc9NNN10UHY6O8CovLz8ugmnbNtFolFQqRTAYRNd1hoeHMQyDgoICuru7AVi1ahW7du06qXADKZ6vvfZa7rzzzuy2oVBoUq+9ybgQzYkVCoVCcflxwYo3v9/Pfffdx8jICHV1dRiGwRtvvJG14ViyZAnd3d3ZN//CwkJmzZrFypUrZzzi5HK5qK2tpa2t7bjn4vE4e/bs4cMPP2RwcJClS5fy/vvvMzo6yo033sjVV189zgR2MkEwa9Ys/H4/+/fvJxAI8Mwzz2QHwufn508qQEOhEA899NDZvdAZpLW1lXnz5lFYWMi7774LkDXs9fl8tLS0UFtbS11dHUePHiUUCrFu3TqWL19OUVERb7/9NnV1dcyZM4cf/ehHHDlyhLKysuzxk8kkPT092LZNYWEho6OjLFiwgNtvv/2smDYrFAqFQnG+uGDFW2Fh4XFmvFu3bqWnp4df//rXzJ8/nyeeeAKPx8PDDz98ztdXV1fHG2+8gWEY2Zqpt99+m+eff55UKkU4HCYQCPDaa69RW1vLvffeOy5SdiI0TWPBggXs3buXQ4cOUVxczL333nvO5onqus7AwACFhYUnjeA5dYKOIDIMg66uLkCKsby8PGKxWLYxwqk1XLBgAQ8++CB5eXnceOONdHV1YZomTzzxBLqu8+lPf5qbbrppSqH1iU98IvvvDRs2sHv37qx4syyLzs5O7rzzTmpqanjqqafYsGEDW7ZsuWCaLBQKhUKhOF0uWPE2GW63mwceeIDR0VE8Hg933303eXl552UMVF1dXXZEUkVFBU899RR79+5l7ty53HzzzdTW1gIwMjJySnYaDsuWLWPXrl3k5eXxhS98YcaEWzKZRNd18vLycLvd9Pb2kkqlWLlyJTt27KCsrIz8/Hw6OzvRdZ38/PysSIpEIvT29uJyudA0Lft72Lx5M/Pnz6etrY2jR4+ybNkyZs2ahWVZDA0N4XK5mD17dnb7srKy7DHnz5+Px+M5JZE1d+7c7KxYj8dDR0cHq1ev5tprr0UIwd/93d+d5bumUCgUCsX546ISb5AzdAWyAul84DjuHz58mBdeeIEjR45w6623ZpshHMaa0Z4KS5Ys4W/+5m8oLi4+I3Fq2zbxeBzDMBgaGkLTNCzLIhgMEolECAQClJaW0tHRga7rNDQ0cN9991FRUcGqVat45JFHADk1YNOmTbzyyiu0tLQghCA/P59vf/vbFBQUMDQ0RHFxMaFQKOttNnfu3HERMmBcanMyTqdWz+PxsHbtWt5++22Ki4vxer0nNNFVKBQKheJi5qITbxcKRUVFhEIhnn/+eQDuvvtuVq9efdaOL4QY5912OjgzOEtLS8nLy2PLli3U1tZmLVNuv/12li1bhtfrJRqNMjIyQnV1dVb0LF68mC9/+cs899xz3H///RQXF7N06VKGhoYAKCkpyUYUT9bkMdOsWrWKV199lVQqxZe//OWTerwpFAqFQnGxosTbNLEsi4GBATweD6FQCCEEX//61zl48CDFxcUsWrTofC8xi2madHV1YRgGS5cu5Utf+tK4NGRFRQWf/OQnx+0TDAYJBoPHHWv58uVcccUV2cYKTdMmnSF6vmloaGDdunWsWbNmSn88hUKhUCguBZR4yxCJRPB4POi6jmmaFBUVYZpmVrS0t7fT1NREOp2mubmZ+vp6wuEwGzZsOO1z6rpOV1cXtbW1GIZBR0dHduqB1+ulqqrqhPtblkVfXx+BQCArvNLpNK2trWzevJl169ZRVVV1xt2VF4NFhsvl4v777z/fy1AoFAqFYsa5LMWbaZoIIbBtm0gkwtDQEOFwmGg0Sl5eHj6fj6NHj2aL8N1uN2VlZXz1q1/F4/Hw61//mm3btmU95U6H4eFhhoaGmDdvHgcOHEDTND772c+ydu1aenp6ePLJJ7MTGCZbf39/P7FYjCuuuILe3l5aWlqwLAvTNLn77ru59tprz/AuKRQKhUKhuBC5LMSbZVnZmi7LsrJF+0IIGhoauP7661m3bh0ejydrZdHW1kZJSQn/+q//isfj4b777ssW4t911100Njby7//+76cs3kzTpL29nZKSEv78z/+c2tpann76aVatWpVNvc6aNYtFixbx9ttvHyfeEokEnZ2drFq1imuuuSZrp6LrOvv27ePYsWNs3LjxLNw1hUKhUCgUFyKXnHjr6+sjHo+PG7pu2zZ1dXWsW7eOxsZGZs2alU0FOpMNxuKM4gL4xje+gaZpx6Ue6+vrT6mb0Zmpads2mzZt4uabb86KwXvvvfe47evq6kin0+Me03Wdzs5OvvSlLx3XHOHxeFi+fDnLly+f9poUCoVCoVBcfFzQ4m1wcJDR0VFcLld2yLxpmjQ2Nk5ax1VVVcVXvvIVAPLy8rKTGQoLC7NC6VSZym+sqKho3CB7B9M06enpwe/3EwgE8Pl8WJZFS0sLGzdu5IYbbpiWfUg4HB4nDkdGRujt7c0Om1coFAqFQnF5csGKN5fLRVlZGZ///OeprKyks7OTwcFBWltb+eCDD8Z5vKVSKYQQhEIh1q9ff07W5/P5CAaDpNNpvF5v9vHu7m6amprQNI3m5mZs28a2bZYuXcqdd9457eL/cDicFawjIyMYhsGDDz44rVmqCoVCoVAoLl0uWPFWU1PDt771rezPjj1FNBpl7969JJPJ7NSBzs7OcanQc0VFRQWDg4NZ8WaaJrqu87nPfY6ysjI6Ozv5+7//e2zb5pZbbjml9RUWFuLxeDAMg8HBQf7kT/5EpUQVCoVCoVBw0U3oDgaDXHfddfT19WEYBseOHWPBggXU1NSc87VUVVWRSCQAGf1rbm5mzZo12SkC1dXVrF27lqVLl57y+oQQVFdX09PTQ1FREUuWLDnr61coFAqFQnHxcdGJN4AVK1ZgGEbWz+xP//RPz8sopKqqKpLJZLaD9LOf/Sxbt24dt82dd97JF7/4xdM6fkNDAy6Xi3vuuWfSxgqFQqFQKBSXHxelIqisrKSxsRGfz8ett9563kxknfFQnZ2dXH311WzevPm4bc5koPwNN9zATTfdRCAQOJNlKhQKhUKhuIS4KMUbSHuN/Pz88+r+X1RURCKRoL6+nttuu+2sH1/N51QoFAqFQjGRi1a8VVRUnO8lUFFRwRe/+EXWrl17RhE2hUKhUCgUiuly0Yq3CwGPx8OmTZvO9zIUCoVCoVBcRlyUDQsKhUKhUCgUlytKvCkUCoVCoVBcRCjxplAoFAqFQnERocSbQqFQKBQKxUWEEm8KhUKhUCgUFxFKvCkUCoVCoVBcRMyIeBNCzBNC/LMQYrcQIiqE6BJC/IcQYtlMnE+hUCgUCoXicmGmIm+fAq4DHgduBb4OlAHvCCFWzdA5FQqFQqFQKC55Zsqk92ng+7Zt284DQojXgGbgm8B9M3RehUKhUCgUikuaGRFvtm33T/LYiBDiEFAzE+dUKBQKhUKhuBw4Zw0LQogSYAlw4FydU6FQKBQKheJSQ4zJbM7siYR4ArgdWGrb9uEptnkIeCjz4xJg7zlZ3IVFKXBc5PIyQF335YW67ssLdd2XF5frdc+3bbvgXJxoWuJNCLEZeHkax3vdtu1rJ9n/r4H/CXzZtu0fT2thQuywbfvK6Wx7KaGu+/JCXfflxf/f3rnH6FGVYfz3FMqlF7VFEORiqxhDm4aSKoEoFFpChQhWqRGClXAVEJAoIpVS0bRBEARFEUzDRapIFbAFC5hCuchFQYGSck3phasU20Jp6VLg+Md7Pjudndmdbzvf7s7s+0u+zH5nnjkzz8555ztn5pwz7rtv4b77Ft3pu2iftweBPQro1qUTJJ2MVdymFq24OY7jOI7jONkUqryFENYBzzSbuaTJwBXAJSGEGc1u7ziO4ziO42xKywYsSPoKcA0wM4RwVhey+G3Jh1QV3Hffwn33Ldx338J99y26zXdLBixI2h/4G/AUcBrwQWJ1WwjhsdJ36jiO4ziO0wdo1SS944Ctgb2AB1LrlgHDWrRfx3Ecx3GcWtNtU4U4juM4juM4m0+pfd4k7SLpckkPSVonKUgalqHbTdJ1kpZH3XOSpksamNLdE/NIf87MyHOipMckrZe0TNJUSVuU6S+PJnwPl/RnSaslrZW0QFK7YcWS+kmaImlp9POEpCNy9n2ipGcktUl6No7u7RZa4HtpzvmemKHtSd+TJN0Uy9k7cf8XSBqc0g2RNFPSG9H3fEmjMvLbRtLPJL0a83sodj1I6wqXi1bQAt9Z5zpIGp3S9XrfkgZLujhes96KPg7Iya8q8V2279rEt6TxkmZJWhw1iyX9RtIOGfnVJr6b9F2n+B4j6Q5JL8djfE3SPEn7ZuRX6nU/lxBCaR/gAOA/wDzgTiAAw1KagcBzwBLgGOwF9mcD7wA3prT3AE8A+6Q+O6Z0E4D3sc6CBwLfBdYDF5bpbzN9bwe8jI3a/TpwGLAAWAPskdLOANqAs6Kfq7B+g4emdCfG9BlRNz1+P6WivpcCd2Sc7yG9zPfDwGzgaGAscCawOqb3ixoB9wMvAUcBXwTuxSau3CWV3+/j9icC44GbYzyM7kq5qJDvgA1qSp/vARX0PQxYCcwHboreDsjJryrxXbbvpdQnvv8E3A4cGzUnYNe5F4BBqfzqFN/N+K5TfI8HfgUcGTVfw6ZQexfYO5FX6df93OMu+Z/QL/H3CWT/mB8c0w9Opf8UeC95YrHK298L7PcxbILgZNq0+I/dsRkPLfQ9NfrbPZE2EKv8zE6k7RAL8o9T298FLEx83xJ4Hbgupbs6FpT+VfId05cCszrZZ2/wvX1G2jej/3Hx+5fj9wMTmg9jP3S/TKTtGXXHpjw+C8xttlxUxXdMD8D0TvZZFd9KrDuInEpMxeK7NN9xfZ3iO0uzf9Qcl0irW3wX8h3TaxPfOdsNjsd+eSKt1Ot+R59SH5uGED7oXMVWcflWKn019hhXzexT0q7AaGBWatX1QH/gkGby6woFfe8DPB8SrwYLIazFaulfktQYPDIB+x+l/cwCRkkaHr/vC2yfobseu9v1haZMdIGSfRelN/hekZH8SFzuHJeHA6+EEBYktnsTuBULcBK6DcCNCd17wB+BCZK2jslFy0XLKNl3USrhO8SrbwGqFN9l+i5KVXwXiQWoWXw34bsolfCdw1qs8rYhkVb2dT+XbnsxfYL5wPPAhZJGSBokaRzwHeDK+MOeZC9Jb0raIGmhpONT60fG5SbvQQ0hLMHe+DCiBR66wvvYncA0bcC2wKfi95ExLf3+10VxOSKhg/bvf03repqivhscJus/1ybp4Yz+ML3V99i4fDouR5L9bt5FwG6SBiV0S4JNhJ3WbQXsntAVKRfdTVd9Nzglnut1ku6WtF9qfVV8F6Xq8d1V3w3qEt9FNXWL72Y1tYrv2D+vv6TdsMeoADMTkrKv+7l0e+UthLAeaz31ww50DXaL9DZsTrgk92HPnw8HJmGVvpmSpiY0Q+NyVcbuViXW9zTPAp+WtF0jQVI/YO/4dWhiuTqjRbsyQwftfad1PU1R32Ctk9Ox1tjRWL/FWyR9I6Hpdb4l7Qz8BJgfQng0cRxZZbJxnEMK6potF93GZvoGa12fij1uOwm7s3K3Nu3wXhXfRalsfG+mb6hXfKc1g4HLsB/7vyRW1S2+05o831DP+J6N3YxYBhyB9ct7KrG+7Ot+Lq2a5y0XSdtgtwp3ACYDy7Ef8mlY36hTGtoQwrTU5nMk3QKcK+myEMLbbHzMmnX7vqlHsC3mSuAM4HeSzsDuCp4LNG4LNx5BimJeOvLdmyjqmxDC6ckN47l+GLiAjbfVe5Xv2JKag5XdY5OrKH4ey9R1CyX4JoQwOfH1fklzsFbrdDY+HquK78JZUMH4LsF33eI7qdkSuAF7xPb5+Pjr/6upV3wnNR35rmt8nw1cCOwKfBu4TdJBiUpet53vnnhsejw2SvHQEMKsEMJ9IYSLge8BJ0vas5PtbwC2ARpDbzuqqX4ksb5HCSG8gLU2x2C3iF/B+ndcGiWvxuVKYIik9EkcklifXKZ9D02t71Ga8J217fvY6KZdJO0Uk3uN79gQmQt8EpgQQngpsXol2WWycR5XFdQlz3eRctFySvLdjhDCGuCvwOdS+VXBd1EqF98l+W5HxeO7oekHXIfdXZoYQliYktQtvhuazny3ow7xHUJ4IYTwSAjhZqw//etYZbRB2df9XHqi8jYKWBVCWJxK/2dc7tHJ9umWWePZ+MhNRDbf2ADsFV29ghDCTVgrZQQ2+nIMMAh4MYSwPMoWYW+nSPcFazz7fyqhg5TvDF2PU9B3HoXON93sW1J/bGqEvbGGyJMpySLaHyPYcS6Pd40buuGSBmTo3mVjX5Ci5aKllOg7dxds2iKtiu+iVCq+S/Sdu4u4rFp8N7gSmwLpyBDCXRnr6xbfDTrznbsLahLfIYR3gYVs2j+t7Ot+hwfQquG3eVNHnB/Td0+lnxTT9+sk3znYo7eBibTHgQUp3VS6aaqQIr5ztB/HathTEmmNodM/SmnnA08mvvcHVgDXpHQzgf8CW1XJd45uS+BRYFlv8o01emZjfXbG52gmxv/H2ETah+IxJoeWj466Y1K+nwZubbZcVMV3zrYfwrpR3JtIq4TvlL7IVCEd+qlKOS/qO0df2fiOukuwbh+TO9DUKr6L+s7ZrhbxndhuADb9zZ2JtFKv+x19Su/zJmlS/HNMXB4iaQWwIoRwL3AtNonuPEkz4sn8LHAe8C/iu1DjqJRzsInrlmJzpRyDDV44J2w6KvWH2LPnq7DHqnthlbdfhBBeK9tjFp35jjX7i7AJ+97CaudTsBr4JY18QgivS7oUmCJpDfBvrIUzjsRQ4xDCBknnAVdIehkr7OOA44DTg7UKWk5ZviUdhfmbB7wIfAzrUzAGm+wQ6DW+f41N0jgDWCtpn8S6l4Ldbp8LPATMkvR97Hb5FKzleVHCz+OSbgQui/+rJVi/z+HY4+aGrlC5aDGl+ZZ0FvAZbMLmV4BPYJN07kg1fSPpEGwOw0aXjrGSPgqsDSHcDpWL79J81y2+Jf0A+x27Gng+pVkR4pOlusV3Ud91i+9Yt1iJNTbewPycBuyE9d1vUOp1v0NaUIsNOZ97EpoRWE33RWxG4eeAi0nMtI3dirwdm725DXgbm9H4qJz9fhV7G0MbViGcBmzRHTX3Ir6xWvVt2OS0bcBi7Fn5gIy8tsAqn8uidiEwKWe/34r/vzZsNO6p3eW5TN/YfHB3R90G4E3swj2ht/nGGhN5vs9P6IZiF7mV2N3iu4A9M/LbFvg58BrW+vsH2XdtCpeL3u4be9PGA9iFcAPWMp1LYrbyCvrO0y3tqp+KlPNOfdctvrEJ5PM016byq018F/VNzeIbazw8GH2sx37H/gCMysiv1Ot+3sdfTO84juM4jlMhemLAguM4juM4jtNFvPLmOI7jOI5TIbzy5jiO4ziOUyG88uY4juM4jlMhvPLmOI7jOI5TIbzy5jiO4ziOUyG88uY4juM4jlMhvPLmOI7jOI5TIf4HLYr+Ep0nfq4AAAAASUVORK5CYII=\n", | |
| "text/plain": [ | |
| "<Figure size 720x432 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": {}, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "if plot_global:\n", | |
| " l = ['historical:1850-1860','historical:1861-2005','rcp45:2006-2099','rcp45:2100-2300','rcp85:2006-2099','rcp85:2100-2300']\n", | |
| " sctypelist = l; markerlist = l.copy(); colorlist = l.copy(); tnamelist = l.copy(); alphalist = l.copy()\n", | |
| " for idx,sctype in enumerate(sctypelist):\n", | |
| " markerlist[idx] = '-'; alphalist[idx] = 0.6\n", | |
| " if 'hist' in sctype:\n", | |
| " colorlist[idx] = 'black'; alphalist[idx] = 0.45; tnamelist[idx] = 'Historical'\n", | |
| " if '1850-1860' in sctype:\n", | |
| " markerlist[idx] = '-.'\n", | |
| " if 'rcp45' in sctype:\n", | |
| " colorlist[idx] = '#99CCFF'; tnamelist[idx] = 'RCP4.5' \n", | |
| " if 'rcp85' in sctype:\n", | |
| " colorlist[idx] = '#CC3333'; tnamelist[idx] = 'RCP8.5' \n", | |
| " if '2100-2300' in sctype:\n", | |
| " markerlist[idx] = '-.'\n", | |
| " \n", | |
| " # find the climatology dataset: \n", | |
| " for ds in ds_master:\n", | |
| " if 'historical:1861-2005' in ds[var].attrs['sctype']:\n", | |
| " hvar = compute_global_mean(ds).compute()\n", | |
| " \n", | |
| " hmodels = hvar.names.values.tolist()\n", | |
| "\n", | |
| " # for each model, compute the time mean from the climatology interval: 1986-2005\n", | |
| " tgm = hvar.sel(time=slice(1986,2005)).mean('time').load()\n", | |
| " \n", | |
| " # calculate the model mean tgm, for use in models which do not have a climatological reference run\n", | |
| " tgm0 = tgm[var].mean('model')\n", | |
| "\n", | |
| " plt.figure(figsize=(10,6))\n", | |
| " rcParams.update({'font.size': 16})\n", | |
| " \n", | |
| " for idx,sctype in enumerate(sctypelist):\n", | |
| " tname = tnamelist[idx];marker=markerlist[idx];color=colorlist[idx];alpha=alphalist[idx]\n", | |
| " data_exists = False\n", | |
| " for ds in ds_master:\n", | |
| " if sctype in ds[var].attrs['sctype']:\n", | |
| " tvar = compute_global_mean(ds).load()\n", | |
| " data_exists = True\n", | |
| " if data_exists:\n", | |
| " year = tvar.time.values\n", | |
| "\n", | |
| " #find the climatology for each models in this scenario:century\n", | |
| " tclimo = 0*tvar[var].mean('time')\n", | |
| " for idx,model in enumerate(tvar.names.values):\n", | |
| " if model in hmodels:\n", | |
| " hidx = hmodels.index(model)\n", | |
| " tclimo[idx] = tgm[var][hidx]\n", | |
| " else:\n", | |
| " #print('using model mean climo data:',model)\n", | |
| " tclimo[idx] = tgm0\n", | |
| "\n", | |
| " num_models = tvar[var].shape[0]\n", | |
| " range5 = 1.64*(tvar - tclimo)[var].std('model') # use std = 1.64 to give 95% and 5% of values\n", | |
| " tvar_mean = (tvar - tclimo)[var].mean('model')\n", | |
| " tvar_95 = tvar_mean + range5\n", | |
| " tvar_05 = tvar_mean - range5\n", | |
| " label = tname+' ('+str(num_models)+' models)'\n", | |
| " plt.plot(year, tvar_mean, marker, color=color, label=label)\n", | |
| " plt.fill_between(year, tvar_05, tvar_95, color=color, alpha=alpha)\n", | |
| "\n", | |
| " plt.plot((1861, 1861), (-2, 6.2), 'k-', linewidth=1.5, alpha=0.75)\n", | |
| " plt.plot((2006, 2006), (-2, 6.2), 'k-', linewidth=1.5, alpha=0.75)\n", | |
| " plt.plot((2100, 2100), (-2, 12), 'k-', linewidth=1.5, alpha=0.75)\n", | |
| " plt.plot((2200, 2200), (-2, 12), 'k-', linewidth=1.5, alpha=0.75)\n", | |
| " plt.ylim(-2,12) \n", | |
| " plt.xlim(1850,2300) \n", | |
| " vtitle = r'CMIP5 global surface air temperature change $^\\degree C$'\n", | |
| " if var == 'ts': \n", | |
| " vtitle = r'CMIP5 global surface temperature change $^\\degree C$'\n", | |
| " plt.title(vtitle,fontsize=16)\n", | |
| " plt.legend(loc='upper left',fontsize='small')\n", | |
| " figfile = 'global_' + var + '.png'\n", | |
| " #plt.savefig(figfile)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": null, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [] | |
| } | |
| ], | |
| "metadata": { | |
| "gist_info": { | |
| "gist_id": null, | |
| "gist_url": null | |
| }, | |
| "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.6.5" | |
| } | |
| }, | |
| "nbformat": 4, | |
| "nbformat_minor": 2 | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment