"
]
},
{
"cell_type": "markdown",
"metadata": {
"lines_to_next_cell": 0
},
"source": [
"# Analyzing precip and temperature data\n",
"\n",
"## Learning objectives\n",
"\n",
"1. Document your project data by capturing important information\n",
" as metadata in a json file\n",
"\n",
"1. Use the groupby method for a pandas dataframe to group observations\n",
" by season\n",
"\n",
"1. Find the probability density function (histogram) of temperature and\n",
" precipitation data and determine possible statistical models (normal distribution,\n",
" exponential distribution) that best describe the variation in the data\n",
"\n",
"1. Find best fit parameters for statistical models that describe the precipitation\n",
" and temperature data, assuming that temperature and precipitation are\n",
" independent of each other.\n",
"\n",
"1. Use groupby to find the median temperature for each decade and each season\n",
" and determine whether the seasonal temperatures have increased at YVR\n",
" since the 1930's\n",
"\n",
"## Introduction\n",
"\n",
"We have about 30,000 days worth of YVR precipitation and temperature\n",
"data generated by the 10-pandas1, 10-pandas2, and 10-pandas3 notebooks, these\n",
"gave you a new csv file: data/processed/weather_YVR.csv. We will use\n",
"that data set to estimate statistical distributions that could produce realistic\n",
"simulated rainfall and temperature inputs to a model. It makes use of several\n",
"new python modules/methods:\n",
"\n",
"1. The **context** module: We will write a new module called context.py that\n",
" will be used to locate folders and libraries associated with a project, such as\n",
" `weather_YVR.csv`.\n",
"1. **The json module**. This reads and writes files written\n",
" in \"javascript object notation\" i.e.\n",
" [json](https://en.wikipedia.org/wiki/JSON). We will use it to\n",
" save \"metadata\" or data about our data.\n",
"\n",
"1. **pandas.DataFrame.groupby**\n",
" We will use groupby to split the dataframe into 4 seasons, each corresponding\n",
" to 3 months of the year: winter (month initials djf), spring (mam),\n",
" summer (jja), and fall (son). This is the tip of the iceberg of what\n",
" you can do with pandas split/apply/combine functions. For more examples\n",
" see [the Python Data Science Handbook](https://github.com/jakevdp/PythonDataScienceHandbook/blob/master/notebooks/03.08-Aggregation-and-Grouping.ipynb)\n",
"\n",
"1. The [scipy.stats](https://docs.scipy.org/doc/scipy/reference/stats.html) module. As you will see\n",
" below, the temperature data is best fit with a a Gaussian (normal) distribution\n",
" using the norm.fit:\n",
" - [scipy.stats.norm.fit](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.norm.html#scipy.stats.norm)\n",
"\n",
" Precipitation, on the other hand, follows an exponential distribution, which\n",
" can be fit with:\n",
" - [scipy.stats.expon.fit](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.expon.html#scipy.stats.expon)"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"in context.py, found /Users/phil/repos/eosc213_students/notebooks/pandas/data/processed\n"
]
}
],
"source": [
"import json\n",
"\n",
"import context\n",
"import pandas as pd\n",
"import scipy.stats as st\n",
"from matplotlib import pyplot as plt\n",
"import pprint\n",
"import numpy as np"
]
},
{
"cell_type": "markdown",
"metadata": {
"lines_to_next_cell": 0
},
"source": [
"## The context module (repeated from 10-pandas1)\n",
"\n",
"As your project grows more complicated, it's good to have a central\n",
"module that keeps track of important files and sets your scripts\n",
"up so that they can import functions and classes from you modules.\n",
"If you were planning to distribute your project using conda, then\n",
"you would need to write an installation script, which is a fair\n",
"amount of work. At this stage, it's easier and more flexible to\n",
"store that information in a file that travels along with your notebook.\n",
"By convention, this file is called:\n",
"[context.py](https://github.com/phaustin/eosc213_students/blob/master/notebooks/pandas/context.py)\n",
"\n",
"Clicking on that link shows you the source -- the code is executed when\n",
"you do `import context` and does the following:\n",
"\n",
"1. Defines Path objects that give the location of the raw and processed data folders\n",
"\n",
"1. Puts the notebooks/pandas folder on python's list of places to look for\n",
" library modules (sys.path). We will use this when we start extracting\n",
" code from notebooks into libraries"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"/Users/phil/repos/eosc213_students/notebooks/pandas/data/processed\n"
]
}
],
"source": [
"#\n",
"# Here is the path to the processed csv file\n",
"#\n",
"print(context.processed_dir)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Reading the processed data\n",
"\n",
"The next cell uses [\"globbing\"](https://en.wikipedia.org/wiki/Glob_(programming)) to\n",
"find every csv file in all folders below `data/processed`. (See\n",
"[this optional tutorial](https://realpython.com/python-pathlib/) for more information)."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"yvr files: \n",
"[PosixPath('/Users/phil/repos/eosc213_students/notebooks/pandas/data/processed/weather_daily_YVR_1938-2017_all.csv'),\n",
" PosixPath('/Users/phil/repos/eosc213_students/notebooks/pandas/data/processed/weather_daily_YVR_1938-2019_all.csv'),\n",
" PosixPath('/Users/phil/repos/eosc213_students/notebooks/pandas/data/processed/weather_YVR.csv')]\n"
]
}
],
"source": [
"yvr_files = list(context.processed_dir.glob(\"**/*YVR*csv\"))\n",
"print(f\"yvr files: \\n{pprint.pformat(yvr_files)}\")"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"yvr_file = context.processed_dir / \"weather_YVR.csv\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Metadata is data about data\n",
"\n",
"Since we care about this dataset, it's a good idea to save\n",
"details for future reference. We can't put this kind of information\n",
"into a csv file, but there are many file formats that can store\n",
"this type of descriptive metadata -- one of the most popular is\n",
"[json](https://en.wikipedia.org/wiki/JSON). In the next cell we\n",
"write the data into a nested dictionary called meta_dict, and\n",
"dump that information into a new json file called metadata.json"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"meta_data = dict()\n",
"\n",
"file_info = dict()\n",
"file_info[\"written_by\"] = \"10-pandas3-process-daily-data.ipynb\"\n",
"file_info[\"contact\"] = \"paustin@eoas.ubc.ca\"\n",
"file_info[\n",
" \"description\"\n",
"] = \"EC canada YVR data downloaded by 10-pandas2-download-daily-data.ipynb\"\n",
"\n",
"meta_data[\"weather_daily_YVR_1938-2017_all.csv\"] = file_info\n",
"\n",
"file_info = dict()\n",
"file_info[\"written_by\"] = \"10-pandas3-process-daily-data.ipynb\"\n",
"file_info[\"contact\"] = \"paustin@eoas.ubc.ca\"\n",
"file_info[\n",
" \"description\"\n",
"] = \"Reduced YVR data subsetted by 10-pandas3-process-daily-data.ipynb\"\n",
"\n",
"meta_data[\"weather_YVR.csv\"] = file_info\n",
"\n",
"history_file = context.processed_dir / \"history_metadata.json\"\n",
"with open(history_file, \"w\") as f:\n",
" json.dump(meta_data, f, indent=4)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Take a minute and look at history_metadata.json with jupyter to see how the metadata is written out."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Read weather_YVR.csv into a dataframe\n",
"\n",
"Here's the dataframe produced by 10-pandas3\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"there are 29632 days in the dataframe\n"
]
},
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
Date
\n",
"
Year
\n",
"
Month
\n",
"
Day
\n",
"
T_mean (C)
\n",
"
T_high (C)
\n",
"
T_low (C)
\n",
"
Rain (mm)
\n",
"
Snow (cm)
\n",
"
Total Precip (mm)
\n",
"
\n",
" \n",
" \n",
"
\n",
"
0
\n",
"
1938-01-01
\n",
"
1938
\n",
"
1
\n",
"
1
\n",
"
4.4
\n",
"
9.4
\n",
"
-0.6
\n",
"
NaN
\n",
"
NaN
\n",
"
0.3
\n",
"
\n",
"
\n",
"
1
\n",
"
1938-01-02
\n",
"
1938
\n",
"
1
\n",
"
2
\n",
"
4.5
\n",
"
7.2
\n",
"
1.7
\n",
"
NaN
\n",
"
NaN
\n",
"
0.5
\n",
"
\n",
"
\n",
"
2
\n",
"
1938-01-03
\n",
"
1938
\n",
"
1
\n",
"
3
\n",
"
1.7
\n",
"
7.2
\n",
"
-3.9
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
\n",
"
\n",
"
3
\n",
"
1938-01-04
\n",
"
1938
\n",
"
1
\n",
"
4
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
\n",
"
\n",
"
4
\n",
"
1938-01-05
\n",
"
1938
\n",
"
1
\n",
"
5
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" Date Year Month Day T_mean (C) T_high (C) T_low (C) Rain (mm) \\\n",
"0 1938-01-01 1938 1 1 4.4 9.4 -0.6 NaN \n",
"1 1938-01-02 1938 1 2 4.5 7.2 1.7 NaN \n",
"2 1938-01-03 1938 1 3 1.7 7.2 -3.9 0.0 \n",
"3 1938-01-04 1938 1 4 2.2 7.2 -2.8 0.0 \n",
"4 1938-01-05 1938 1 5 2.2 7.2 -2.8 0.0 \n",
"\n",
" Snow (cm) Total Precip (mm) \n",
"0 NaN 0.3 \n",
"1 NaN 0.5 \n",
"2 0.0 0.0 \n",
"3 0.0 0.0 \n",
"4 0.0 0.0 "
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"yvr_df = pd.read_csv(yvr_file)\n",
"print(f\"there are {len(yvr_df)} days in the dataframe\")\n",
"yvr_df.head()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Handling missing data\n",
"\n",
"What should we do about missing vlaues like the Rain and Snow totals for Jan./Feb. 1938?\n",
"This is a judgement call -- we could through away the whole row, or do something\n",
"more complicated like masking the missing values in calculations where they are referenced.\n",
"For this notebook we'll fill the NaNs with zeros -- not a great idea since 0 is an\n",
"actual temperature, but we can see if that makes a difference later"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"yvr_df_zeros=yvr_df.fillna(0.)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Use apply to tag the 29,190 days with their season\n",
"\n",
"We know that the seasons are quite different, and if we are interested in\n",
"generating daily data we have to take that into account. In the next cell\n",
"we set up a lookup table (dictionary) called season_table that maps the\n",
"month number to one of four seasons. With this table we can write\n",
"a function called find_season that takes a row of the dataframe and\n",
"returns the season for that row"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"season_table = dict()\n",
"for i in range(12):\n",
" m = i + 1\n",
" if m in [12, 1, 2]:\n",
" # winter\n",
" season_table[m] = \"djf\"\n",
" elif m in [3, 4, 5]:\n",
" # spring\n",
" season_table[m] = \"mam\"\n",
" elif m in [6, 7, 8]:\n",
" # summer\n",
" season_table[m] = \"jja\"\n",
" else:\n",
" # fall\n",
" season_table[m] = \"son\"\n",
"\n",
"\n",
"def find_season(row, season_table):\n",
" month = row[\"Month\"]\n",
" return season_table[month]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### creating a new column with apply\n",
"\n",
"One of the most useful pandas dataframe methods is [apply](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.apply.html). This method has lots of bells and whistles, all we need here is to apply the `find_season` function to every\n",
"row of the dataframe, which means that we want to return a new column with the same number of rows\n",
"as the other columns. We do that by tell apply that we want to return `axis=1` (rows are axis=0)"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"season = yvr_df.apply(find_season, args=(season_table,), axis=1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now add that column to the dataframe"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
Date
\n",
"
Year
\n",
"
Month
\n",
"
Day
\n",
"
T_mean (C)
\n",
"
T_high (C)
\n",
"
T_low (C)
\n",
"
Rain (mm)
\n",
"
Snow (cm)
\n",
"
Total Precip (mm)
\n",
"
season
\n",
"
\n",
" \n",
" \n",
"
\n",
"
0
\n",
"
1938-01-01
\n",
"
1938
\n",
"
1
\n",
"
1
\n",
"
4.4
\n",
"
9.4
\n",
"
-0.6
\n",
"
NaN
\n",
"
NaN
\n",
"
0.3
\n",
"
djf
\n",
"
\n",
"
\n",
"
1
\n",
"
1938-01-02
\n",
"
1938
\n",
"
1
\n",
"
2
\n",
"
4.5
\n",
"
7.2
\n",
"
1.7
\n",
"
NaN
\n",
"
NaN
\n",
"
0.5
\n",
"
djf
\n",
"
\n",
"
\n",
"
2
\n",
"
1938-01-03
\n",
"
1938
\n",
"
1
\n",
"
3
\n",
"
1.7
\n",
"
7.2
\n",
"
-3.9
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
"
\n",
"
3
\n",
"
1938-01-04
\n",
"
1938
\n",
"
1
\n",
"
4
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
"
\n",
"
4
\n",
"
1938-01-05
\n",
"
1938
\n",
"
1
\n",
"
5
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" Date Year Month Day T_mean (C) T_high (C) T_low (C) Rain (mm) \\\n",
"0 1938-01-01 1938 1 1 4.4 9.4 -0.6 NaN \n",
"1 1938-01-02 1938 1 2 4.5 7.2 1.7 NaN \n",
"2 1938-01-03 1938 1 3 1.7 7.2 -3.9 0.0 \n",
"3 1938-01-04 1938 1 4 2.2 7.2 -2.8 0.0 \n",
"4 1938-01-05 1938 1 5 2.2 7.2 -2.8 0.0 \n",
"\n",
" Snow (cm) Total Precip (mm) season \n",
"0 NaN 0.3 djf \n",
"1 NaN 0.5 djf \n",
"2 0.0 0.0 djf \n",
"3 0.0 0.0 djf \n",
"4 0.0 0.0 djf "
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"yvr_df[\"season\"] = season\n",
"yvr_df.head()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Grouping the data by season\n",
"\n",
"The next cell uses the groupby method to sort all of the\n",
"rows into 4 seasons. The resulting group object (`seasons`) has\n",
"4 dataframes inside it, keyed by the season marker djf, mam, etc."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"seasons = yvr_df.groupby(\"season\")"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n"
]
},
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
Date
\n",
"
Year
\n",
"
Month
\n",
"
Day
\n",
"
T_mean (C)
\n",
"
T_high (C)
\n",
"
T_low (C)
\n",
"
Rain (mm)
\n",
"
Snow (cm)
\n",
"
Total Precip (mm)
\n",
"
season
\n",
"
\n",
" \n",
" \n",
"
\n",
"
0
\n",
"
1938-01-01
\n",
"
1938
\n",
"
1
\n",
"
1
\n",
"
4.4
\n",
"
9.4
\n",
"
-0.6
\n",
"
NaN
\n",
"
NaN
\n",
"
0.3
\n",
"
djf
\n",
"
\n",
"
\n",
"
1
\n",
"
1938-01-02
\n",
"
1938
\n",
"
1
\n",
"
2
\n",
"
4.5
\n",
"
7.2
\n",
"
1.7
\n",
"
NaN
\n",
"
NaN
\n",
"
0.5
\n",
"
djf
\n",
"
\n",
"
\n",
"
2
\n",
"
1938-01-03
\n",
"
1938
\n",
"
1
\n",
"
3
\n",
"
1.7
\n",
"
7.2
\n",
"
-3.9
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
"
\n",
"
3
\n",
"
1938-01-04
\n",
"
1938
\n",
"
1
\n",
"
4
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
"
\n",
"
4
\n",
"
1938-01-05
\n",
"
1938
\n",
"
1
\n",
"
5
\n",
"
2.2
\n",
"
7.2
\n",
"
-2.8
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
djf
\n",
"
\n",
"
\n",
"
59
\n",
"
1938-03-01
\n",
"
1938
\n",
"
3
\n",
"
1
\n",
"
5.0
\n",
"
9.4
\n",
"
0.6
\n",
"
3.6
\n",
"
0.0
\n",
"
3.6
\n",
"
mam
\n",
"
\n",
"
\n",
"
60
\n",
"
1938-03-02
\n",
"
1938
\n",
"
3
\n",
"
2
\n",
"
4.7
\n",
"
8.3
\n",
"
1.1
\n",
"
0.5
\n",
"
0.0
\n",
"
0.5
\n",
"
mam
\n",
"
\n",
"
\n",
"
61
\n",
"
1938-03-03
\n",
"
1938
\n",
"
3
\n",
"
3
\n",
"
8.1
\n",
"
12.8
\n",
"
3.3
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
mam
\n",
"
\n",
"
\n",
"
62
\n",
"
1938-03-04
\n",
"
1938
\n",
"
3
\n",
"
4
\n",
"
8.6
\n",
"
13.3
\n",
"
3.9
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
mam
\n",
"
\n",
"
\n",
"
63
\n",
"
1938-03-05
\n",
"
1938
\n",
"
3
\n",
"
5
\n",
"
6.2
\n",
"
10.6
\n",
"
1.7
\n",
"
1.3
\n",
"
0.0
\n",
"
1.3
\n",
"
mam
\n",
"
\n",
"
\n",
"
151
\n",
"
1938-06-01
\n",
"
1938
\n",
"
6
\n",
"
1
\n",
"
11.7
\n",
"
18.9
\n",
"
4.4
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
jja
\n",
"
\n",
"
\n",
"
152
\n",
"
1938-06-02
\n",
"
1938
\n",
"
6
\n",
"
2
\n",
"
13.6
\n",
"
20.0
\n",
"
7.2
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
jja
\n",
"
\n",
"
\n",
"
153
\n",
"
1938-06-03
\n",
"
1938
\n",
"
6
\n",
"
3
\n",
"
14.7
\n",
"
21.1
\n",
"
8.3
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
jja
\n",
"
\n",
"
\n",
"
154
\n",
"
1938-06-04
\n",
"
1938
\n",
"
6
\n",
"
4
\n",
"
16.7
\n",
"
25.0
\n",
"
8.3
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
jja
\n",
"
\n",
"
\n",
"
155
\n",
"
1938-06-05
\n",
"
1938
\n",
"
6
\n",
"
5
\n",
"
16.4
\n",
"
22.2
\n",
"
10.6
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
jja
\n",
"
\n",
"
\n",
"
243
\n",
"
1938-09-01
\n",
"
1938
\n",
"
9
\n",
"
1
\n",
"
16.2
\n",
"
20.6
\n",
"
11.7
\n",
"
1.3
\n",
"
0.0
\n",
"
1.3
\n",
"
son
\n",
"
\n",
"
\n",
"
244
\n",
"
1938-09-02
\n",
"
1938
\n",
"
9
\n",
"
2
\n",
"
15.3
\n",
"
18.9
\n",
"
11.7
\n",
"
2.3
\n",
"
0.0
\n",
"
2.3
\n",
"
son
\n",
"
\n",
"
\n",
"
245
\n",
"
1938-09-03
\n",
"
1938
\n",
"
9
\n",
"
3
\n",
"
15.0
\n",
"
19.4
\n",
"
10.6
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
son
\n",
"
\n",
"
\n",
"
246
\n",
"
1938-09-04
\n",
"
1938
\n",
"
9
\n",
"
4
\n",
"
14.8
\n",
"
17.8
\n",
"
11.7
\n",
"
3.3
\n",
"
0.0
\n",
"
3.3
\n",
"
son
\n",
"
\n",
"
\n",
"
247
\n",
"
1938-09-05
\n",
"
1938
\n",
"
9
\n",
"
5
\n",
"
15.6
\n",
"
19.4
\n",
"
11.7
\n",
"
1.3
\n",
"
0.0
\n",
"
1.3
\n",
"
son
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" Date Year Month Day T_mean (C) T_high (C) T_low (C) \\\n",
"0 1938-01-01 1938 1 1 4.4 9.4 -0.6 \n",
"1 1938-01-02 1938 1 2 4.5 7.2 1.7 \n",
"2 1938-01-03 1938 1 3 1.7 7.2 -3.9 \n",
"3 1938-01-04 1938 1 4 2.2 7.2 -2.8 \n",
"4 1938-01-05 1938 1 5 2.2 7.2 -2.8 \n",
"59 1938-03-01 1938 3 1 5.0 9.4 0.6 \n",
"60 1938-03-02 1938 3 2 4.7 8.3 1.1 \n",
"61 1938-03-03 1938 3 3 8.1 12.8 3.3 \n",
"62 1938-03-04 1938 3 4 8.6 13.3 3.9 \n",
"63 1938-03-05 1938 3 5 6.2 10.6 1.7 \n",
"151 1938-06-01 1938 6 1 11.7 18.9 4.4 \n",
"152 1938-06-02 1938 6 2 13.6 20.0 7.2 \n",
"153 1938-06-03 1938 6 3 14.7 21.1 8.3 \n",
"154 1938-06-04 1938 6 4 16.7 25.0 8.3 \n",
"155 1938-06-05 1938 6 5 16.4 22.2 10.6 \n",
"243 1938-09-01 1938 9 1 16.2 20.6 11.7 \n",
"244 1938-09-02 1938 9 2 15.3 18.9 11.7 \n",
"245 1938-09-03 1938 9 3 15.0 19.4 10.6 \n",
"246 1938-09-04 1938 9 4 14.8 17.8 11.7 \n",
"247 1938-09-05 1938 9 5 15.6 19.4 11.7 \n",
"\n",
" Rain (mm) Snow (cm) Total Precip (mm) season \n",
"0 NaN NaN 0.3 djf \n",
"1 NaN NaN 0.5 djf \n",
"2 0.0 0.0 0.0 djf \n",
"3 0.0 0.0 0.0 djf \n",
"4 0.0 0.0 0.0 djf \n",
"59 3.6 0.0 3.6 mam \n",
"60 0.5 0.0 0.5 mam \n",
"61 0.0 0.0 0.0 mam \n",
"62 0.0 0.0 0.0 mam \n",
"63 1.3 0.0 1.3 mam \n",
"151 0.0 0.0 0.0 jja \n",
"152 0.0 0.0 0.0 jja \n",
"153 0.0 0.0 0.0 jja \n",
"154 0.0 0.0 0.0 jja \n",
"155 0.0 0.0 0.0 jja \n",
"243 1.3 0.0 1.3 son \n",
"244 2.3 0.0 2.3 son \n",
"245 0.0 0.0 0.0 son \n",
"246 3.3 0.0 3.3 son \n",
"247 1.3 0.0 1.3 son "
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#\n",
"# this creates a new object of type DataFrameGroupBy\n",
"#\n",
"print(type(seasons))\n",
"seasons.head()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Use a dictionary comprehension put the dataframes into a dictionary\n",
"\n",
"In the next cell we create a new dictionary by iterating over the DataFrameGroupBy object. We\n",
"will use a [dictionary comprehension](https://jakevdp.github.io/WhirlwindTourOfPython/11-list-comprehensions.html),\n",
"to simplify the loop -- note that we wind up with 4 dataframes in the dictionary,\n",
"one for each season. We do this because a dictionary of dataframes is simpler to work with than\n",
"the DataFrameGroupBy object, which somewhat obscures the season keys."
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"dict_keys(['djf', 'jja', 'mam', 'son'])"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"season_dict = {key: value for key, value in seasons}\n",
"season_dict.keys()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here is how you get the fall (son) dataframe from the dictionary"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
Date
\n",
"
Year
\n",
"
Month
\n",
"
Day
\n",
"
T_mean (C)
\n",
"
T_high (C)
\n",
"
T_low (C)
\n",
"
Rain (mm)
\n",
"
Snow (cm)
\n",
"
Total Precip (mm)
\n",
"
season
\n",
"
\n",
" \n",
" \n",
"
\n",
"
243
\n",
"
1938-09-01
\n",
"
1938
\n",
"
9
\n",
"
1
\n",
"
16.2
\n",
"
20.6
\n",
"
11.7
\n",
"
1.3
\n",
"
0.0
\n",
"
1.3
\n",
"
son
\n",
"
\n",
"
\n",
"
244
\n",
"
1938-09-02
\n",
"
1938
\n",
"
9
\n",
"
2
\n",
"
15.3
\n",
"
18.9
\n",
"
11.7
\n",
"
2.3
\n",
"
0.0
\n",
"
2.3
\n",
"
son
\n",
"
\n",
"
\n",
"
245
\n",
"
1938-09-03
\n",
"
1938
\n",
"
9
\n",
"
3
\n",
"
15.0
\n",
"
19.4
\n",
"
10.6
\n",
"
0.0
\n",
"
0.0
\n",
"
0.0
\n",
"
son
\n",
"
\n",
"
\n",
"
246
\n",
"
1938-09-04
\n",
"
1938
\n",
"
9
\n",
"
4
\n",
"
14.8
\n",
"
17.8
\n",
"
11.7
\n",
"
3.3
\n",
"
0.0
\n",
"
3.3
\n",
"
son
\n",
"
\n",
"
\n",
"
247
\n",
"
1938-09-05
\n",
"
1938
\n",
"
9
\n",
"
5
\n",
"
15.6
\n",
"
19.4
\n",
"
11.7
\n",
"
1.3
\n",
"
0.0
\n",
"
1.3
\n",
"
son
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" Date Year Month Day T_mean (C) T_high (C) T_low (C) \\\n",
"243 1938-09-01 1938 9 1 16.2 20.6 11.7 \n",
"244 1938-09-02 1938 9 2 15.3 18.9 11.7 \n",
"245 1938-09-03 1938 9 3 15.0 19.4 10.6 \n",
"246 1938-09-04 1938 9 4 14.8 17.8 11.7 \n",
"247 1938-09-05 1938 9 5 15.6 19.4 11.7 \n",
"\n",
" Rain (mm) Snow (cm) Total Precip (mm) season \n",
"243 1.3 0.0 1.3 son \n",
"244 2.3 0.0 2.3 son \n",
"245 0.0 0.0 0.0 son \n",
"246 3.3 0.0 3.3 son \n",
"247 1.3 0.0 1.3 son "
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"season_dict[\"son\"].head()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Fitting the distributions\n",
"\n",
"What does it mean to \"fit\" a distribution? We want to be able to generate a series of random\n",
"numbers that resemble temperatures and precip amounts that we might see at YVR during a \n",
"particular season. We need to find two things:\n",
"\n",
"1. A probability density function (pdf) that has a similar shape to the histogram\n",
" we get from the data.\n",
" \n",
"1. The best choice of parameters for that pdf that getting it \"closest\" in some statistical\n",
" sense to the data.\n",
" \n",
"You'll see below that temperature is best fit with a Gaussin (normal) distribution, whilc\n",
"precipitation more closely follows an exponential distribution. Both the normal and\n",
"exponential distributions are functions of two parameters. For the normal distribution,\n",
"the best estimate (maximum likelihood values) of those parameters is provided by the \n",
"mean and the standard deviation of the temperature data. For the exponential distribution,\n",
"they are given by the minimum value and the mean.\n",
"\n",
"### Daily average temperature\n",
"\n",
"The next cell shows histograms temperature for each of the seasons. We use\n",
"[scipy.stats.norm.fit](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.norm.html#scipy.stats.norm) to\n",
"find the mean and standard deviation (called the `loc` and `scale` here) of the data\n",
"and then shows that fitted distribution as a red line. Specifically, for temperatures x the red\n",
"line is a plot of:\n",
"\n",
"\\begin{align*}\n",
"f(y) &= \\frac{\\exp(-y^2/2)}{\\sqrt{2\\pi}} \\\\\n",
"y &= (x - loc) / scale\n",
"\\end{align*}\n",
"\n",
"The four plots show that a normal distribution give as pretty good representation of each of the seasons."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Missing values again\n",
"\n",
"How many of our temperature measurements are missing? Here's the count for the whole dataframe"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"We are missing 24 temps and 25 precips out of 29632 records\n"
]
}
],
"source": [
"missing_temps=np.sum(np.isnan(yvr_df['T_mean (C)'].values))\n",
"missing_precip=np.sum(np.isnan(yvr_df['Total Precip (mm)'].values))\n",
"print(f\"We are missing {missing_temps} temps and {missing_precip} precips out of {len(yvr_df)} records\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"So not a showstopper, but still unprofessional/incorrect to set these to 0. Instead\n",
"we'll remove them with the [dropna](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.dropna.html) method\n",
"below"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"lines_to_next_cell": 2
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"6 nan values dropped\n",
"8 nan values dropped\n",
"5 nan values dropped\n",
"5 nan values dropped\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Users/phil/mini37/envs/e213/lib/python3.6/site-packages/numpy/lib/histograms.py:824: RuntimeWarning: invalid value encountered in greater_equal\n",
" keep = (tmp_a >= first_edge)\n",
"/Users/phil/mini37/envs/e213/lib/python3.6/site-packages/numpy/lib/histograms.py:825: RuntimeWarning: invalid value encountered in less_equal\n",
" keep &= (tmp_a <= last_edge)\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAgEAAAHxCAYAAAD0nG/lAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xd8VfX5wPHPQ0iAsEcII0DClg1CQBlGVBy4ah1Y2qJVaau22tIqba2rarWtq2qHVX91FhBFUZEhMzLDhrBHkBCQPcIIGc/vj3MvhpCQcce543m/Xvd1b849594nN+d889zvFFXFGGOMMdGnmtsBGGOMMcYdlgQYY4wxUcqSAGOMMSZKWRJgjDHGRClLAowxxpgoZUmAMcYYE6UsCTBBJSL/FZGnRGSwiGwstr2TiKwQkWMi8ks3YzTGmGhR3e0ATHRS1XSgU7FNDwFzVLW3SyEZY0zUsZoAEyraAJluB2GMMdHEkgATUCLSW0SWe6r5xwM1PdvTRCTb83gWcCnwqojkikhHF0M2JiqJSJaI/FZEVovIcRF5U0QSReRLz/X7lYg09Oz7oYjsEZEjIjJPRLoWe53/isg/PMflish8EWkmIi+JyCER2SAiVuMXIiwJMAEjInHAJ8C7QCPgQ+D7JfdT1aFAOnC/qtZR1U1BDdQY4/V94AqgI3Ad8CXwe6AJzv8Lb3+dL4EOQFNgOfB+ide5FXjEc1wesNCzXxNgIvBCIH8JU3GWBJhAGgDEAi+par6qTgQyXI7JGFO2V1T1W1XdhZOYL1bVFaqaB0wCegOo6luqesyz/XGgp4jUL/Y6k1R1maqe8hx3SlXfUdVCYLz3dYz7LAkwgdQC2KVnr1K1w61gjDHl+rbY45Ol/FxHRGJE5FkR2SoiR4Esz/NNKvM6/gvZ+MKSABNIu4GWIiLFtrV2KxhjjF/8ALgBuByoDyR7tktZB5jQZUmACaSFQAHwSxGpLiI3Aakux2SM8U1dnHb+A0A88Iy74RhfWBJgAkZVTwM3AXcAh4DbgI/djMkY47N3cJr1dgHrgEXuhmN8IWc31xoTHCIyFHhDVdu6HYsxxkQrqwkwbukGbHc7CGOMiWY2bbAJOhF5GbgeGOV2LMYYE82sOcAYY4yJUtYcYIwxxkQpSwKMMcaYKBUVfQKaNGmiycnJbodRquPHj1O7dm23wwgJ9ll8x63PYtmyZftVNSHob1xBoXwtl8fOb4d9DsH5DCp6LUdFEpCcnMzSpUvdDqNUc+bMIS0tze0wQoJ9Ft9x67MQkZCe1jmUr+Xy2PntsM8hOJ9BRa9law4wxhhjopQlAcYYY0yUsiTAGGOMiVKWBBhjjDFRypIAY4wxJkpZEmCMMcZEKUsCjDHGmCgVFfMEGFMRyWO/KHefrGeHByESY8KLXTvhy2oCjDE+EZGrRGSjiGwRkbGlPD9ERJaLSIGI3Fxsey8RWSgimSKyWkRuC27kxhhXkgArNIyJDCISA7wGXA10AW4XkS4ldvsGuAP4oMT2E8CPVbUrcBXwkog0CGzExpjigt4cUKzQuALIBjJEZLKqriu2m7fQ+E2Jw72FxmYRaQEsE5Fpqno4CKEbY86VCmxR1W0AIjIOuAE4cz2rapbnuaLiB6rqpmKPc0RkL5AA2PVsTJC40SfACg1jIkdLYGexn7OB/pV9ERFJBeKArX6KyxhTAW4kAVZomKCrSMclUyVSyjat1AuINAfeBUapalEpz48GRgMkJiYyZ86cKoTpvtzc3LCNvTxjuheUu4/3d4/kz6GiQukzcCMJCHih4dknLAqOUDoZ3BbIz6IihVRFBOtvFUbnRTbQqtjPSUBORQ8WkXrAF8AjqrqotH1U9XXgdYC+fftquK5AF8mr591RkdEBI9OAyP4cKiqUPgM3koCAFxoQPgVHKJ0MbgvkZ1GRQqoivAVZoIXReZEBdBCRFGAXMAL4QUUOFJE4YBLwjqp+GLgQjTFlcWN0wJlCw1MIjAAmV+RAKzSMCS2qWgDcD0wD1gMTVDVTRJ4UkesBRKSfiGQDtwD/FpFMz+G3AkOAO0RkpefWy4Vfw5ioFfSaAFUtEBFvoREDvOUtNIClqjpZRPrh/LNvCFwnIk94hhF5C43GInKH5yXvUNWVwf49jDEOVZ0CTCmx7dFijzNwavxKHvce8F7AAzTGlMmVGQOt0DDGGGPcZzMGGmOMMVHKkgBjjDF+0+DkUX6S8Sm18064HYqpAEsCjDHG+M2DX3/Ao7P+w+R3fkXnvdvdDseUw5IAY4wxflHvVC63rpnBolbdqHP6JJ+8O4YRK6eCVmoqGBNEtpSwMcYYvxi58kvi8/N4/PKfsq92Q178/HmenfYqA3augT+mQZ06bodoSrCaAGOMMT6LLcxn1LLPSG/Tiw1NUzhQuwGjbn2Cvw7+EdetT4e+fWH1arfDNCVYEmCMMcZn165Pp1nuQd5I/d6ZbSrVeO3i2xg54ik4cgT696f5F19Y80AIsSTAGGOMb1S5O+MTNjVuzdyUPuc8vah1D1i5EgYNotPf/gY//jHk5roQqCnJkgATdVoc3cvT016lxdG9bodiTES46Js1dN27jTf63QhS2hpxQGIiTJ3K9jvvhA8+gH79IKfCy8aYALEkwESVznu38/G7v2Hkyqn8afo/3Q7HmIhwd8Yk9sU34NOuaeffMSaGHT/+McyYARs3whtvBCU+UzZLAkzUuGjHKia8/zCK8Haf4Vy2NYOhW5a4HZYx4W3DBi7bmsG7fYaTVz2uYscMHQr9+8Pnnwc2NlMuSwJMVLhu3VzenvAYu+s14aYf/Y0/Db2HTY1b89jM16lRcNrt8IwJXy++yKnqcbzX+5rKHXfttZCRAd9+G5i4TIVYEmAi3wsv8Mpnf2VFy87cMvIv7K6XQEFMdR674qe0ObyHny7+yO0IjQlP+/bBO+/wcdehHIyvX7ljhw937r/80v9xmQqzJMBErqIi+PWvYcwYvug0kB/f+iRHa343WcnCNj35vPNg7l30IUlH7NuIMZX2z3/CqVO82feGyh/bsye0bGlNAi6zJMBEprw8GDkSXnwRfvELfnH9Q6W2Vz516V0USjX+OPM/LgRpTBg7dQpeew2uuYatTVpV/ngRpzZg+nQ4bU1ybrEkwESeI0fg6qth3Dh47jl4+WWKqsWUuuueek145eIRXLl5EZdsWxbkQI0JY++/D3v3wpgxVX+N4cPh2DFIT/dfXKZSLAkwkWXvXhgyxClU3nkHHnqo7HHLHm/2u4GtjVry+Ff/Iq4gP0iBGhPGVOGFF5wq/UsvrfrrXHYZ1KgBX3zhv9hMpVgSYCLLSy9BZqZTqPzoRxU6JD8mlscv/ykph3Zzd8akAAdoTASYNg3WrXNqAcpJss+rdm0nibAkwDWWBJjIoQoffeQUKsOGVerQ9JQ+fNnxYu5fON5mEjSmPC+8AC1awG23+f5aw4fDpk2webPvr2UqzZIAE/aSx35B8tgvuOLuf8KmTTxSreOZbd5bRTw19G5E4Q+z3gxwxMaEsdWrnRn/fvELiKvg5EDn4x0qaLUBrrAkwESMqzYtoAhhWocBVTp+V/2mvHbRLQzfOJ+BWSv9HJ0xEeLFFyE+HkaP9s/rpaRAly6WBLjEkgATMa7etIBlLS9gX51GVX6N/6TeRFaD5jwx41/EFlonwYoQkatEZKOIbBGRsaU8P0RElotIgYjcXOK5USKy2XMbFbyoTZXs3u2MCvjJT6BR1a+zcwwfDnPnOiMFTFC5kgRYoWH8rfWh3XTZu52pnS726XXyqsfxxOWjaX8wmzuXTvZTdJFLRGKA14CrgS7A7SLSpcRu3wB3AB+UOLYR8BjQH0gFHhORhoGO2fjg3/+GggJ44AH/vu6110J+vtPMYIIq6EmAFRomEK7atACAqR19SwIAZrfrx4z2qTww/38kHtvv8+tFuFRgi6puU9XTwDjgrOnjVDVLVVcDRSWOvRKYoaoHVfUQMAO4KhhBmyqaMQMGDID27f37uhdfDA0aWJOAC6q78J5nCg0AEfEWGuu8O6hqlue5MgsNz/PeQuN/gQ/bhLKrNy5gdbP27Krf1C+v9+Rlo/nqjZ9xz5JJPHXZPX55zQjVEthZ7OdsnCS9qse2LLmTiIwGRgMkJiYyZ86cKgXqttzc3LCNHaBaXh6Dliwh++ab2Vbi9xjTvaDc472/e1mfQ5fevWkwaRILRo6EapHdUh1K54IbSUDACw0In4IjlE4Gt1X1s3i0xR56797Iglt/WKHCqGKasGdBd27KWcbJ7nee2Rqsv1UYnRelDRJXfx6rqq8DrwP07dtX09LSKhxcKJkzZw7hGLt3dE3qzrUMKSjg8fyOzJp6vMReFfhXssY5Zkz3Qp7/uuTxcGPtPrx0aDZp9epB376+hh3SQulccCMJCHihAeFTcITSyeC2qn4Ws58eB8Aj9QaxbY3/TunDCRfyxzVvMO7rg2dqGLJGVj6+qgij8yIbKD5xfBKQU4lj00ocO8cvURm/67czE4BlLS8IyOvPTenjTDz0+ecRnwSEEjfqXHwtNKp6rIlQV2+cz8YmrdnWOMmvrzsn5UIAhmxf7tfXjTAZQAcRSRGROGAEUNEeldOAYSLS0NO3Z5hnmwlB/bLXsaFJG47UqhuQ1z8UX9/pb2D9AoLKjSTACg3jP3v30i97nV86BJa0tXESu+omWBJwHqpaANyPcx2uByaoaqaIPCki1wOISD8RyQZuAf4tIpmeYw8Cf8IpEzKAJ739fUxoqVZUSJ9d61maVLIPt59dey0sXQp79gT2fcwZQU8CrNAwfvXpp8Rokc9DA0slwty2fRiYtZLqhf7qaxB5VHWKqnZU1Xaq+rRn26OqOtnzOENVk1S1tqo2VtWuxY59S1Xbe27/59bvYM6v874d1Dt9giWtupa/sy+8swdOmRLY9zFnuNIF0woN4zcff0xWg+asT0gJyMvPS+lDvdMn6J2zISCvb0w46Jft9AcIeE1Ajx6QlGRNAkEU2eMwTGQ7fBhmzmRqx4t8W8nsPBa06UmBVGPI9hUBeX1jwkG/7HVk10sgp55/huCWScSpDZg+HU6fDux7GcCSABPOPvsM8vOZ2mlgwN7iaM06LG/ZmUu2LwvYexgT0lTpl50Z+FoAr+HDITcX5s0LzvtFOUsCTPj6+GNo2ZJVzTsE9G3mpfShx54tND5+OKDvY0woanXkWxJzD7I0KcD9Abwuuwxq1rQmgSDxKQkQkW7+CsSYSsnNhalT4aabUAlsLjvXM1RwUJY1CZjok+qZH2BJsGoC4uPh0kstCQgSX0vPf4nIEhG5V0Qa+CUiYypi6lQ4dQpuuingb7W2WTsO1KrHJTZU0EShvtmZHK5Zh81NWgfvTYcPh82bYdOm4L1nlPIpCVDVQcBInAl8lorIByJyhV8iM+Z8PvoIEhJg8OCAv5VKNdJTejN4+wooKrmchTGRLTV7HUtbXhDwGrezeIcKWm1AwPn8V1XVzcAjwMPAJcDfRWSDiAT+K5qJTqdOOVOL3ngjxMQE5S3npfQh4cRhWLkyKO9nTEjYt492B7OD1x/AKzkZunZ1rnMTUL72CeghIi/iTPozFLhOVS/wPH7RD/EZc66vvnL6BAShKcArPbmP82CaTVBposjXXwOwJNhJADi1AfPmwdGjwX/vKOJrTcCrwHKgp6rep6rLAVQ1B6d2wBj/++gjqF8fhg4N2lvuq9OQtYntnL4IxkSLr78mLyaWtc3aB/+9r70WCgpgxozgv3cU8TUJuAb4QFVPAohINRGJB1DVd30Nzphz5OfD5Mlw3XUQFxfUt56X0hsWLLBvJiZ6pKezskUnTlePDf57X3QRNGhgtW8B5msS8BVQq9jP8Z5txgTG3Llw8CB8//vBf+uUC51vJrNmBf29jQm648dh+XIygjU0sKTq1Z1VBRctcuf9o4SvSUBNVc31/uB5HO/jaxpTto8/dsYRDxsW9Lde3rIz1KljTQImOixaBIWFZLjRH8Crf3/IzIRjx9yLIcJV9/H44yLSx9sXQEQuBE76HpYxpSgqgkmT4JprnEQgyPJjYpnevBsXjPuEwfWvPe96BVnPDg9iZMYEwNdfg4iT/LplwADnul+61JlAyPidrzUBDwIfiki6iKQD43GWCTbG/xYscNYZd6EpwGtuSh9aHfmWtgd3uRaDMUGRng49e3KsRm33YkhNde4XL3YvhgjnU02AqmaISGegEyDABlXN90tkxpQ0aZLTGfCaa1wLYW6KM1RwyPblbGuc5FocxgRUfr7THHDnne7G0agRdOxo/QICyB9TQPUDegC9gdtF5Md+eE1jzjVzJgwaBPXquRZCdoNmbG3U0lYVNJFt5UqnY2AQZuQsV//+Tk2AqtuRRCSfagJE5F2gHbASKPRsVuAdH+My5myHD8Pq1fD4425HwryUPoxYNZ0aBafJqx7cYYrG+Evy2LKn5L0r4xP+CKTOPQV1XWwOAKdfwLvvwjffQJs27sYSgXztGNgX6KJqKZoJsPnznW8CQ4a4HQlzU/pw57LP6Lczk69TersdjjF+1y87kx0NmrG3bmO3Q3FqAsCpDbAkwO98bQ5YCzTzRyDGnNe8eRAb+12B4KLFrbqTF1PdmgRMZFKlb/Y6lro1P0BJPXpAzZrWLyBAfK0JaAKsE5ElQJ53o6pe7+PrGnO29HTo1w9q1Sp/3wA7GVeTJUndGLJ9OU+7HYwxftb24C6anDhCRssQSQJiY+HCC22EQID4WhPwOHAj8AzwfLGbMf5z4gRkZIREU4DX3LZ96LT/G5of3ed2KK4TkatEZKOIbBGRsaU8X0NExnueXywiyZ7tsSLytoisEZH1IvK7YMduztU3ex0AGa1cnCSopAEDYNkyOH3a7Ugijq9DBOeKSBugg6p+5Vk3oNy1XUXkKuBlz75vqOqzJZ6vgdO58ELgAHCbqmaJSCzwBtDHE/s7qvpnX34HEwYWLXKm6w2hJGBeSh+Y/RaDt69gQs/gz14YKkQkBngNuALIBjJEZLKqriu2213AIVVtLyIjgOeA24BbgBqq2t1TdqwTkf+palZwfwtTXGp2Jgdq1WNro+APgS2rs+I1W2L5R14e193zDz57+8EgRxXZfF1K+B5gIvBvz6aWwCflHOMtNK4GuuAMKyxZ73Sm0MBZkvg5z/YzhQZOgvBT77cKE8HmzXNm57v4YrcjOWNTkzbsrtPY+gVAKrBFVbep6mlgHHBDiX1uAN72PJ4IXCYigjOSqLaIVMdZg+Q0YKszuexMf4DzzIgZbCtadAKgd84GlyOJPL42B9wHDMRz4arqZqBpOcdYoWEqJz0devVylg8OFSLMS+nD4KyVxBQVlr9/5GoJ7Cz2c7ZnW6n7qGoBcARojHNtHwd2A98Af1PVg4EO2JQtIfcgyYd3u7doUBl2123Ct3Ua0Wv3JrdDiTi+dgzMU9XT4skYPf+cyxsuWFqhUbLL91mFhogULzRuwCk04oFfWaER4U6fhoULYfRotyM5x9y2F3Lbmhn0zNnE8qQL3A7HLaV9XSxZBpS1TyrO/CItgIZAuoh8parbzjpYZDQwGiAxMZE5c+b4GrMrcnNzQyr2Md0LztnWfvEaADqmdWZMu3Of94fEWqW/d3lyO3dk6M4NIfUZVlUonQu+JgFzReT3QC0RuQK4F/isnGMCXmhA+BQcoXQyuK20z6JeZiZ9Tp5kbaNG7C/jc6pKgeIPNdp2o2hyNcYez2Bx9w5nPefr3zSMzotsoFWxn5OAnDL2yfZ8UagPHAR+AEz1TDW+V0Tm48w9ctb1rKqvA68D9O3bV9PS0gLwawTenDlzCKXY7yil/f2xBRs4EVuD3x/rSMEaX/89lG5M9wKer8JrH6/dmbHfLiKte3doHALzF/gglM4FX//KY3Ha79cAPwWm4HTcO5+AFxoQPgVHKJ0Mbiv1s/AMC+r2859DQkKpx5VWmAVHAy5s3pHqS1bxfKezZ8vOGpnm0yuH0XmRAXQQkRRgFzAC5zotbjIwClgI3AzMUlUVkW+AoSLyHk7N3gDgpaBFbs7RL3sdK1p0oiAmMAmAL7z9Ali82NX1QyKNT30CVLVIVf+jqreo6s2ex+U1B5wpNEQkDqfQmFxiH2+hAcUKDZx2w6HiqI1TaFhPkUiWng6dO5eZALhtcetudN+zhVqnT7kdiis8bfz3A9OA9cAEVc0UkSdFxDtfyJtAYxHZAvwa58sDOB2E6+BMOpYB/J+qrg7qL2DOqJN3ggv2bmdpyxAaGljMmmbtKZRqNl+An/m6dsB2SukDoKptyzrG08bvLTRigLe8hQawVFUn4xQa73oKjYM4iQI4hcb/4RQaghUaka2w0FnT/Lbb3I6kTEuSunLvoon02r2RhW16uh2OK1R1Ck4tYPFtjxZ7fApnZE/J43JL227c0WfXemK0iCWhND9AMSfiarExoQ1dbOZAv/LH2gFeNXEu6EblHWSFhqmQNWvgyJGQmh+gpOUtL6AIIXVnZtQmASYy9MteR4FU+67aPQStbN6JLksWQVERVPPHIrjG1+aAA8Vuu1T1JWCon2Iz0W7ePOc+hJOAozXrsKFpMv2yM90OxRif9MvOJDOxHSfi3J+auywrWnRyVhTdZEMF/cXXyYL6FLv1FZGfAXX9FJuJdunpzqphrVqVv6+LliR1pU/OBqoXujNKwRhfxRQV0mPPZpa37Ox2KOd1VudA4xe+1qcUXy/gzziz+N3qa1DGoOrUBIRwLYBXRlJX4vPz6PrtVrdDMaZKOu3bQXx+HitahHYSsLVxEtSrZysK+pGvawdc6q9AjDnLpk2wd29YJAHejlT9sjNZFcLtqcaUxTsdbyj3BwBQqQapqVYT4Ee+jg749fmeV9UXfHl9E8XCoD+A1746jchq0JzU7HW8kXqT2+EYU2m9czayP74+O+snuh1K+QYMgD//2VldND7e7WjCnq/NAX2Bn+NM89sS+BnOokB1sb4Bxhfp6dC0KXToUP6+ISAjqSt9s9chWuR2KMZUWu+cjU4tQAgtGlSm/v2d4cPLon7xLr/wNQloAvRR1TGqOganT0CSqj6hqk/4Hp6JWt7+AOFQKAFLWnWh0cmjtDuQ7XYoxlRKvVO5tDuYzcrmod0UcEZ/z1Iz1i/AL3xNAlrjrOTndRpI9vE1TbTbscO5hUFTgFdGkrdfwDqXIzGmcnrlbARCvz/AGQkJ0Lat9QvwE1+TgHeBJSLyuIg8BiwG3vE9LBPV0tOd+zBKArIatmBf7QY2X4AJO71zNlKEsLp5R7dDqbgBA6wmwE98nSzoaeBO4BBwGLhTVZ/xR2AmiqWnQ/360K2b25FUnAhLkrqSutOSABNeeu3eyKYmrcmtEUad7Pr3h127INua33zlj3kX44Gjqvoyzqp/KX54TRPN5s2DQYMgJsbtSColI6krSUf30vzoPrdDMaZiVOmds5GV4dIU4DVggHNvTQI+83XGwMeAh4HfeTbFAu/5GpSJYnv3woYNYdUU4JVRbL4AY8JByqEcGpzKDZ/+AF49e0JcnCUBfuDrAkLfA3oDywFUNUdEbGigqbqvv3buwzAJWJ+QzNG4eFJ3ZpI89oty9896dngQojKmbGHXKdCrRg3o08f6BfiBr80Bp1VV8SwnLCK1fQ/JRLV586BWLecCDzNF1WJY3vICqwkwYaN3zkZy42qxpXFor89Rqv79YelSKLA1O3zha03ABBH5N9BARO4BfgL8x/ewTNSaNw8uusip6oMKfaMOJUtadeWhee/Q4ORRDteq53Y4xpxX75wNrGregaJq4dX/BnD6Bbz8srPkeO/ebkcTtnwdHfA3YCLwEdAJeFRVX/FHYCb6xOTmwsqVYdkU4JWR1AWAvtnrXY7EmPOrmX+KzvuyQn7RoDJ5Jw2yfgE+qXISICIxIvKVqs5Q1d+q6m9UdYY/gzPRpX5mprN6YBgnAaubdyQvpro1CZiQ1+3brcQWFYbPTIElJSc7U4tbvwCfVLk5QFULReSEiNRX1SP+DMpEp/qrV0Ns7HcZfhjKqx7HquYdbb4AE/J673I6Ba5sEUaTBHF2E+F/6iWT8tlMLi+l2dA63laMrx0DTwFrRORNEfm79+aPwEz0abBqFfTtG/Yrg2UkdaXbt1uodfqU26EYU6ZeuzfyTf1E9tdu6HYoVbaiRWfaH8ym3qlct0MJW74mAV8AfwTmAcuK3YypnJMnqbtxY1g3BXhlJHUltqiQXrs3uh2KMWVyVg4M0/4AHis9Ux333L3J5UjCV5WaA0RkpqpeBnRR1Yf9HJOJRosXU62gICKSgGVJF1CEkLozk4VterodjjHn2rWLFsf2h11TQEmrm3ekCKF3zkbSU8JvWHEoqGpNQHMRuQS4XkR6i0if4jd/BmiixLx5qAhcfLHbkfjsWI3arG+aEjWdA0XkKhHZKCJbRGRsKc/XEJHxnucXi0hysed6iMhCEckUkTUiUjOYsUctT4/6cK8JyK0Rz+Ymreids8HtUMJWVZOAR4GxQBLwAvB8sdvfyjvYCg1zjnnzyG3XDho0cDsSv1jSqit9cjZQvTCyJzIRkRjgNeBqoAtwu4h0KbHbXcAhVW0PvAg85zm2Os404z9T1a5AGpAfpNCj26JF5MVUZ13Ttm5H4rNlLS+gd85GRIvcDiUsVSkJUNWJqno18BdVvbTEbej5jrVCw5wjLw8WLOBIjx5uR+I3GUldic/Po+u3W90OJdBSgS2quk1VTwPjgBtK7HMD8Lbn8UTgMhERYBiwWlVXAajqAVUtDFLc0W3RItY1bcfp6rFuR+KzpUldaHAql/b7d7odSljydbKgP1XhMCs0zNmWLIGTJzkUQbN+ZSRFzWJCLYHipW+2Z1up+6hqAXAEaAx0BFREponIchF5KAjxmvx8WLo0/NYLKIP3WkuN/GstIHydNrgqSis0Sg4MP6vQEJFzCg0gARinqn8p7U1EZDQwGiAxMZE5c+b483fwm9zc3JCNLVjavP02ySLktG/PgRKfxZju4VqdXpfDE5vzg6Nrqd/9+lL3ON/fPYzOCyllm1Zwn+rAIKAfcAKYKSLLVHXmWQeHybVcnlD5m9bZvJm+J0/SYkAHV66vxFp+vq61MbnjG/Kj3EyadR+DeZ2mAAAgAElEQVR2ZnMofNZlCZVzAdxJAgJeaACo6uvA6wB9+/bVtLQ0X2IOmDlz5hCqsQXN449D797UaNbsnM/ijjBbO6C4xKZduXzdEl5YXQ2VcyvdskamlXlsGJ0X2UDx1WeSgJwy9sn2NOnVBw56ts9V1f0AIjIF6AOcdT2Hy7VcnpD5m653prR+quACstcE/1/AmO4FPO/n901J7EqvtevPet3zXV9uC5lzAd/nCUBEBonInZ7HCSKSUs4hlSk0KKvQUNUTgLfQMOHq5ElYuBAuvdTtSPwuI6krjU4epd2BbLdDCaQMoIOIpIhIHDACmFxin8nAKM/jm4FZntVHpwE9RCTec51fAqwLUtzRa9EiaNqU7PqJbkfiN0uTupB0dC/Nj+5zO5Sw41MSICKPAQ8Dv/NsisXpuHc+VmiY7yxcCKdPw9Dz9icNS0taRX5bpaeN/36ca3M9MEFVM0XkSRHxtoO8CTQWkS3Ar3FGFqGqh3BGF2UAK4Hlqhq+VT/hYvFiZwU+Ka3CNTx9t3CX/TuoLF9rAr4HXA8cB1DVHKDu+Q6wQsOcZfZsiImBwYPdjsTvdjRozr7aDegX4esIqOoUVe2oqu1U9WnPtkdVdbLn8SlVvUVV26tqqqpuK3bse6raVVW7qap1DAy0gwdh48awXp+jNBuappAbV4u+uywJqCxfG2ZOq6qKiAKISO2KHKSqU3Cq8otve7TY41PALWUc+x7l1zaYcDFrlrNeQN3z5o7hSYQlSV2jYYSACRdLljj3AwbA9JPuxuJHhdViWN6iM/2sJqDSfK0JmCAi/wYaiMg9wFfAf3wPy0SF3FynUIrA/gBeGUldSTq6jxZH97odijFOU4AI9OvndiR+tzSpC533ZlE377jboYQVX+cJ+BvOOP6PgE7Ao6r6ij8CM1Fg/nwoKIjI/gBeGZ5+AZHeJGDCxKJF0K1bRNa8ZSR1oRpKn102hXBl+Nox8FfAelX9rar+RlVn+CkuEw1mz4bYWBg40O1IAmZ9QjJH4+IjunOgCROqTk1AhPUH8FrZvBMFUs06B1aSr80B9YBpIpIuIveJSOSMOTGBN2uWUyDFx7sdScAUVYthecsL6LfTCibjss2b4dAhpz9ABDoZV5O1zdpZH5xK8rU54AnPHP73AS2AuSLylV8iM5HtyBFYtiyi+wN4LW7djY4HviEh95DboZhotmiRcx+hSQDA0pZd6LV7E7GFtqRMRfk8WZDHXmAPcABo6qfXNJEsPR2KiiK6P4DXvGRnTYSBO1a6HImJaosWOX0BOof38sHnk5HUlZoFp+m2J+IX7vIbX/sE/FxE5uBM89kEuEdVI2cpOBM4s2dDjRoR/a3Ea11iWw7Uqsfg7cvdDsVEs8WLITXVmZcjQi1LugCIioW7/MbXeQLaAA+qqn3FMZUzaxZcfDHUrOl2JAGnUo2vk3szJGuF0zkrgmZqM2HixAlYtQrGjnU7koDaX7sh2xq2sPkCKqFKNQEiUs/z8C/ANyLSqPjNf+GZiHTwoFMgRUF/AK/0lN4kHD9M531ZbodiotGyZVBYGLEjA4rLSOrKhbvWO82NplxVbQ74wHO/DFjquV9W7GdjyjZ3rvONOJqSgOReAAzevsLlSExUWrzYuY+CJGBpUhcanTzqTI9sylWlJEBVr/Xcp6hqW8+999bWvyGaiDN7tjMsMDXV7UiC5tu6TdjUuDWDsywJMC5YtAjatoWmkd9v27uYEF9/7W4gYcLXjoEzK7LNmLPMmgWDBkFcnNuRBFV6Sm/671xLjfw8t0Mx0WbRoqioBQDIatiCffENLAmooCp1DBSRmkA80EREGgLenk71cOYLMKZ0e/dCZib88IduRxJ06cm9uWvpp6RmZ5Ke0sftcEwESx773eKqSYf38PWuXTy6rx7vjI2CRVdFWJrUhastCaiQqtYE/BSn/b8zZ/cH+BR4zT+hmYg0Z45zH0X9AbwWt+pGXkx1BmXZYBoTPEM8TVDzPf1SosHSpC6wbRvk5LgdSsirap+Al1U1BfhNiT4BPVX1VT/HaCLJ7NnOhCUXXuh2JEF3Mq4mS5O6MMTmCzBBNGj7CnbVTWBroyS3QwmaM/0C5s93N5Aw4Ou0wa+ISDcRuVVEfuy9+Ss4E4FmzYIhQ6C6r1NUhKf05D5csC+LhNyDbodiokC1okIG7ljF18m9omp+inVN2zqdj61JoFy+dgx8DHjFc7sUZ96A6/0Ql4lEOTmwaVNUNgV4pac4Uwhbk4AJhh57tlA/7/iZ8y5aFMRUd2YjtSSgXL5+HbsZ6AmsUNU7PasIvuF7WCYizZ7t3F966Vkdl7zGdC/gjgjvuLSuaQr74+vbUEETFIO3L6cIYX6bnm6HEnyDBsFTT8GxY04TpCmVrwsInVTVIqDAM4vgXsDmCTClmz0bGjaEnlFYIHk4Uwj3cpIAm9HMBNigrJWsbdaOQ/H13Q4l+AYNcq4x7+qJplS+JgFLRaQB8B+c0QHLgSU+R2Ui06xZcMklEb2ASUWkJ/ch4fhhWLPG7VBMBKudd4I+ORuc/gDRaMAAqFbNmgTK4WvHwHtV9bCq/gu4Ahilqnf6JzQTUXbsgO3bo7o/gJd3CmGmT3c3EBPRBuxcQ2xRIenJUTonRd260KuXs2y5KVNVFxDqU/IGNAKqex6Xd/xVIrJRRLaIyDnLWolIDREZ73l+sYgkl3i+tYjkishvqhK/cUGx/gDRbm/dxmxo0iZikgC7nkPT4O0rOBFbg2UtL3A7FPcMGuQ0B+Tnux1JyKpqx8Dnz/OcAkPLelJEYnAmFLoCyAYyRGSyqhZf+/Eu4JCqtheREcBzwG3Fnn8R+LKKsRs3zJ4NCQnQtavbkYSE9JTedE7/0lniNT7e7XCqzK7n0DU4ayWLW3XjdPVYt0Nxz6BB8Pe/w4oVUbVWSWVUdbKgS89zKzMB8EgFtqjqNlU9DYwDbiixzw3A257HE4HLRJxBriJyI7ANyKxK7MYFqk5/gLQ0p43OkJ7cG/LyIqGq0q7nENTi6F7aHczm6+ToGhp4joEDnXvrF1Amn4YIikg88GugtaqOFpEOQCdV/fw8h7UEdhb7ORsoubLFmX1UtUBEjgCNReQk8DDOtw6rOgwXW7dCdrY1BRSzpFVXqFHDaRK48kq3w/FFwK9nERkNjAZITExkjnfq6TCTm5sbtNj/lL8MgLbDejAmqSAo71lRibWc4cCB5v2s+7doQe4nn5DZJ3T6RgTzXCiPr/ME/B/OqICLPT9nAx8C50sCSpu2Siu4zxPAi6qaK+XMfhUuBUconQyB0vzzz+kELKldmxOe37W0QiBYhUNoqM7Bbt2ImzSJpdddd86zYXReBPx6VtXXgdcB+vbtq2lpaVWL1GVz5swhWLF//vNH2VOnEX84mAKHQmumwDHdC3h+TeBnDM0ameY8uPxyan35JWmXXBIysyYG81woj69/iXaqepuI3A6gqielvP/OTqLQqtjPSUDJVR68+2SLSHWgPnAQ5xvGzSLyF6ABUCQip0pbryBcCo5QOhkC5vXXoXlzUn/0ozMXYWmTAgWrcAgVv7j1Vnj4YdI6doQWZy++GUbnRVCuZ1MJhYUMzFrJzPapIfNPz1WDBsE778DmzdCxo9vRhBxfG2hPi0gtPJm/iLQDylssPQPoICIpIhIHjAAml9hnMjDK8/hmYJY6BqtqsqomAy8Bz1iBEeK8/QEuvdQKpJKGDXPuZ8xwNw7f2PUcalasoOGpY98NRY12gwY599YvoFS+JgGPAVOBViLyPjATeOh8B6hqAXA/MA1YD0xQ1UwReVJEvOsOvInTZrgFp8/BOcOOTJhYsQK+/RYuu8ztSEJPjx7QtGlYDxW06zkEeZLKaFo6+Lw6d4bGjS0JKEOV61491f4bgJuAATjtfg+o6v7yjlXVKcCUEtseLfb4FHBLOa/xeOWjNkE3fryzYuCNN7odSeipVg2uuMJJAoqKwnbkhF3PIWb6dGeNitoN3Y4kNIg4tQGWBJSqykmAqqqIfKKqFwKRveqLqRpVGDfOqfZu1MjtaELTsGHw/vuwahX0jvLhXMZ3x4/D/PnM631uZ9NoU3yRsruOJfDHzZsZ9LM3yW7Q7Mz2rGeHuxFaSPH1q8ciEennl0hM5Fm0CL75Bm67rfx9o9UVVzj3YdwkYELI3LmQn2/zA5QwrcMAAK7ZON/lSEKPr0nApcBCEdkqIqtFZI2IrPZHYCYCjB/vjIW/oeTcMeaM5s2he3dLAox/zJgBNWuSkdTF7UhCSnaDZqxu1t6SgFL4mgRcDbTDmSb4OuBaz72JdoWFMGECXH011I/CZUwrY9gwp73y+HG3IzHhbvp0GDyYvNgabkcScqZ0GkSv3ZtoeWSv26GEFF9XEdxR2s1fwZkw9vXXsHs3jBjhdiShb9gwOH0a5s1zOxITznbtgnXrvht6as4ypZMzhfBVVhtwlvDsjmxC37hxzsI4117rdiShb/Bgp9kkvOcLMG7znj/efibmLN80bM7axHYM32ijBIqzJMD4X0EBTJwI110HtWu7HU3oq1ULhgyxfgHGNzNmQGKi08fElGpKp4H0ydlI86P73A4lZFgSYPxv1izYv99GBVTGsGGQmQk7d5a/rzElFRU5ScDll4ftfBPB4G0SuHrjApcjCR12thj/Gz8e6tZ1OgWairneM7ne+++7G4cJT6tXw7591h+gHFmNWrKuaQrXWJPAGZYEGP86fRo+/tiZIbBmTbejCR8dOzp9A954w5lkyZjK8DYlXX65u3GEgSmdBtJ313oSj5U7uW1UsCTA+Nf06XD4sI0KqIp77oGtW50JX4ypjBkzoFu3c1ajNOea0slZUMiaBByWBBj/GjcOGja0byRV8f3vO3MqvPGG25GYcHLyJKSn26iACtrWOIkNTdpYk4CHJQHGf06ehE8/df6ZxcW5HU34iY+HkSNh4kSqHzvmdjQmXKSnQ16e9QeohCmdB9E3ez3k5LgdiussCTD+M2UK5ObaqABf3H035OWR+NVXbkdiwsX06U7SPWSI25GEjSmdBlINdfovRbkqryJozDnGj4emTWk3/QSFX9nCkuUpvspZcZ8ltqP5J5NJrnkFWc/ZZEumHDNmOEvlxse7HUnY2NKkNZsat6bjxIlw//1uh+Mqqwkw/pGbC59/DjffTGG1GLejCWvje15Jk2+y6L5ni9uhmFC3Z48zPND6A1Tal50GOlN179njdiiusiTA+Mdnnzl9AmxUgM8+7XIJ+XFx3L5qmtuhmFA3bpxzP3y4u3GEoS86D3SG406a5HYorrIkwPjHuHHQsiUMHOh2JGHvWI3abOk/kOvXz3VqWIwpTWEhvPKKc83ZVMGVtqlJG+jcGT780O1QXGV9AozvDh+GqVPhvvtsylI/yUy7ggvSZzsF1J13uh2OCSHeviRDtyzhrW3buK/bLXxRRv8Scx4icPPN8MwzsHcvNG3qdkSusBLb+O6TT5yZAm1UgN/s7ngBWxsl2ZwBpkx3LPuM3XUaM63jRW6HEr5uucVZdyGKmwQsCTC+Gz8ekpMhNdXtSCKHCON6DIMFC5w14o0ppv3+bxiStYJ3+wynIMYqdKuse3fo0CGqmwRcSQJE5CoR2SgiW0RkbCnP1xCR8Z7nF4tIsmf7FSKyTETWeO6HBjt2U8L+/c4Qpdtuc6rXjN983G0oVK8Ob77pdijnZddz8I1a/jl5MbGM63ml26GENxGnNmDOHGcBpigU9CRARGKA14CrgS7A7SLSpcRudwGHVLU98CLwnGf7fuA6Ve0OjALeDU7Upkwff+x0ULJRAX53oHYDuOEGeOcdZ0a4EGTXc/DVO5XL99fO5NMul3Awvr7b4YS/m292yrBPPnE7Ele4UROQCmxR1W2qehoYB9xQYp8bgLc9jycCl4mIqOoKVfXO85gJ1BSRGkGJ2pRu3DhnBbyePd2OJDLdc49T2zJ5stuRlMWu5yC7ZfUM4vPzePvC69wOJTL06gXt2sHEiW5H4go3koCWwM5iP2d7tpW6j6oWAEeAxiX2+T6wQlVD8ytSNNizx1nxbsQIawoIlMsvh9at4T//cTuSstj1HEyFhYxa/jmLk7qSmdjO7Wgig3eUwMyZcOCA29EEnRs9Skr7b1FyAfXz7iMiXXGqFMtcMUNERgOjARITE5kzZ06lAw2G3NzckI2tPG3/9S9aqbKkXTtOFvsdxnQvqNLrJdaq+rGRxvtZzElPp83QoaT8978s+t//ONW8uduhlRTw6zlcruXy+ONabzx/Pt2PfMvaUaPC9loJpev8lfc/BSChTgIjCguZOfYx1qWdPfti95b+b3IJpXLfjSQgG2hV7OckoORSTt59skWkOlAfOAggIknAJODHqrq1rDdR1deB1wH69u2raWlp/orfr+bMmUOoxnZe27c7w2pGjaL/j3981lN3VHHM8pjuBTy/xno6w3efRdbINGjbFt5+mwHr18Ptt7sdWkkBv57D5Vouj1+u9aeeIqduE35RYyCFa8Jzeu6QvM61EwPrJxIzcxHPN776rKeyRqb5/e1Cqdx3ozkgA+ggIikiEgeMAEo2eE7G6SgEcDMwS1VVRBoAXwC/U9X5QYvYnGvsWKfn+lNPuR1J5GvdGq68Et56y+nAFFrseg6WzEyYOZP3el9j63P4mwhfdBrIwB0rqX8yupbxDnoS4GkTvB+YBqwHJqhqpog8KSLXe3Z7E2gsIluAXwPeYUf3A+2BP4rISs8tOqd5ctPChTBhAvz2t85UwSbw7r4bdu2CaaG1noBdz0H0yitQowb/s2GBATG5SxqxRYXckxFdEwe5UiejqlOAKSW2PVrs8SngllKOewqwr55uUoVf/xqaN3eSABMc110HCQnODILXXON2NGex6zkIDh2Cd9+FkSM5ZMMCA2JdYlsmdUnjniWTGN9jGDsbNHM7pKCwGQNN5UyYAIsWOc0AtWu7HU30iIuDO+5wVmuM8qVPo9Jbb8GJE/CLX7gdSUR77pI7KKxWjd/PfsvtUIImxHpnmJB26hQ7f/oAuQnJDN/QhCJbtCS47roL/vpXeO45ePFFt6MxwVJYCK++CkOGOGPax+1yO6KItadeE14bcCu/TX+Xi3asYmGbyJ//xGoCTMW98gqtjnzLU0Pvpsg6JgVfp07w85/Dyy/DvHluR2OC5fPPISsLfvlLtyOJCm/0u5Gd9RN57KvXiSkKuY64fmdJgKmY/fvh6aeZ2a4f85N7uR1N9PrLXyAlxWkayM11OxoTDH//O7Rq5UwhbQIuL7YGT1/6Ezrv38HtK6e6HU7AWRJgKuaJJyA3l2fSfuJ2JNGtTh3473+db4bWMTPyrV0Ls2bBffc5Q3JNUEzteDELW3dnTPp7cPCg2+EElCUBpnwbNsA//wmjR7O1Savy9zeBNXiwM0LjX/+C6dPdjsYE0iuvQM2azhBREzwiPHHZaOrlHYfHHnM7moCyJMCU76GHID4eHn/c7UiM15/+BJ07O50FDx92OxoTCAcPnhkWSOOSSy2YQNvQNIUPel3lfAFau9btcALGkgBzfrNnO8PS/vAHaGrzuISMWrXg7bdh92548EG3ozGB8PTTcPKkDQt00QuDRkLdus41piWXxIgMlgSYshUWOtXObdrAAw+4HY0pKTUVfvc7JxkI3aWGTVW88w688AL87Ge2TLeLDsXXhyefdFYY/PRTt8MJCEsCTNnefRdWroQ//9lplzSh549/dP5JjB7tjOAw4W/+fLjnHhg61BkZYNz1s59Bly4wZgycOuV2NH5nSYAp3fHjThNAaiqMGOF2NKYscXHOt8aDB50e5Ca8ZWXB977n1L59+CHExrodkYmNhZdegm3bInKSLksCzLkKC52JSXJynCpJKW05eBMyevRwOm1OmADjx7sdjamqY8ecNSLy851+OI0auR2R8briCmeehqefdsrFCGIDT83Z8vKc3sgffeTUBAwc6HZEpiIeeshps7z3XrjkEmgWHYufRIzCQrj9dgoy1zHq1ieZ/39bgC1uR2WK+9vfoGvX7/rhRAhLAswZ3X71Ia9//BQXf7OaPw29mzcLLgJbH8BVyRX4/LOeHe5MJPP229C7t9OePHmy1eCEk4cfhi++4PFh99qMnKGqfXv41a+ctTtGjXL6bEQASwKiQEX+kTQ5fohxHz5Op31ZPHjtGD7pemkQIjN+1bkzPPOMM6Lj7bedqYVN6HvzTXj+ebj/ft6rfZXb0ZgSipeftfMu5POGzWkx7EoeGXYfH/a4AvAk4mHK+gQYWh3ew8T3HqLdgWzuuemPlgCEswcecFabW7rU7UhMRcyd6ywKNWxYRHY6izTHa8TzvR89T0ZSV/765cs8PuNfVC8scDssn1hNQJS7YO823p7wGHGF+Ywc8RTLW17gdkjGF9WqwdSpzmRCJrRt3Qrf/z60bet06LS1AcLC4Vr1GHXrkzw857+MzpjEBfuy4LcXhe1kanbWhbmKVPWXJXXnWt6Y+CS5NeL5wYin2dKktR8jM66xBCD0HTnijAQoKnKWCm7QwO2ITCUUVovhmaF3sS6xLc9OfQX69oVJk+DCC90OrdKsOSBKXbF5Ee+O/yN76zTi5h/+xRIAY4JBFT75BAYNgs2bnVE47du7HZWpok+6XsrNI//i/DBoELz/vrsBVYElAVGm/slj3LtwAv+a9Azrm6Zwy8jnyKkXntVYxoQNVRovWOB8Y/ze95yZ5z7+GC61/jfhbm2z9k4fnNRU+OEP4Te/gYLw6SdgzQHRQJWLvlnDiFXTuGrTAmoU5jOjfSoPXPdbTsRZ1bExAaMKX37JqrsepOeezexo0Iy/X/MrPumaRuH8ajDfhuBGhKZN4auvnCGEzz8Pq1bBBx9AQoLbkZXLkoAIlpB7iJvXfsWtq6eTcmg3R2rU5n89r2R8z2Gsb9rW7fCMiVyqMH26sxb94sU0qp/IV/fcz8/qX05BjBW7ESk2Fl591Zmr4957oWVLp6bnhhvg+ushKcntCEvlytkoIlcBLwMxwBuq+myJ52sA7wAXAgeA21Q1y/Pc74C7gELgl6o6LYihB1VVOv1VKypkyPbl3L5qGpdtWUJ1LWJxq268PPAHfNnxYvJiawQgUhPN7Hr2UIXsbMjIcL4NLlgArVvD668zdFMCv+wlFKyxBCASnV1WN6PjD1/gprWzGJaxiLbTp8N997GqWQd63j/KSQpCaFnioJ+RIhIDvAZcAWQDGSIyWVXXFdvtLuCQqrYXkRHAc8BtItIFGAF0BVoAX4lIR1UtDO5v4TvvSTOmewF3VOGffZ28EyQfyiHl4C7aHN5NyqEckg/m0O5gNg1O5bI/vj5v9LuRCT2Gsa1xaGagxj8qmiwGYkKTqL2ejx2DtWthzRpYvdq5rVkDhw87zyclwT//CT/5CcTFkT/2CyB82omNbzYlJPPspT/h2bQ7aXcgmys3L+SKzYvhkUfgkUfo36KFszDbZZdBq1bQvDk0buzKLJ9upKWpwBZV3QYgIuOAG4DihcYNwOOexxOBV0VEPNvHqWoesF1Etnheb6FfIywqcobwFFda5ubdpnrm8YVPTgdAUMTztGgRMZ5btSLnvl1RETFaSJMd+XTbI1QvKqRW/ilqnz5FfP5J4k+fonb+qTPbauWfom7ecVof/pbkQzkknDh8Vii76zQmq1ELpnQaSHpyb77q0J/8GFuBzARc6F/PR49Cbq4zP7/3VlBw7uP8fOef+7FjzjHHjvHXiUupe/oEtU+fpPbpkzQ4eYz2B3bS+si3Z17+WFwtNiYksyH5IjYkJLMxoQ2rmnfidFYsPDrDr7+KCTMibG3Sin80acU/LrqVrF/2hs8+48Rbb1HrlVecBdq84uKcZKB5c2jRwrl5H994Y8CGkbqRBLQEdhb7ORvoX9Y+qlogIkeAxp7ti0oc29IfQRX/NtXk+CGWvvqjKr3Osiocc3s5z+fFVOdEbC2Ox9VkV/1EZrZPJathC7Y3bEFWoxbsaNCMU7E1qxKuMb4K+ev5mamv8INVVWtl+C2QFxNLblwtcmvEkxsXz+rmHRnfYxgbmiazMSGZ7HpNbZ0GUyHJf18BJDHml7/n30tP03nfdprmHqJp7kEScw/SNPcA368TAxs2wKxZ39UsDRkSUUlAaVdLya/ZZe1TkWOdFxAZDYz2/JgrIhsrGuCOMt4oQJoA+8+7R2EBFB6DU8fg6D7YuTY4kQXZLyvyWUSJQHwW8lyFdmtT2ZctZZtfr2dfrmWAkZ5blRXmw8l8OHnU+XnvNtiQXumXsfPbYZ/Dd59BZmlPlraxXbuqvE2FrmU3koBsoFWxn5OAkgs0e/fJFpHqQH3gYAWPBUBVXwde91PMASMiS1W1r9txhAL7LL4TRp9FwK/ncLmWyxNGf9OAss8htD4DNyYLygA6iEiKiMThdAyaXGKfycAoz+ObgVmqqp7tI0SkhoikAB2AJUGK2xhzLruejQljQa8J8LQJ3g9MwxlS9JaqZorIk8BSVZ0MvAm86+kodBCnYMGz3wScTkcFwH1h0ZPYmAhl17Mx4U00hMYrRiMRGe2p7ox69ll8xz6LyGN/U4d9DqH1GVgSYIwxxkQpW0DIGGOMiVKWBLhERG4RkUwRKRKRviWe+52IbBGRjSJypVsxBpOIXOX5fbeIyFi34wkmEXlLRPaKyNpi2xqJyAwR2ey5b+hmjMY30Xp+27kNItJKRGaLyHpPmf+AZ3tIfA6WBLhnLXATMK/4xhJTqV4F/MMzNWvEKjb17NVAF+B2z+cQLf6L87cubiwwU1U7ADM9P5swFOXn93+xc7sAGKOqFwADgPs8f/+Q+BwsCXCJqq5X1dImPTkzlaqqbge8U6lGsjNTz6rqacA79WxUUNV5OL3mi7sBeNvz+G3gxqAGZfwpas9vO7dBVXer6nLP42PAepyZMUPic7AkIPSUNg2rX6ZSDWHR+DuXJ1FVd4NTiABNXYvz+yQAACAASURBVI7HVJ2d32eL2nNbRJKB3sBiQuRzsHUtA0hEvgKalfLUH1T107IOK2VbpA/hiMbf2UQPO78NIlIH+Ah4UFWPSoisN2FJQACp6uVVOKzCUyNHkGj8ncvzrYg0V9XdItIc2Ot2QKbK7Pw+W9Sd2yISi5MAvK+qH3s2h8TnYM0BoScap1KtyNSz0ab4VLujgLJqjkzos/P7bFF1bnuWzX4TWK+qxdYODo3PwSYLcomIfA94BUgADgMrVfVKz3N/AH6C06v0QVX90rVAg0RErgFe4rupZ592OaSgEZH/AWk4K4t9CzwGfAJMAFoD3wC3qGrJDlYmTETr+W3nNojIICAdWAMUeTb/HqdfgOufgyUBxhhjTJSy5gBjjDEmSlkSYIwxxkQpSwKMMcaYKGVJgDHGGBOlLAkwxhhjopQlAcYYY0yUsiTAGGOMiVKWBJiA86yhnSYivxeRN9yOxxhjjMMmCzLGGGOilNUEGGOMMVHKkgATcCKSJSKXi8jjIvJese0fisgeETkiIvNEpKubcRpjziUiD4vILhE5JiIbReQyzwJnL4lIjuf2kojU8OyfJiLZIjJGRPaKyG4RudPt38OUzpIA46YvcVZJbAosB953NxxjTHEi0gm4H+inqnWBK4Es4A/AAKAX0BNIBR4pdmgzoD7QErgLeE1EGgYvclNRlgQY16jqW6p6TFXzgMeBniJS3+WwjDHfKQRqAF1EJFZVs1R1KzASeFJV96rqPuAJ4EfFjsv3PJ+vqlOAXKBTsIM35bMkwLhCRGJE5FkR2SoiR3G+XYCz5KgxJgSo6hbgQZwkfa+IjBORFkALYEexXXd4tnkdUNWCYj+fAOoEOFxTBZYEGLf8ALgBuByn2jDZs13cCsgYcy5V/UBVBwFtAAWeA3I8P3u19mwzYcaSAOOWukAecACIB55xNxxjTEki0klEhno6/Z0CTuI0EfwPeEREEkSkCfAo8N55XsqEKEsCjFvewalC3AWsAxa5G44xphQ1gGeB/cAenE68vweeApYCq4E1OB17n3IpRuMDmyzIBJyIfAP8EKfqP0lVf+JySMYYY7CaABNgIpIAJOB86+8CbHc3ImOMMV6WBJiAEZF+wGbgFeATIAn4j6tBGWOMOcOaA4wxxpgoZTUBxhhjTJQKaBIgIld55preIiJjS3n+1yKyTkRWi8hMEWlT7LlRIrLZcxtVbPuFIrLG85p/FxEbV26MMcZUQcCaA0QkBvh/9u47Pur6fuD4652EhL0hkIAQGWHLng4cUFyIq6C2irX1p1Xb/uTXSn/1Z6tddlhrW1vrRuvCjYIioEFlrwCywzSDITuMhJD374/v9+AIGZfk7r433s/H4x539x2fe98ld/e5z3h/NgKjgFxgCXCTqq71O+ZiYJGqHhWRu4GRqjpeRJrjTD8ZiJOcYhkwQFX3i8hi4Mc4U8pmAH9T1Y8qi6Vly5basWPHWj2fI0eO0KBBg1qV4SWL31vREv+yZcu+UdVWXsdRkWC8lyF6/h7VFavPC+y5VVeg7+WkoD7qmQYDOaq6BUBEXsfJEHeqEqCqn/kdvxBnGhk4i1TMUtV97rmzgDEikgU0VtUF7vaXgHE4C9FUqGPHjixdurRWTyYrK4uRI0fWqgwvWfzeipb4RWR71Ud5JxjvZYiev0d1xerzAntu1RXoezmU3QHpwNd+93PdbRW5g9Nf5hWdm+7eDrRMY4wxxlQglC0B5fXVl9v3ICLfwWn6v6iKc6tT5p3AnQCpqalkZWVVEW7lCgsLa12Glyx+b0V7/MaY2BTKSkAu0N7vfjvKWWBCRC7DWZv6IndJWd+5I8ucm+Vub1dVmQCq+jTwNMDAgQO1tk0t0d4UZfF7K9rjN8bEplB2BywBuohIhogkAxOAaf4HiEg/4N/AWFXd7bdrJjBaRJqJSDNgNDBTVQuAwyIy1J0VcCvwfgifgzHGGBOzQtYSoKolInIvzhd6IvC8qq4RkUeApao6DfgTzhrTb7oz/Xao6lhV3Sciv8apSAA84hskCNwNvAjUwxlDUOmgQGOMMcaUL5TdAajqDJxpfP7bHvK7fVkl5z4PPF/O9qVAryCGaYwxxsQlyxhojDEmMh09CsePex1FTAtpS4AxxhhTIydOQPfukJsLGRn0at0azj/f2ea7NGnidZRRzyoBplo6Tp5e5THbHr0yDJEYY2LazJmwYwfcdhscPUrdpUvhiSeguPj0MW3bOpWB73wHbr/du1ijmFUCjDHGRJ4XX4RWreCZZ6BOHZZmZTHy/PNh61ZYvx7WrXMuc+fCj38MN90Edet6HXXUsUqAMcaYyLJvH3zwAdx9N9Spc3p7UhJ06eJcrr7a2TZrFoweDR99BNde6028UcwGBhpjjIksr7/uNPtPnFj1sRdf7LQYvPFGyMOKRVYJMMYYE1lefBH69IG+fas+NikJrr/eaTk4ciTkocUaqwQYY4yJHOvWwZIlzoDAQI0f70wnnF71wGVzJqsEGGOMiRxTpkBiItxyS+DnXHABtGljXQI1YJUAY4wxkeHkSXj5Zbj8ckhNDfy8xES48UaYMQMOHw5dfDHIKgHGGGMiw+zZkJ9fva4An/HjneyC06ZVfaw5xSoBxhhjIsOUKdCs2enpf9UxbBi0a2ddAtVklQBjjDHeO3gQ3n3XSfqTklL98xMS4Nvfho8/hv37gx9fjLJKgDHGGO9Nneo059ekK8Bn/HhnzYH33gteXDHOKgHGGGO8N2WKsw7AoEE1L2PQIMjIsC6BaghpJUBExojIBhHJEZHJ5ey/UESWi0iJiNzgt/1iEcn2uxwXkXHuvhdFZKvfvgCySRhjjIlYOTkwb57TCiBS83JEnC6B2bPhm2+CF18MC1klQEQSgSeBy4EewE0i0qPMYTuAicCr/htV9TNV7auqfYFLgKPAJ36H/NS3X1WzQ/UcjDHGhMGUKU6f/ne+U/uyxo93phq+807ty4oDoWwJGAzkqOoWVS0GXgeu8T9AVbep6iqgtJJybgA+UtWjoQvVGGOMJ0pL4aWXYNQoSE+vfXl9+0LXrtYlEKBQVgLSga/97ue626prAvBamW2/FZFVIvK4iNRgGKkxxpiIMHcu7NhRuwGB/kSc1oCsLNi1KzhlxrBQLiVcXseOVqsAkbZAb2Cm3+afAzuBZOBp4AHgkXLOvRO4EyA1NZWsrKzqPPRZCgsLa12Gl4IV/6TeJVUeE4rXyV7/yCUiY4AngETgWVV9tMz+FOAlYACwFxivqttEpA7wLNAf57PoJVX9fViDN9578UVo3BjGjaPj5Ipz/0/qXcLEydPZ9uiVVZc5fjz8+tfw1ltwzz3BizUGhbISkAu097vfDsivZhnfBt5V1RO+Dapa4N4sEpEXgP8p70RVfRqnksDAgQN15MiR1XzoM2VlZVHbMrwUrPgnVvIm9dl2S+0fpyx7/SOT39ifUTjv+SUiMk1V1/oddgewX1U7i8gE4A/AeOBGIEVVe4tIfWCtiLymqtvC+yyMZwoL4e234eaboV694JXbs6dzeeMNqwRUIZTdAUuALiKSISLJOM361c3neBNlugLc1gFERIBxwFdBiNUYUzNVjv1x709xb78FXOq+fxVoICJJQD2gGDgUnrBNRHj7bWf532B1BfgbPx6+/BLy8oJfdgwJWUuAqpaIyL04TfmJwPOqukZEHgGWquo0ERkEvAs0A64WkYdVtSeAiHTEaUmYW6boV0SkFU53QzZwV6ieQ1wrKoKlS53MWwcPOpcDB5ictZxGRUdoVHSURkVHOCkJ3H/VJA7Vbeh1xMYb5Y39GVLRMe7nwkGgBU6F4BqgAKgP/Leq7gt5xCZyvPgidO4Mw4cHv+zx4+Ghh+DNN+EnPwl++TEilN0BqOoMYEaZbQ/53V6C001Q3rnbKGcgoapeEtwoTbkefBD+/OezNt+emMThlAYcTqnPkeT69Nq1mbFr5/Kf/gH005lYFMjYn4qOGQycBNJwfgh8ISKzVXXLGScHeXwPxO4YjWh6XnV37mRoVhZbv/c9ts91futVNuYotZ6zP9DntzrvIBM6ZFDy5FO81Sqj3GN6pzepdtyh4OXfLaSVABOlSkvhtdfg0kvhd7+DJk2gaVNo0oTMX80549CZz/2QseusEhDHAhn74zsm1236bwLsA24GPnbH/OwWkXnAQOCMSkCwx/dA7I7RiKrn9etfA5Dx0ENkdOgAVD7maFLvEh5bnRTwmKOJk6dztOOFPDB3Cq9/uY+8Jq3POiYU45dqwsu/m6UNNmebP9/pR7vjDhg8GDIznbW969Y969Bp3S9icO5a0g7t9iBQEwECGfszDfB1+t4AfKqqipMs7BJxNACGAuvDFLfxkqqTIOjii8GtAITCh90uAODK9V+E7DGinVUCzNmmTnW+8K+6qspDP+h+IQBXrbM3WTxS1RLAN/ZnHTDVN/ZHRMa6hz0HtBCRHOB+wJdC/EmgIc7g3iXAC27yMBPrVq2CzZuDkyGwEl83bUN22y5cve7zkD5ONLPuAHOmkyedubVXXAGNGlV5+I5mbclu25Wx6z7n6SHXhyFAE2kCGPtzHGc6YNnzCsvbbuLA/PnO9SWhH+L1YbcLePCz5+mwP5/tzdJC/njRxloCzJnmzYOCAmcRjgBN634RvXZt5ty9uSEMzBgTMxYsgDZtQtoV4DPd7RKw1sryWUuAOdPUqU7SjgC6Anw+7HY+D376LGPXzeWv598SwuCMMTFh/nwYNqxGKwZWllWwPAWNW7E0vTtXrf+CJ4ePr/bjxTprCTCn+boCrroKGjQI+LTdjVqw6JxeXL3uC2fAjzHGVGT3bmc8QChyA1Tgo67D6b5nG6mHbXnhsqwlwJz2+efOghvV6Arwmdb9In4/8x/03L2l6oONMXHH9wt+1KaFPANcn60sq+av+ppa1s5Zxb5v/kZmZrYMy2NGC2sJMKdNnQr16zuDAqvpo8zhnEhI5Oq1ZRM8GmPMaf3z1lOckMRXbTqH7THXtj6XosQk+uXbDNSyrBJgHCUlTh7vq692KgLVdKBeYz7P6O90CZSWhiBAY0ws6J+/njWpnShKSg7bYxYn1WFNaif65W8I22NGC6sEGMfcubBnT426Anymdb+Q9MN7Tk//McYYP0knS+hTsInl6d3C/tgr0rrRZ2cOSSerXg49nlglwDimToWGDeHyy2tcxKwuQzmWlOKkHDbGmDK6795KvZIilqV3D/tjr0jLpF5JEd32bAv7Y0cyqwQYOHHC6QoYO7ZWa3ofTa7HnM6DnVW7Sqy2bYw5U3+3T355mjctAYCNCyjDKgEGPvsM9u6tVVeAz7TuFzrdCp9+GoTAjDGxZEDeOvIbtWRn4/CP0M9r3IpdDZvbuIAybIqgcboCGjUi84tSihbUbsrO3HMHQOPGTpfA6NFBCtAYEwv6561nuQddAQCIsCIt01oCyrCWgHh34gS88w5cc01QRusWJSXDddc5ZR4/HoQAjTGxoPXhvbQ7tNuTrgCfFWmZZOwvoNnRg57FEGlCWgkQkTEiskFEckRkcjn7LxSR5SJSIiI3lNl3UkSy3cs0v+0ZIrJIRDaJyBvu8qWmpubMgf37g9IVcMpNN8GhQ/DRR8Er0xgT1U6NB/BgZoCPb1xA34KNnsUQaULWHSAiiThLhY4CcoElIjJNVdf6HbYDmAj8TzlFHFPVvuVs/wPwuKq+LiJPAXcA/wpq8PHkjTegSROn6X7e7OCUeckl0Lq10yVw7bXBKdMYE9X6562nKLEOa1LP9SyG1amdKZEE+uWt57NOgwJah2Dbo1eGITLvhLIlYDCQo6pbVLUYeB24xv8AVd3mrh8eUHYZERHgEuAtd9MUYFzwQo4zxcXw7rswbhykpASv3KQkuPFG+OADOHw4eOUaY6LWgLx1rGrThROJdTyL4VhyXda3zjjVKmFCWwlIB772u5/rbgtUXRFZKiILRcT3Rd8COKCqvvln1S3T+Js1Cw4eDG5XgM9NNzljAt5/P/hlG2OiS1ERvXbleNoV4LMiLZPzCjaSUHrS61AiQihnB5S3RmR1lpg7R1XzReRc4FMRWQ0cCrRMEbkTuBMgNTWVrKysajz02QoLC2tdhldW5x0ktR78/ZUzv5Ave+oJMuo34Lk9Ryl95X0m9Q7O42VlZUFpKUNTUznyj3+wul27WpcZza8/RH/8xtTKihWknCyJmErAd1fMoPPer9nYqqPX4XgulJWAXKC93/12QH6gJ6tqvnu9RUSygH7A20BTEUlyWwMqLFNVnwaeBhg4cKCOHDmyBk/htKysLGpbhlcmTp7OpN4lPLb69J87ueQEty1exLTM4fxpXc0TBJVn2y0jnRu33Ubdv/yFkb17Q4sWtSozml9/iP74jakVN5X48jSPpgf6WX4qadAGqwQQQHeAiPSqYdlLgC7uaP5kYAIwrYpzfI/ZTERS3NstgRHAWlVV4DPAN5PgNsDam2vgwq3LaVx8lA+7XRC6B5kw4fTCRMaY+LVgATuapLKnYTOvI2FbszT2121kSYNcgYwJeEpEFovID0WkaaAFu7/U7wVmAuuAqaq6RkQeEZGxACIySERygRuBf4vIGvf07sBSEVmJ86X/qN+sggeA+0UkB2eMwHOBxmROu3L9F+yv24j5Hc4L3YP07QuZmbaWgDHxTBXmz4+IrgDAkgaVUWV3gKqeLyJdgO/hfDEvBl5Q1VkBnDsDmFFm20N+t5fgNOmXPW8+UG4PtapuwZl5YGoo5UQRo3IW8UG3CyhJDGGPkIgzQPDhhyEvD9JtDKcxcefrryE/n+U9rvY6klNWpGUycssyGhUd4XBKA6/D8VRAswNUdRPwIM6v8IuAv4nIehG5LpTBmdC4aOtyGhYfY3oouwJ8JkxwfgnYLAFj4tOCBQCerBxYkRVp3UhA6VOwyetQPBfImIA+IvI4TpP+JcDVqtrdvf14iOMzIXDFhi/ZV68xCzr0Cf2Dde0KaWkwb17oH8sYE3nmz4f69VnfOsPrSE5ZmdaVUoT+eeu8DsVzgbQF/wN4BvhfVT3m2+hO33swZJGZ0FBl+PZVzM3oz8mExNA/nggMG3bq14AxJs4sWACDBoXn8yZAh1MakNOivQ0OJLDugCuAV30VABFJEJH6AKr6ciiDM8GXfmgPrY/sD+8gnWHDYOtW2LUrfI9pjPHesWOwYgUMH+51JGdxBgducLor41gglYDZgP9E8vruNhOFfM1fK8K5ktewYc61tQYYE1+WLnWmCfs+AyLI8vRuNDt+mI77A05fE5MCqQTUVdVC3x33dv3QhWRCqV/+Bo4lpbA+nEky+veHOnWsEmBMvPG954cO9TaOcqxIywSI+y6BQCoBR0Skv++OiAwAjlVyvIlg/fPXs6ptl9BODSyrbl2nImCVAGPiy/z50KULtGrldSRnyWnRnsPJ9awSEMAxPwHeFJEvROQL4A2cJEAmyqSUFNNj15ZTNeCwGjbMaRo8cSL8j22MCT9Vp+IfgV0BAKUJiaxs2zXukwZVWQlwE/p0A+4Gfgh0V9VloQ7MBF/PnZtJLi05lTs7rIYNcwYJrVwZ/sc2xoTf1q2we3dEDgr0WZHWje67t1L3xHGvQ/FMoEsJDwL64Czic5OI3Bq6kEyo9M/3YFCgjw0ONCa+uIsGRWpLADjjApK0lN47c7wOxTOBJAt6GfgzcD5OZWAQMDDEcZkQ6J+3nq+9WsSjfXsnbbBVAoyJDwsWQKNG0LOn15FUKNvtGu0fx10CgYwOGwj0cFfwM1GsX/4GFrev6aKQQWBJg4yJH/Pnw5AhkBg5SYLK2le/CVubtY3rwYGBdAd8BbQJdSAmtBru3UPbwr3eruQ1bBhs2wY7d3oXgwk6ERkjIhtEJEdEJpezP0VE3nD3LxKRjn77+ojIAhFZIyKrRaRuOGM3IVJYCKtWRXRXgM+KtG70z1sft0mDAqkEtATWishMEZnmu4Q6MBNcbTY5NV1PZgb42LiAmCMiicCTwOVAD5wxQz3KHHYHsF9VO+OsN/IH99wk4D/AXaraExgJ2PSRWLBkCZSWRkklIJPWR/aTdniP16F4IpDugF+FOggTem1yNnI8KZl1Xi7i0b8/JCc7lYBrr/UuDhNMg4Ecd4lvROR14Bpgrd8x13D6c+Qt4B8iIsBoYJWqrgRQ1b3hCtqEmG9QYAQmCSrLN1C6X94G8hu39jia8KuyEqCqc0WkA9BFVWe76wYE1MkjImOAJ9zjn1XVR8vsvxD4K87Mgwmq+pa7vS/wL6AxcBL4raq+4e57EWc544NuMRNVNTuQeOJZm5wNrGrTmROJdbwLIiXFkgbFnnTga7/7ucCQio5R1RIROQi0ALoCKiIzgVbA66r6x7IPICJ3AncCpKamkpWVVeugCwsLg1JOpImU59X7ww+p26EDS/ymBE/qXVKrMlPr1b6M8iR0b0fJq8l8v3gd3Xqf3XIRjtfTy79blZUAEfkBzhuwOdAJ5w39FHBpFef5mglH4XwwLBGRaarq/wthBzAR+J8ypx8FblXVTSKSBiwTkZmqesDd/1NfhcFULbnkBK23bWZa/7Feh+I0D/7rX1Bc7LQKmGgn5Wwr27la0TFJnJ51dBSYIyLLVHXOGQeqPg08DTBw4EAdOXJkbWMmKyuLYJQTaSLieanCxo0wbtwZsUycPL1WxU7qXcJjq0OR6TSJQamdSVy1sdzyt90yMgSPeSYv/26BvKL34DT5LQJwv5gDaTOpsplQVbe5+0r9T1TVjX6380VkN84vhQOYauu5azOJJeFLEtSxkjf7FVtS+Ofx407SoEGDwhKPCalcoL3f/XZA2RVZfMfkuuMAmgD73O1zVfUbABGZAfQH5mCi18aNsG9fRCcJKmtFWjcmLvuA5JITFCd52FrqgUAGBhaparHvjvsmDmQYZXnNhOnVCw9EZDCQDGz22/xbEVklIo+LSEp1y4w3vjmwng4KdJ2qiFiXQKxYAnQRkQwRSQYmAGUHDk8DbnNv3wB86k45ngn0EZH67ufKRZw5lsBEI997OwoGBfqsSMsk5eQJeuze4nUoYRdIS8BcEflfoJ6IjMJJHfxBAOcF0kxYeQEibYGXgdtU1dda8HNgJ07F4GngAeCRcs4Naj9ipPS11cSEwrUcadmK7w5vAgS/T616mnL4jRYcff991vXpE/BZ0fz6Q/THXxG3j/9enC/0ROB5VV0jIo8AS1V1GvAc8LKI5OC0AExwz90vIn/BqUgoMENVa9dmbLy3eDE0bgzdPJyOXE2+Hyf98tefSiAULwKpBEzGmeKzGvgvYAbwbADnBdJMWCERaQxMBx5U1YW+7apa4N4sEpEXOHs8ge+4oPYjRkRfWw3lrdtEXs/MEPWnVV9Gq+5ctXkzqdV4PaP59Yfoj78yqjoD53PBf9tDfrePAzdWcO5/cKYJmlixeLHT1ZcQaFZ67+1u1IK8Rq3ol7+BF7wOJswCmR1QCjzjXqrjVDMhkIdT+785kBPdZsV3gZdU9c0y+9qqaoE7xWgcTjIjU5G8PNIP7+Hzzld7HckpK9K7cdWnX0JBAbRt63U4xpha8o0DSjlRxFcrsnl68HX8qZYDAcNtRVpmXGYODGTtgK0isqXsparzVLUEZ8nhmcA6YKqvmVBExrplDxKRXJxfCf8WkTXu6d8GLgQmiki2e+nr7ntFRFbjtEy0BH5TzeccXxY6jSg7O0dO05yNCzAmNvXcvYU6pSdZ2bar16FU24q0TNof3EWrwv1ehxJWga4d4FMX5wu7eSCFB9BMuASnm6DseRU2EarqJYE8tnEtWEBRYh32dMxwqmIRYE1qp9NJg667zutwjDFBcl6BM7ErEgYhV5cvaVDfgg3M6hL5SY6CpcqWAFXd63fJU9W/AvZFHC0WLOCr1E6URtC0l+KkOjBggLUEGBNj+uZvJL9RS/Y0DOh3YkRZk3ouJxIS6RtnXQKBJAvq73c3AadloFHIIjLBU1wMy5axvM/lXkdytmHD4MknLWmQMTHkvIKNUdkVAFBUJ4X1rTrStyC+KgGBDN98zO/ye2AATp+9iXTZ2VBUFLYkQdUybBgUFTlJg4wxUa/psUN0PFAQtZUAgOy0TPoUbCKh9KTXoYRNIN0BF/tdRqnqD1Q1vqpK0cptbl8RqZUAsC4BY2JE33xnPEB2WhRXAtpm0qj4GJ325nodStgE0h1wf2X7VfUvwQvHBNXChdCuHTsbt8T7JEFlpKdD+/ZOJeBHP/I6GmNMLZ1XsJFShNWpnb0OpcZ8FZi+BRvY1KqDx9GERyDdAQOBu3FS/qYDd+GsG94IGxsQ2RYsiOzUncOGWUuAMTHivIKNbGrZniMp9b0Opca2NE/nUEoD+uVvrPrgGBFIJaAl0F9VJ6nqJJwxAe1U9WFVfTi04ZkaKyiA7dsjez3voUOdGAsKqj7WGBO5VDmvYCPZbaNvaqA/lQSy23aNq8GBgVQCzgGK/e4XAx1DEo0JHjdJUMS3BIC1BhgT5dof3EWLY4dYGcXjAXyy23Ylc8926hUf9zqUsAikEvAysFhEfiUiv8RZUvil0IZlam3BAmfqXf/+VR/rlX79TicNMsZELd/c+mieGeCTnZZJopbSe1eO16GERSCzA34L3A7sBw4At6vq70IdmKmlBQucL9mUCF5pOSXFkgYZEwPOK9jIsaQUNrSM/sF0vlUE4yVpUKDLPNUHDqnqE0CuuyiQiVQnTsDSpZHdFeAzbJgTa3Fx1ccaYyJS34KNfJXaiZLEyFiptDb21W/CjiapVgnwcbsAHgB+7m6qgy39GdlWroTjx6OnElBU5CQ2MsZEnxMn6LVrMyvbdvE6kqDJTsvkvIJNXocRFoG0BFwLjAWOAKhqPjY1MLL5mtcjeWaAjw0ONCa6ffUVdUuKY2I8gE9220zSD++hVeE+r0MJuUAqAcWqqoACiEiD0IZkam3hQkhLc5LxRDr/pEHGmOizaBEQt56vuwAAIABJREFUnSsHVsSXNKhfHHQJBFIJmCoi/waaisgPgNnAM6ENy9TKggVOK4CI15EExpIGGRO9Fi9mb73G5DZJ9TqSoFmT2onihKS4yBdQ5SgOVf2ziIwCDgGZwEOqOivkkZma2bULtm6FH/7Q60gCN2wYTJ0K+flOC4YxJnosXux0BUTLj44AFCUls651Bn3zN9Jx8vQqj9/26JVhiCo0Kq0EiEgiMFNVLwOq/cUvImOAJ4BE4FlVfbTM/guBvwJ9gAmq+pbfvtuAB927v1HVKe72AcCLQD1gBvBjt7vCQHQkCYIz3lh980t5D7jrnn/wceaIM46L5jeXMTHv0CFYu5bsETd7HUnQZad15fqvPiWh9CSlCYlehxMylXYHqOpJ4KiINKluwW4F4kngcpy1Bm4SkR5lDtsBTAReLXNuc+CXwBBgMPBLEWnm7v4XcCfQxb2MqW5sMW3BAkhKiuwkQWWsST2XosQ69M9b73UoxpjqWLYMVGNqUKBPdttMGhYfo/Per70OJaQCmdR5HFgtIrNwZwgAqGpVS78NBnJUdQuAiLwOXAOs9Stjm7uvtMy53wJmqeo+d/8sYIyIZAGNVXWBu/0lYBzwUQDPIz4sXOgkCapXz+tIAnYisQ6r2nRhQN46r0MxxlTH4sUAMTU90Od00qCNbGzV0dtgQiiQSsB091Jd6YB/FSoX55d9Tc/1rWKYW872s4jInTgtBqSmppKVlRXgQ5evsLCw1mWEwuq8g6duJ5SUcOfChay96DI+f+X9U9sn9YbUejCpd4QtJ+wnoW9X+s78kJ9lHuVkcvKp7b7XPFJf/0BFe/zGlGvRIujcmQP1GnsdSdBtbZbGwZQG9C3YwNTzRnsdTshUWAkQkTmqeinQQ1UfqEHZ5Y0SCbTvvqJzAy5TVZ8GngYYOHCgjhw5MsCHLl9WVha1LSMUJvr1rfcu2Mo9RUU8U7cnH64+8087qXcJj62O3Gxeq1J68UzJe8yes53l7bqf2r7tlpFA5L7+gYr2+I0p1+LFcNFFXkcRGiKsbNs15qcJVjYmoK2IXASMFZF+ItLf/xJA2bmA/0T1dkB+gHFVdG6ue7smZca8gXlOT8vS9LJDLyLf8rRuwOnnYIyJcHl5zmXwYK8jCZkVaZl0/WYH9YuPeR1KyFRWCXgImIzzRfsX4DG/y58DKHsJ0EVEMkQkGZgATAswrpnAaBFp5g4IHI0zS6EAOCwiQ0VEgFuB9ysrKJ4MyFtPXqNW7Gzc0utQqm1vg6ZsbdbWxgUYEy2WLHGuY7gScGpFwZ2xu6JghZUAVX1LVS8H/qiqF5e5XFJVwapaAtyL84W+DpiqqmtE5BERGQsgIoNEJBe4Efi3iKxxz90H/BqnIrEEeMQ3SBC4G3gWyAE2Y4MCHaoMyF3LMr+m9GizLL2HM0PAZnwaE/kWLXJmIvXr53UkIeOb9RDLSYMCSRb065oWrqozcOby+297yO/2Es5s3vc/7nng+XK2LwV61TSmWJV+aA9tC/eyND2aKwHdueGrOXQ4UMD2ZpY0yJiItngxnHce1K3rdSQhs69+E7Y3bUPf/I1ehxIygS4lbCLcALcvfVlUVwKccQHWJWBMhCstdboDYrgrwCe7bWZMLytslYAYMSBvHUfq1GV96wyvQ6mxTS3P4VBKA6sEGBPp1q+Hw4dhSKCzvqNXdlpX2hbuJfXwN16HEhIBVQJE5HwRud293UpEovebJkYNzF3HirRMTkZxekuVBJandWNArlUCjIlobpKgeGkJAGK2S6DKSoCI/BJ4APi5u6kO8J9QBmWqp0HRUbrt2cayKJwaWNay9G50/WYHjY8Xeh2KMaYiixdD48aQGTvLB1dkbeq5FCck0S9GBwcG0hJwLTAWN2WwquYDjUIZlKmevgUbSdTSU33q0WxZencS0Jjug4s1IjJGRDaISI6ITC5nf4qIvOHuXyQiHcvsP0dECkXkf8IVs6mlRYtg0CBIiP0e5aKkZNamZnBeQZy2BADF7ip9CiAiDUIbkqmugblrKUVYEQOVgOy0TE5Kgo0LiBIBLhR2B7BfVTsDjwN/KLP/cWyqb/Q4dgxWrYqLrgCf7LaZ9N6ZQ0LpSa9DCbpAKgFTReTfQFMR+QEwG3gmtGGZ6hiQt44NrTpwOCX662dHk+uxrnWGVQKix6mFwlS1GPAtFObvGmCKe/st4FI32RciMg7YAqwJU7ymtrKzoaQkvioBac6Kgl1icEXBKisBqvpnnDfu20Am8JCq/j3UgZnAJJSepF/+epa2i/7xAD7L0rs5XRwxWOuOQRUt9lXuMW4SsYNAC7dV8QHg4TDEaYJl0SLnOp4qAb6kQTHYTVllsiAR+W/gTVWdFYZ4TDVlfrOdRsXHojo/QFnL0ntw2/LpdNuzzetQTNUCWdSromMeBh5X1UK3YaD8BwjyiqAQu6s6huN5df/gA5q0asXCjRtho9NPHo4VSj1dCVVbc+y1Rtx6fD3tel961u5oXqU2kGXlGgMzRWQfTlPfW6q6K7RhmUD5ptNFc6bAsnwVmv7WJRANAlkozHdMrogkAU2AfThLi98gIn8EmgKlInJcVf/hf3KwVwSF2F3VMSzP6/vfhwsuOONx/FczDRWvV0Lt3aorqWs2lRuDb7XTmvLy/zGQtMEPAw+LSB9gPDBXRHJV9bKQR2eqNDBvLbsaNie3SarXoQRNXuNW7GzY3MYFRIdTC4UBeTgLhd1c5phpwG3AAuAG4FN3sPEFvgNE5FdAYdkKgPFeR78v+GZHD7Ji82YebXcBT4Xhiz+SZKd15b75b1C/+BhHk+t5HU7QVGd+x25gJ7AXaB2acEx1Dchbz7K0blBJc2rUEWFpeg8GWtKgiBfIQmHAczhjAHKA+3FWJzVR6LyCTYDzhRhvsts6Kwr22bnJ61CCKpAxAXfjtAC0whkg+ANVtUXfI0F+Pu0P7uLFAVd7HUnQLU/vzlUbvnTWK08vO87MRJIAFgo7jrNSaGVl/CokwZmgOq9gI6UIq1M7ex1K2PkqPn3zN7LwnD4eRxM8gXSwdAB+oqrZoQ7GVNO8eUBsjQfwOZX4aMECuOEGb4MxxgAw5OuvWN+6I0dS6nsdStgdqNeYbU3bxtyywhV2B4hIY/fmH4EdItLc/xKe8Eyl5s3jWFIKa1I7eR1J0K1J7cSxpJRTFR1jjLfqnjjOgLy1zOtwnteheCY7rWvMTROsbEzAq+71MmCpe73M736VappOVERuEZFsv0upiPR192W5Zfr2xe/4hPnzWdW2CyWJ3o2YDZWSxCRWte0C8+d7HYoxBmeRspSTJczr0NfrUDyT3TaTNoX7aHModlYUrLASoKpXudcZqnque+27nFtVwbVJJ6qqr6hqX1XtC3wX2FamO+IW335V3V2N5xs7jh6FFStisivAZ1l6N1i+3ElTaozx1IjtKzmRkMji9j29DsUz2WnuioIx1CUQyCqCcwLZVo5apRP1cxPwWgCPF18WL4aSkpjKFFjW0vQeTnrSpQE1PBljQmj49pWsSMuMqelx1bW29bkcT0pm8Nexk+W6sjEBdd2+/5Yi0sxvPEBHIC2AsmucTrTMMeM5uxLwgtsV8H/lVBrig9tXvjwt+hcNqsgKt9Zt4wKM8VaTY4fpvTMnrrsCAIqT6rC4XU/O3xY74+Qr60z+L+AnOF/4yzid+vMQTjN/VWqTTtTZKTIEOKqqX/ntv0VV80SkEc56Bt8FXjrrwYOcajTS0oz2/uAD6nbowPcH1wOqTqXpacrNGmvA0ffac/SDDyjs1SuiXv/qirT/H2OqY+iO1SSgzOsYv4MCfT7P6MeDnz1Pm0PfsLNxS6/DqbUKKwGq+gTwhIjcV8MFg2qTTtRnAmVaAVQ1z70+LCKv4nQ7nFUJCHaq0YhKM1paChs2wA03BJxG0+uUmzV136hR1H//fRo2aBA5r38NRNT/jzHVNGL7SgqT65HdNtPrUDz3RUZ/+Ox5Lti2gjf7jPI6nFoLZBXBv4tILxH5tojc6rsEUPapdKIikozzhT6tzDG+dKJwZjpRRCQBJ8HI676DRSRJRFq6t+sAVwFfEW/WrYMDB2DECK8jCb3hw2HvXup9HXtLeBoTLUZsz2Zxu54xOROpuja07MDuBs24YNsKr0MJikAGBv4S+Lt7uRgnb8DYSk8iKOlELwRyVXWL37YUnMWMVgHZOLnKn6kqlpjjmzY3fLi3cYSD+xybrImdgTjGRJM2h76h0768uM4PcAYRvujYlxHbshEt9TqaWgukWncDcB6wQlVvF5FU4NlACq9NOlFVzQKGltl2BBgQyGPHtHnzoFUr6NIFiK081mfJzIRmzWhslQBjPDFi+0oA5nWM70GB/r7I6M/1az6jx64trGkT3SmUA1lA6JiqlgIlbhbB3UCVeQJMCM2b5/xCjoeJEQkJMHw4Tb6Kv14fYyLBiO3ZfFO/CRtadfA6lIjhmyVxYQx0CQRSCVgqIk1xmt2XAcuBxSGNylRs1y7IyYmP8QA+w4fTYPt22Lev6mONMcGjyojtK1lwTh9UqrPobGzb07AZa1tncMHWOKgEqOoPVfWAqj4FjAJuU9XbQx+aKdeCBc51PIwH8PE914ULvY3DmHizbh2phfv40roCzvJ5x34MzF1LveLjXodSK5UlC+pf9gI0B5Lc28YL8+ZBcjIMiKOhEYMHowkJljTImHCb4ySHtUGBZ/uyYz+SS0sY8nV0d1VWNjDwsUr2KXBJkGMxgZg3DwYOhLp1vY4kfOrX53CXLjS2xYSMCa85c9jRJJXcpm28jiTiLGnXg+NJyVE/VbCyZEEXhzMQE4Djx2HZMvjRj7yOJOwO9exJ448/hhMnoE4dr8MxJvaVlMBnn/FlxjCvI4lIRXVSWNyuZ9SPC6hyiqCI1MeZw3+Oqt4pIl2ATFX9MOTRxbmOk6efcX9A7lreLi7mzs0pfFJmX6w72KsX7d55B1audFpCjDGhtWwZHDrEfOsKqJAvhTC5udCundfh1Eggwz1fAIoB30i0XOA3IYvIVGhg3loAlsXw8sEVOdSrl3PDugSMCQ93PIBVAir2RYY7PG7WLG8DqYVAKgGdVPWPwAkAVT1G+Qv/mBAbmLuOLc3S2NugqdehhF1Rq1bQvr1VAowJl9mz4bzz2Fe/ideRRCxfCuFYrwQUi0g93NX9RKQTUBTSqMzZVOmft47lcdgKcMr558Pcuc4CSsaY0Dl2zKlwX3qp15FENjeFMLNmRe3nUiCVgF8CHwPtReQVYA7ws5BGZc7SaV8uLY4dYmk8VwKuvBJ27oQlS7yOxJjYNm8eFBXBZZd5HUnE+yKjP3zzDWRnex1KjVRaCRARAdYD1wETcZb1Hejm9TdhdGmOk6Txy4x+HkfioSuugKQkeO89ryMxJrbNnu281y64wOtIIp4vhTCffOJtIDVU6ewAVVUReU9VBwDxNRw9wnxr4wK+Su1EbpNUr0MJK98MiUm9S5j48XxeTu9F26f/w2V6/qljtj16pVfhGROb5syBoUOhYUOvI4l4exo2g/POcyoBkydXfUKECaQ7YKGIDAp5JKZCrQr3MSB/PR93tfm6n3QdSud9uXTa+7XXoRgTm/bvd6YHWldA4EaNcrpQjhzxOpJqC6QScDGwQEQ2i8gqEVktIqtCHZg57VubnJz5M7tYJWB25yEAjNq0yONIjIlRn30GqjYosDpGj4biYvj8c68jqbZAKgGXA51w0gRfDVzlXldJRMaIyAYRyRGRs9pJRCRFRN5w9y8SkY7u9o4ickxEst3LU37nDHArIjki8jd33EJMG71xAZubp7Op5Tleh+K5gsatWNmmC6M3LfA6FGNi05w5TjfAkCFeRxI9zj/fSeUeheMCAllFcHt5l6rOE5FE4EmcSkQP4CYR6VHmsDuA/araGXgc+IPfvs2q2te93OW3/V/AnUAX9zKmqliiWePjhQzbsYpPugyD2K/vBOSTLkPpn7+BVoW2tLAxQTdnDlx4oaXnro569ZzXLBYrAbUwGMhR1S2qWgy8DlxT5phrgCnu7beASyv7ZS8ibYHGqrpAVRV4CRgX/NAjx6U5i6lTepKZNh7glE+6DAVgVI51CRgTVLm5sGGDdQXUxKhRsHYt5OV5HUm1hLISkA74j97KdbeVe4yqlgAHgRbuvgwRWSEic0XkAr/jc6soM6Z8a9MCChq2YGXbLl6HEjE2tTyHrc3aMnrjQq9DMSa2uKmCbVBgDYwe7VxHWfbAKhcQqoXyftFrgMcU4CxYtFdEBgDviUjPAMt0Cha5E6fbgNTUVLKysgKNu1yFhYW1LqO6Huh6hEu3LWPdRZdxf59SoOYZqVLrOdPsolXZ+PcPG8IFMz9kcqdDYf+71IQX/z/GVNucOdCqFfjW6jCB690bUlOdLoGJE72OJmChrATkAu397rcD8is4JldEkoAmwD63qb8IQFWXichmoKt7vP9STeWViXve08DTAAMHDtSRI0fW6slkZWVR2zKq67W//YI6xcX8vsUIFqyu3Z9qUu8SHqtlGV4qG39Ws+G8ffI91nywgrum3eRhZIHx4v/HmMqUXaUUVRa9O53F7Xpy3/9+5E1Q0UzEaQ346CMnhXBCKBvagyeUUS4BuohIhogkAxOAaWWOmQbc5t6+AfjUTVDUyh1YiIicizMAcIuqFgCHRWSoO3bgVuD9ED4HT31r4wL2123E4vZWKy9rRVome+o3ZfQm6xIwJhg67c0ltXAfX3bs63Uo0Wv06KhLIRyySoDbx38vMBNYB0xV1TUi8oiIjHUPew5oISI5wP2AbxrhhcAqEVmJM2DwLlX1DQW/G3gWyAE2A7FZZT1xgktzFjOn82BOJiR6HU3EKU1IZHbnwYzcstTJcW48U4upwKNEZJk75XeZiFwS7tjNaedvd7645tnSwTXnG0sRReMCQto+rKozgBlltj3kd/s4cGM5570NvF1BmUuB2P9pnJVFk6IjfNx1uNeRRKxPug7jplWfOMlNxsT0TNGI5TcVeBROd90SEZmmqmv9Djs1FVhEJuBMBR4PfANcrar5ItIL5wdDTA/0jWQjtq9ke9M25DZt43Uo0atNG+jTxxkX8MADXkcTkOjotIhH77zDkTp1nWUqTbnmdziPI3Xq2oJC3qrxVGBVXaGqvjE9a4C6IpISlqjNGZoeO8RFW5bx2bkDvQ4l+o0eDV9+GTUphK0SEIlKS+G995ib0Z+iOvaZWJGipGSyzh0A778ftWt5x4DaTgX2uR5YoarWt+OB61fPIeXkCV7ray1qtRZlKYSjd7h4LFu4EHbu5ONB3/E6koj3SZehXPnhPFi82Fn1zIRbbaYCOzud6b9/AEaX+wBBnu4LsTtlszrP69SUW1W+M+UjCrpkctUl7biKyJxKHMnTnP1f84TSUkYkJ1Pw3HPk1KsX0Ple/j9aJSASvfsu1KnDZ51s8caqfNZpkLPu+XvvWSXAGzWeCgwgIu2Ad4FbVXVzeQ8Q7Om+ELtTNqvzvCa6UwSHbV/FfTvzuX/Af/NOBE8jjuRpzttuGXnmhjFjaPfFF7R75RVIqbo118v/R+sOiDSqTiXg0ks5nNLA62gi3qG6DeHii21cgHdqMxW4KTAd+LmqzgtbxOYMN2d/xIG6DZmeeb7XocSOe+6B3bvhzTe9jqRKVgmINKtXw+bNcO21XkcSPcaNc/Kdr1/vdSRxp5ZTge8FOgP/57diaOswP4W41uLIAb61cQFv97rUxh8F06hR0K0b/P3vXkdSJasERJp333UyT11TdoC1qdBY97vGWgM8oaozVLWrqnZS1d+62x5S1Wnu7eOqeqOqdlbVwaq6xd3+G1Vt4LdaaF9V3e3lc4k3N66eTXJpCa+eZwMCg0oE7r3XGau0KLIXOrNKQKR5910YMcLJQW0C064dDBpklQBjqkG0lJtWfsyi9r3Y3LJ91SeY6rn1VmjcOOJbAyJzlEW82rIFVq6Exx7zOpLoM24c/OIXkJ8PaWleR2NMxBuxbSUdDuzksQu+63UoUe+sdRhcD3UZyXdee4MRTcewp2Eztj16ZZgjq5q1BESSd991rm08QPWNG+dcTys7Js0YU55bsj9ib73GlpU0hF7qfyXJpSXcnB252e2tEhBJ3nkH+vaFjAyvI4k+3btDly7WJWBMIPLzGbVpIW/2vozipDpeRxOztjVP59NzB3JL9kfUOXnC63DKZZWASLFzJyxYYK0ANSXitAZ8+ikcPOh1NMZEtuefJ0lLLUNgGEwZcDWtj+zn8g2ROQvWKgGR4v33nRwB113ndSTRa9w4OHHCWc/bGFO+kyfhmWf4okNftjez8TOh9nlGPzY3T+f2pR94HUq5rBIQKd55Bzp3hp49vY4keg0Z4syqsC4BYyo2cybs2MGr1goQFioJvNT/KvoVbIAlS7wO5yxWCYgEBw44zdjXXec0a5uaSUx0cgbMmAFHj3odjTGR6amnIDWVWV0szXa4vN3rUgqT60XkdEGrBESCd96BkhIbDxAMt90Ghw/DE094HYkxkefrr2H6dLjjDkoSbYZ4uBSm1OfN3pfB66/Drl1eh3OGkFYCRGSMiGwQkRwRmVzO/hQRecPdv0hEOrrbR4nIMhFZ7V5f4ndOlltmbKQZPXYMHn4Y+vWDwYO9jib6jRjhtAY8+ih8843X0RgTWZ57zhl79IMfeB1J3Hmp/1XOmKWnn/Y6lDOErBIgIonAk8DlQA/gJhHpUeawO4D9qtoZeBxnOVGAb4CrVbU3zsIjL5c575aYSTP6xBOwY4eTICjBGmaC4ve/h8JC+M1vvI7EmMhRUgLPPANjxkDHjl5HE3e2Nk93Xvt//QuKi70O55RQfusMBnJUdYuqFgOvA2UT4l8DTHFvvwVcKiKiqitU1bcc6RqgrojE3uoWu3fD734HV1/trIRngqNHD7jjDvjnP50sjMYYpxsgPx/+67+8jiR+3XcfFBQ4XcARIpSdQunA1373c4EhFR2jqiUichBogdMS4HM9sEJVi/y2vSAiJ4G3gd+oqgY7+FDrOHk6v/7kn9xUeIRvpV7B5grSTpqqlZeys3Wd85mrLzHrqtv50difRWS6TmPC6qmnID0drrT3gmfGjHFmgf3tbzBhgtfRAKGtBJQ3zL3sl3Wlx4hIT5wugtF++29R1TwRaYRTCfgu8NJZDy5yJ3AnQGpqKllZWdUKvqzCwsJal+HvN823cnP2x6y+bAzjRrYFSoJWdnlS68Gk3qF9jFCqfvxNWJ03lrHvv8mxm8aSldUgZLEFItj/P8ZUy7ZtztTAhx6CJBsQ6JmEBGd1wZ/8BJYtgwEDvI4opJWAXMB/aap2QH4Fx+SKSBLQBNgHICLtgHeBW1V1s+8EVc1zrw+LyKs43Q5nVQJU9WngaYCBAwfqyJEja/VksrKyqG0Z/j6946cUJtdjYrdb2L869G/KSb1LeCwMjxMqNYm/YcaNzK03k3OencKwB+/3dPplsP9/jKmWZ55x/v+//32vIzETJzqLnf397/Dii15HE9IxAUuALiKSISLJwASg7Oou03AG/gHcAHyqqioiTYHpwM9V9VSuRRFJEpGW7u06wFXAVyF8DqExezaXbFnKP4Z9m/31m3gdTcwqTKnPEyNuYtiO1ZZF0MSvQ4ecWQFXXeUsu2281aSJUxF47TVnXJjHQlYJUNUS4F5gJrAOmKqqa0TkEREZ6x72HNBCRHKA+wHfNMJ7gc7A/5WZCpgCzBSRVUA2kAc8E6rnEBInT8KkSXzdJJUpA672OpqY91rfMWxt1hYeeMB57Y2JJ6pw553OdNnJZ83SNl65915nhsAz3n99hbR9WFVnADPKbHvI7/Zx4MZyzvsNUNH8Lu87UWpjyhRYtYo/jP0ZRUnJXkcT804k1uFPF97GP99/FF56CW6/3euQjAk532DZm7M/4ncz3+CPF97KP9/fB+/bAOSI0K0bjBrlTBH/3vc8DcUmpodTYSE8+CAMHcqH3S7wOpq4MSNzhJOI6f/+z0nOZEwc6L57C7+c/TRzM/rzr6E3eB2OKeuxx5z05tddh3iYN8AqAeH05z87c0Qfe8zWCAgnEfjTnyAvz9IJm7jQoOgoT773KPvrNeK/r5qEin3UR5zevZ2BgQsX0vWJJ5yuGw/Yf0a45OU5X0Q33gjDh3sdTfy58EInKdPvf2/phE1sU+V3M5+kw4Gd/Gjsz9hng48j1w03wC9+QdsZM5xMgh6I3jlj0ebBB520nY8+6nUkcanj5Ol0bn05Mw9P58XLv8+vLz07d7olFDIx4dlnuWbdXP544a0sbt/L62hMVR55hG8+/ZSWP/4x9Orl/GAJI2sJCIfsbGdA4H33wbnneh1N3MppeQ5Te1/Gd5dPp92BnV6HY0zQNcjJgfvu4/OO/WwcQLRISGDd//4vdOrktAzs2BHehw/ro8UjVZg0CZo1cxJEGE89fv4tnExI5Kefl12Typgod/gwPR95BJo3t3EAUeZkw4bw/vtQVATjxjkDBsPE/ktC7YMP4NNP4Ve/cioCxlO7G7Xg2UHjuGbdXEZtWuh1OMYEhyrcdRf18vLgtdfY26Cp1xGZ6srMhFdfdVqOf/CDsA0UtEpAKE2d6iwS0b073HWX19EY11NDrmdVm8489e7v+N6S9z0blWtM0Dz3HLz6KtsmToSLLvI6GlNTV17pLIH+6qvOLLIwsEpAKJSWOgt1jB8P/frBZ59BnTpeR2VcR1LqM/6mR5nVZQgPffoMv/nknySdjN7FlUycW7XKGW80ahTbb77Z62hMbf38587YgAcecBZ9CjGbHRBkPe5/i8c/fIxvbVrI1N6X8eCweyh+fKnXYZkyjiXX5e5xP+dnc1/i7kVvcc6BnfCLi6CpNaOaKLJzJ3z7205X43/+A2vXeh2RqS0ReOEF2LDBaUlessRZfjhErBIQTFu38vZ/fkrXb3bw8KU/4IUBYy0pUARTSeAPIyeypXm3ArgHAAALsElEQVQav5v5pJO/4cMPbQaHiXwlJfDPfzpZMI8fd34xtm5tlYAI50vnXNak3iVMdPdte/RK8A0UHDgQrrkGFi6ERo1CEpN1BwTL3LkwaBBph/Yw8cZf8cLAa6wCECXe7DOa747/tfOrasgQmDev6pOM8cr8+TBoEPz4xzB0KKxeDbZMdezJyHDGlR075iSbCxGrBATDv/8Nl10GLVsy7ta/8EVGf68jMtW08Jw+Tm27aVO45BJ45RWvQzLmTHv2OIvNjBhB/uZc7hr3czr2vY+Oz2+i4+TpdJw8ndV5Byv8tWmi0KWXwvr1zoJDIWKVgNo4cQLuuccZ+T9qFCxaxNbm6V5HZWqqa1enIjB0KHznO/DLX9rMAeO9kyedlLJdu8LLL/PUkOu57Pv/4uPMEdbaGA+SQ7varI0JqIn1653+mtdeg5Ur4ac/dXLSJyZ6HZmpBd8vqDqD/5vf7UvmxkceIecfz/NR5gg+7jqMNamd2PaHqzyO0sSVJUvghz+EpUvh4ovhySd5dMoWr6MyMSSkLQEiMkZENohIjohMLmd/ioi84e5fJCId/fb93N2+QUS+FWiZIXHyJI1Xr+apIdezuUU7Z97/5Mms3lnIj67+KR0TLqLjLz62ZrgYcSKxDj+94sf8zxU/YVej5ty98E2mT/kJXz51B9x/P3z5pfPrzACheZ/HrT174L33nB8Ww4c7Y1Ty8px543PmOJ89xgRRyFoCRCQReBIYBeQCS0Rkmqr6D1+9A9ivqp1FZALwB2C8iPQAJgA9gTRgtoh0dc+pqszgOHoUZs1yfvF/+CH99+yhV0ISC8/pzYsDrmZ25yEUNG4V9Ic1EUKEt3pfxlu9L6PZ0YNclrOYMRvn0+7JJ+HxxyE11Unved11zoIfdet6HbEnQvE+V9X4qGGpOtPA5s1zLl9+CZs2AVCUmMTqNl34cvgEnh10LYUr68PKGR4HbMIhkB+TwVzsLJTdAYOBHFXdAiAirwPXAP4fDtcAv3JvvwX8Q0TE3f66qhYBW0Ukxy2PAMqsvc2bnbWejx2DJk3giitY07kzEw734nBKg6A+lIl8++s34c0+o3izzygaFh3l4s1L+NbGBVz8/BQa/PvfAByo25BdDZuzu0FzdjdqzvVXDoK0NGjbFtq2pd6OHbBtG6SkOH18vus6daK9XzcU7/MFYYq99lSd6XonTkBxsXN94gQcOAB79566/P4/X9Ls2GGaHjtM82OHaHrsEJ335tL82CEA9tVrzLL07iwdOZGl6T34qk1nipJC2xdsDIS2EpAOfO13PxcYUtExqloiIgeBFu72hWXO9Y24q6rMGjmj9qXKz867knkd+rK4fU9OJNZhUmYJh1fbEIp4V5hSnw96XMQHPS4i5UQRF2zLJnPPNlIL99H6yD5SD+8jY8dqiv/0Ocmlp7MQVvpP6l8pSEhwKgUip2+XvS5bafC/77t99dVOi0Xohep9HjwPPgivv+7cdgd6Djl+nO1HSxG/gZ8CiKpzwblO4PR9FBK0lDqlJ0kqLSHp5Mkz/saV+TnOr/sD9Rqzr15jDtRrxOzOg1ma3oNl7bqzuXm7aK8MmigVym+18v6jyw61ruiYiraXN4ah3OHbInIncKd7t1BENlQQZ7nuAVj41qn7P4KWwDfVKSOSWPyhsTHwQyuOv7jYuQTTX//qXKqvQzWPD8X7/MyTa/lerkD4/59OlkDhPucCwGpYPTuoDxGp75NgsOd2mvwhoMMCei+HshKQC7T3u98OyK/gmFwRSQKaAPuqOLeqMgFQ1aeBp2safFkislRVBwarvHCz+L0V7fFXIlTv81OC/V6G2P17xOrzAntuoRLK2QFLgC4ikiEiyTgDgKaVOWYacJt7+wbgU1VVd/sEd1RxBtAFWBxgmcaY8AnF+9wYEyYhawlw+/7uBWYCicDzqrpGRB4BlqrqNOA54GV3QNA+nA8Q3OOm4gwuKgHu8Y0YLq/MUD0HY0zlQvU+N8aEh6hlRAuIiNzpNktGJYvfW9Eef6yJ1b9HrD4vsOcWsse2SoAxxhgTn2ztAGOMMSZOWSWgHCLyvIjsFpGv/LY1F5FZIrLJvW7mZYyVqSD+X4lInohku5crvIyxIiLSXkQ+E5F1IrJGRH7sbo+K17+S+KPi9Y9lInKj+zcpFZGBZfZFffpiT1Kqh0i0fwZXJBI/36wSUL4XgTFltk0G5qhqF2COez9SvcjZ8QM8rqp93Uuk5iAtASapandgKHCPm142Wl7/iuKH6Hj9Y9lXwHXA5/4by6QvHgP8U5x0yFFDTqdvvhzoAdzk938XjV4kuj+DKxJxn29WCSiHqn6OM4rZ3zXAFPf2FGBcWIOqhgrijwqqWqCqy93bh4F1OFnkouL1ryR+4zFVXaeq5SUaOpW+WFW3Av5pyqPFqfTNqloM+NI3R6Vo/wyuSCR+vlklIHCpqloAzh8SaO1xPDVxr4iscpvaIr4pTZzV5voBi4jC179M/BBlr38cKS/1cbRV3GLhOVQl6j4DKhMpn29WCYgf/wI6AX2BAuAxb8OpnIg0BN4GfqKqh7yOp7rKiT+qXv9oJSKzReSrci6V/SoOKH1xhIuF5xA3IunzzVbECdwuEWmrqgUi0hbY7XVA1aGqu3y3ReQZ4EMPw6mUiNTBeYO8oqrvuJuj5vUvL/5oev2jmapeVoPTAkpfHOFi4TlUJWo+AyoTaZ9v1hIQOP/Up7cB73sYS7W5/1g+1+IMkoo4IiI4GebWqepf/HZFxetfUfzR8vrHqVhIXxwPKdWj4jOgMpH4+WbJgsohIq8BI3FWdtoF/BJ4D5gKnAPsAG5U1YgcfFdB/CNxmqIV2Ab8l68PKpKIyPnAF8BqoNTd/L84/WYR//pXEv9NRMHrH8tE5Frg70Ar4ACQrarfcvf9Avgezujtn6jqR54FWkPutNO/cjp98289DqnGov0zuCKR+PlmlQBjjDEmTll3gDHGGBOnrBJgzP+3dz+hUpVhHMe/v4UlSQS5sBsSEhLRIrMoQqM/lEEgbQqCSFzVTghyV4ugjVALlXDV1o20CKFI0BA0EMOKvKVBIEhUFlFQahfyPi3Oe3EY73CZAbvee74fGDhn5rzveQfmGR7ec877SFJPmQRIktRTJgGSJPWUSYAkST1lEiBJUk+ZBAiAJKsHytz+MlT29qbFHt+wJLuTPN62VyTZ1cpwTic5meS59tlh1+lXnxjLGofrBOgaSd4G/q6q9xZ7LPNJcjvwSVU92vZ3AVPAa1U1k2QN8ERVHUiyHVi7lBdOkSZlLGshzgRobEnWJTmb5IOWre9P8kySz1sG/0g7blWrmPdFkq/miri09seSfNlem9r7TyY5muTD1v/+tszmsBeBT1ubW4BXgR1VNQPdOv1VdaAde5ButT5JQ4xlmQRoUuuBPcD9wL3Ay8BjwE66ZTAB3gQ+q6qHgaeAd5OsoiuOsaWqHgReAvYO9LsReB24D7gb2DzPuTcDpwbGcX5UJa6q+gO4OcnqCb+ntNwZyz1mFUFN6lxVnQZI8i1wpKoqyWlgXTvmWeD5JDvb/kq6tbF/At5P8gBwBbhnoN+TVfVj6/fr1tfxoXNPAb+NMdZfgTuB38doI/WFsdxjJgGa1MzA9uzA/ixXf1cBXqiq7wcbtuuUF4ANdLNR/4zo9wrz/0Yv0/0JAfwA3JXk1qr6a8RYV7Y2kq5lLPeYlwN0PR0CdsxdC0yysb1/G/BzVc0C2+iqno3jDN3UIVV1ia405965O5+TTCV5pW0HuIOucp+kyRjLy5RJgK6nd4AVwDdJpts+wD5ge5ITdNOHF8fs92O6MqNz3qKbUvyunecjrk4xPgScqKp/J/oGksBYXrZ8RFBLUpLjwNaq+nOB4/YAB6vqyP8zMknjMJYXlzMBWqreoLsxaSHT/mlINzRjeRE5E6CR2qM48wXd01Xl3bnSEmEsaxSTAEmSesrLAZIk9ZRJgCRJPWUSIElST5kESJLUUyYBkiT11H8MUxDt14gkvgAAAABJRU5ErkJggg==\n",
"text/plain": [
"