Skip to content

Instantly share code, notes, and snippets.

@acherrera
Created August 17, 2019 19:40
Show Gist options
  • Save acherrera/cb03b52e1b4026251100a2d1c9ed4127 to your computer and use it in GitHub Desktop.
Save acherrera/cb03b52e1b4026251100a2d1c9ed4127 to your computer and use it in GitHub Desktop.
Analyzing some weather data to predict soil moisture
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Soil Moisture Analysis\n",
"\n",
"## Problem\n",
"\n",
"I mountain bike. A lot. But when it rains I can no longer mountain bike because I will tear up the trails and the mud is not good for my bike. This is a very common issue and there are many, many posts on Facebook asking if trails are dry, when they'll be dry, etc. People get tired of answering these questions, but this needs to be done or people might go out and ride trails that are too wet. \n",
"\n",
"## Solution\n",
"\n",
"Create a program that can predict how dry the trails are using weather input data. If there is a lot of moisture, we want that to evaporate faster than small amounts of moisture. Initially I'll just be taking whatever moisture the previous reading had, divided that by a number and add any rain the occured. Side note - taking the previous value and dividing by a value to get the current value is a first order differential equation in case you want sound like a stuck up prick. Bonus stuck up points - that is also (basically) Newton's Law of Cooling. Okay, let's get started\n",
"\n",
"\n",
"## Data Source\n",
"The trails I ride often are right next to the Des Moines aiport (KDSM) so I will take the weather readings from there. In case you want to follow along at home: https://mesonet.agron.iastate.edu/request/download.phtml\n",
"\n",
"If you want to use a different airport, um, good luck? The prodcedure should be very similar and issues will be similar as well, but I don't know of a good source for long term historical data for every airport on Earth off the top of my head. ASOS reporting is standarized across all airports so a different location should be pretty simple to integrate once you find the raw data"
]
},
{
"cell_type": "code",
"execution_count": 69,
"metadata": {},
"outputs": [],
"source": [
"# Imports\n",
"import pandas as pd \n",
"import matplotlib.pyplot as plt\n",
"from datetime import datetime, timedelta"
]
},
{
"cell_type": "code",
"execution_count": 70,
"metadata": {},
"outputs": [],
"source": [
"# Bring in that data\n",
"file_name = './raw_data/kdsm_2019_partial.txt'\n",
"asdf = pd.read_csv(file_name)"
]
},
{
"cell_type": "code",
"execution_count": 71,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Index(['station', 'valid', 'tmpf', 'dwpf', 'relh', 'drct', 'sknt', 'p01i',\n",
" 'alti', 'mslp', 'vsby', 'gust', 'skyc1', 'skyc2', 'skyc3', 'skyc4',\n",
" 'skyl1', 'skyl2', 'skyl3', 'skyl4', 'wxcodes', 'ice_accretion_1hr',\n",
" 'ice_accretion_3hr', 'ice_accretion_6hr', 'peak_wind_gust',\n",
" 'peak_wind_drct', 'peak_wind_time', 'feel', 'metar'],\n",
" dtype='object')"
]
},
"execution_count": 71,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"asdf.columns # Check out the data "
]
},
{
"cell_type": "code",
"execution_count": 72,
"metadata": {},
"outputs": [],
"source": [
"asdf.valid = pd.to_datetime(asdf.valid) # convert to datetime object - MUCH better for plotting"
]
},
{
"cell_type": "code",
"execution_count": 73,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"11 0.0001\n",
"17 0.0001\n",
"20 0.0001\n",
"23 0.0000\n",
"25 0.0001\n",
" ... \n",
"69472 0.0000\n",
"69485 0.0000\n",
"69497 0.0000\n",
"69510 0.0000\n",
"69522 0.0000\n",
"Name: p01i, Length: 9329, dtype: float64"
]
},
"execution_count": 73,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"asdf[asdf.p01i.notnull()].p01i # Get rid of all the 'NaN' - just checking"
]
},
{
"cell_type": "code",
"execution_count": 74,
"metadata": {},
"outputs": [],
"source": [
"rain_vals = asdf[asdf.p01i.notnull()] # Actually get ride of 'NaN"
]
},
{
"cell_type": "code",
"execution_count": 75,
"metadata": {},
"outputs": [],
"source": [
"# Because pandas will not to import so much crap in the future. Future proofing code\n",
"# This is based off the warning that comes if you plot dates\n",
"from pandas.plotting import register_matplotlib_converters\n",
"register_matplotlib_converters()"
]
},
{
"cell_type": "code",
"execution_count": 76,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Text(0, 0.5, 'Rain in Inches')"
]
},
"execution_count": 76,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAE8CAYAAAAMvfwgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3de5hcVZnv8e+vu5OA3CFBISEEBFEEBAnoceQiCoRBCCpoELl5iZdhdMZREYdBRXR0nCPqAZWIoCAICAphiKIMgjojmgARTBCJgCRcwzXhFtLd7/ljrUoqle5OdVO7qnbX7/M89XTVvr1rVXfXW2uvtddWRGBmZtZIXa0ugJmZjT5OLmZm1nBOLmZm1nBOLmZm1nBOLmZm1nBOLmZm1nBOLjZikiZLelpSd6vLYmbtxcmlw0m6V9JzOUk8JOn7kjasZ9+IuC8iNoyIvhHE3V9Sf45beVw9/BoMeOwTJP22EceqOua9kt4yyLpjqurwXG29GlyO9SSFpGfy8R+V9EtJb29knGGW6SZJj0jqqVl+iaRTh9ivUpdJVctOlbRE0k6SptW8l4sl/UjSHjXHOVLSbZKWSVoq6brKMSV9Ocf4YM0+J+fln27Mu2C1nFwM4LCI2BDYHdgDOKVJcR/IyanyOGy4B6j9QGuFiLioUgfgEGrqVVDYnfKxXwX8CDhX0skFxRqUpJ2AvYAxpLq/mGOdAbwf2Dci7syL78713Bh4A3AP8L+S9sn77AycC5wEbAK8HJgF9Fcd+i/AcTXhjsvLrSBOLrZKRDwEXEtKMgBIOlTSrflb4WJJn6taNyV/++vJr2+Q9AVJ/yNpuaRfSBo/3HJIGifp65IeyI+vSxqX1+2fv9meLOkh4PxhHvtESXfk8t1d/Y1W0nhJ/yXpSUmPS/qNpC5JFwKTgavzN+hPjaBO20i6Krc07pb0oap1X5Z0Uf5Wvjx/C999qONVRMTSiDgP+CjwWUkb52N+UNKf8/EWSXpvVbxFkg6ser2epKckvUrSBrnF8Xh+H34vabMhinA8cAMpwR1fdcyPAu8A/i2/Zz9ex/vzn8AMUmK5e4B69kfE4oj4DHAR8O951WuBP0fEryNZFhGXRcQDVbv/FthK0stzrKnAC8DtQ5XJXhwnF1sln0o4BFhUtfgZ0re8TYFDgQ9LOmKIw7wbOBHYEhgLfGIERflX4PWkJPcaYG+g+vTKy4DNgW2BmcM89iPAW0nfhE8EzpT02rzuX4AlwATgpcBngIiIY4H7yC28iPiP4QRU6pOaA/wvsDUwDfiMpP2qNnsbcB7pff5v4OvDrNdPgfWBPfPrB0m/y42BDwFnS3p1XncB8J6qfacDf4mIO0gthx5gIjCe1CJ4YZB6dQHHkj7sLwIOqySiiPgmcAXwhfyeHTVE2b9O+tvaNyLuq6OuPwFeL2kMMA/YQ9JX8xePDQbYPoAfsrr1chzpPbACObkYwJWSlgOLSR++n62siIgbIuL2/M3xNtI31P0GOQ7A+RHxl4h4DriMqlbQALbO344rj3fm5ccAp0fEIxGxFPg86UOsoh/4bESsyHHqFhHXRMRf87fcG4FfAPvk1SuBrYBtI2JlRPwmGjP53huB9SLiKxHxQkT8hdTimlG1zfUR8cvcf3UhQ79va4mIZ4CnSEmXiJgdEffkel4H3JjLAemD9QhJ6+fXx+aYkN6DCcDLI6I3IubmYw/kgLztFaTE+UBNnep1EPBfNa2NoTwAdAMbR8SfgTcD2wOXA49KOreqbhUXAO+RNBY4Erh4BOW0YXByMYAjImIjYH/glaRvrABIep2kX+WO0qdI34KHOtX1UNXzZ4Gh+hweiIhNqx6X5eVbA3+r2u5veVnF0oh4fp21GoCkQ5Q6oB+X9CTw96yuz1dJrbZf5FNXjers3RaYUp1IgY+TWmAVw3nf1pK/sW8CPJ5fHy7pD1X1PIBcz4i4F7iVlGAm5HWX5EN9j5SILs+nH7+kwUcDHg9cExFP5SS8xqmxYTgKOGEY7/dEoA9Yluvz24h4R0SMz3U5GFjj1GVELAIeBr4I3BoRD4+gnDYMLe8MtfYRETdK+j7wn0Dl1NfFwFnAIRHxvKSvM3RyaYQHSB/IC/LryXnZqqKO5KC53+YK0mmRqyJipaQrAQFExHLSqbF/kbQLcL2kuRHx3yONmS0m9Qvs+iKOsS5vA54Dbs6J5sekb+g/i4heST8n1zP7AenU2HhSq+kRgIhYAZwGnCZpe1If3ALSaa9VJG0EvB3oz31fAOOATSXtlDvk633PFgBvIb3fz0fEuk4Jvg24KSJW1q6IiN9JugrYZYD9LgC+BRxdZ7nsRXDLxWp9HThQ0mvy642Ax3Ni2ZvUp1K0HwGnSpqQBwScRjpnPhzKHdWrHqQ+oHHAUqBX0iGkUzKVHd4qaQdJIp1i6mP1qKOHSadeRuK3+fj/lMvSI2m3qr6eEZO0haTjSb+3MyJiGanvZQzpFGe/pMNJrdJql5NOk32Yqv4HSW+RtHPuT1kG9LLmyKuKI4GnSS3d3fPjVcAfWN23Ufd7FhF/JLU4Plc92KGqXJI0SdIXSEnxX/PyN0l6b26BkfuVDgVuGiDMhaTf95X1lMleHCcXW0Pu47iA9IEO8BHg9NwncxqpH6VoZ5A6am8jjei5JS8bjjeQvsnXPj5KqsMTpEQ5u2qfHYHrSB+avwO+FRG/yuv+nZTwnpQ0rEEK+Rv23+cy/Y2U3L7NME991bhT6RqayjDbD0fEl3K8R0kDKa4GHiO1QufUlGl5Xr81a74HE4GrgOXAn/J+lw4Q/3jg3Ii4PyIeqjyAs4Fjc3KaBeyV37NLBjjGGiJiHmkQwlcknZgXb5/r+TTwe2An4I25vwzS7/FIYEHe7mpSK2ut1k9EPBMR1+XWmRVMjemvNLOykfQlYMuIeH+ry2Kjj/tczDpQPo10Aqv71swayqfFzDqMpJOAe4EfR8QfWlwcG6V8WszMzBrOLRczM2u4juhzGT9+fEyZMqXVxTAzK5Wbb7750YiYMJJ9C00ukqYB3yBN1XBuRHy5Zv2+pCGDuwEzIuLyvPxNwJlVm74yr78yX+S3H+k6BIATImL+UOWYMmUK8+bNa0CNzMw6h6S/rXurgRWWXPKUEWcDB5ImA5wraXZELKza7D7SiJU1rhvI1xbsno+zOXlKjqpNPllJRGZm1n6KbLnsDSyqTJ+dL6KaDqxKLnmOIyQNdAVwRWUKi2eLK6qZmTVSkR36E0lzKlUsycuGawZpOpBqX1S658WZeb6otUiaKWmepHlLly4dQVgzMxupth4tJmkrYFfS5HkVp5D6YPYiTS8+4N33ImJWREyNiKkTJoyoP8rMzEaoyORyP7BN1etJedlwvBP4afXspxHxYL5HxQrSPTH2ftElNTOzhioyucwFdpS0Xb5BzwzWnCCvHkdTc0ost2bIM9ceQZpcz8zM2khhySUiekm3SL0WuAO4LCIWSDo9TwGOpL0kLSHdLOgcSZX7dyBpCqnlc2PNoS+SdDtpttzxDH+2XDMzK1hHTP8yderU8HUuVkbLnl/Jlbfez4y9JjO2p627SG0UknRzREwdyb7+azVrY+f++m5Ou2oB/32H78pr5eLkYtbG7n/yeQCWr+htcUnMhsfJxayNRb4NvVpcDrPhcnIxa2e5S7RLTi9WLk4uZm2sPw+4cW6xsnFyMSsBJxcrGycXszZWuVBA7nWxknFyMWtj/Tm7uOViZePkYtbGOuEiZxudnFzM2lgltXi0mJWNk4tZGwuPFrOScnIxKwF36FvZOLmYtbFwh76VlJOLWRtbdRFli8thNlxOLmZtzC0XKysnF7M2tuoiSmcXKxknF7M2tqrl0tpimA2bk4tZW/NFlFZOTi5mbSw85b6VlJOLWRvzlPtWVoUmF0nTJN0paZGkTw+wfl9Jt0jqlXRkzbo+SfPzY3bV8u0k/T4f81JJY4usg1krre7Qb2kxzIatsOQiqRs4GzgE2Bk4WtLONZvdB5wAXDzAIZ6LiN3z4/Cq5V8BzoyIHYAngPc1vPBmbWL1UGRnFyuXIlsuewOLIuLuiHgBuASYXr1BRNwbEbcB/fUcUOk/7ADg8rzoB8ARjSuyWXtZfT8Xs3IpMrlMBBZXvV6Sl9VrPUnzJN0kqZJAtgCejIjeER7TrFQ85b6VVU+rCzCEbSPifknbA9dLuh14qt6dJc0EZgJMnjy5oCKaFcujxaysimy53A9sU/V6Ul5Wl4i4P/+8G7gB2AN4DNhUUiUpDnrMiJgVEVMjYuqECROGX3qzNhB4tJiVU5HJZS6wYx7dNRaYAcxexz4ASNpM0rj8fDzwd8DCSOcIfgVURpYdD1zV8JKbtQm3XKysCksuuV/kJOBa4A7gsohYIOl0SYcDSNpL0hLgKOAcSQvy7q8C5kn6IymZfDkiFuZ1JwMfl7SI1AfzvaLqYNZqnv7FyqrQPpeImAPMqVl2WtXzuaRTW7X7/S+w6yDHvJs0Es1s1AtP/2Il5Sv0zdqYB4tZWTm5mLUxJxcrKycXszbm02JWVk4uZm3MLRcrKycXszYWNT/NysLJxayNefoXKysnF7M25tRiZeXkYtbG+p1drKScXMzamU+LWUk5uZi1sVUd+s4xVjJOLmZtzEnFysrJxayN+SJKKysnF7M25paLlZWTi1kb82gxKysnF7M25osoraycXMxKwH0vVjZOLmZtzA0XKysnF7M25haLlZWTi1kbc8vFysrJxayN9Tu7WEkVmlwkTZN0p6RFkj49wPp9Jd0iqVfSkVXLd5f0O0kLJN0m6V1V674v6R5J8/Nj9yLrYNZKnv7FyqqnqANL6gbOBg4ElgBzJc2OiIVVm90HnAB8omb3Z4HjIuIuSVsDN0u6NiKezOs/GRGXF1V2s7bhpGIlVVhyAfYGFkXE3QCSLgGmA6uSS0Tcm9f1V+8YEX+pev6ApEeACcCTmHUQ5xYrqyJPi00EFle9XpKXDYukvYGxwF+rFn8xny47U9K4F1dMs/bliyitrNq6Q1/SVsCFwIkRUWndnAK8EtgL2Bw4eZB9Z0qaJ2ne0qVLm1Jes0bz9C9WVkUml/uBbapeT8rL6iJpY+Aa4F8j4qbK8oh4MJIVwPmk029riYhZETE1IqZOmDBhRBUwaxfOMVY2RSaXucCOkraTNBaYAcyuZ8e8/U+BC2o77nNrBkkCjgD+1NBSm7URX0RpZVVYcomIXuAk4FrgDuCyiFgg6XRJhwNI2kvSEuAo4BxJC/Lu7wT2BU4YYMjxRZJuB24HxgNnFFUHs1br71/3NmbtqMjRYkTEHGBOzbLTqp7PJZ0uq93vh8APBznmAQ0uppmZNVhbd+ibdTqPFrOycnIxa2MeLWZl5eRiVgJuwVjZOLmYtTGPFrOycnIxa2M+LWZl5eRi1sZ8NszKysnFrK05u1g5ObmYtbFKy8UpxspmnclF0lGSNsrPT5X0E0mvLb5oZuakYmVVT8vl3yJiuaQ3Am8Bvgd8u9himRl4CLKVVz3JpS//PBSYFRHXkO6vYmYF82gxK6t6ksv9ks4B3gXMyTfncl+NWRO45WJlVU+SeCdpZuOD8z3sNwc+WWipzAxwn4uV1zqTS0Q8CzwCvDEv6gXuKrJQZpZFzU+zkqhntNhnSbcSPiUvGsMg0+GbWWM5p1hZ1XNa7G3A4cAzABHxALBRkYUys6TffS5WUvUklxci9SoGgKQNii2SmVU4t1hZ1ZNcLsujxTaV9AHgOuC7xRbLzMCzIlt5rfM2xxHxn5IOBJYBOwGnRcQvCy+ZmVVN/+IkY+WyzuQCkJOJE4pZkzmlWFnVM1rs7ZLukvSUpGWSlktaVs/BJU2TdKekRZI+PcD6fSXdIqlX0pE1647Pce+SdHzV8j0l3Z6P+U1JqqcsZqXk7GIlVU+fy38Ah0fEJhGxcURsFBEbr2snSd3A2cAhwM7A0ZJ2rtnsPuAE4OKafTcHPgu8Dtgb+KykzfLqbwMfAHbMj2l11MGslDxazMqqnuTycETcMYJj7w0sioi7I+IF4BJgevUGEXFvRNwG9NfsezDwy4h4PCKeIJ2SmyZpK2DjiLgpj2C7ADhiBGUzKwWnFiurQftcJL09P50n6VLgSmBFZX1E/GQdx54ILK56vYTUEqnHQPtOzI8lAyw3G5U8t5iV1VAd+odVPX8WOKjqdQDrSi4tJWkmMBNg8uTJLS6N2cismv3FOcZKZtDkEhEnvshj3w9sU/V6Ul5W77771+x7Q14+qZ5jRsQsYBbA1KlT/a9ppeSkYmVVz2ixH0jatOr1ZpLOq+PYc4EdJW0naSwwA5hdZ7muBQ7KsTYjtZqujYgHgWWSXp9HiR0HXFXnMc3MrEnq6dDfLU+1D0DuYN9jXTtFRC9wEilR3AFcFhELJJ0u6XAASXtJWgIcBZwjaUHe93HgC6QENRc4PS8D+AhwLrAI+Cvws7pqamZmTVPPRZRdkjbLSaUyTLjeiy/nAHNqlp1W9Xwua57mqt7uPGCtFlJEzAN2qSe+mZm1Rj1J4v8Cv5P0Y0DAkcAXCy2Vma3BfS9WNvXMLXaBpJuBN+VFb4+IhcUWy8zMyqyu01vAn4EnKttLmhwR9xVWKjMzK7V1JhdJ/0iaiuVhoI90aiyA3YotmpmZlVU9LZePATtFxGNFF8bMzEaHeoYiLwaeKrogZramxY8/2+oimI1YPS2Xu4EbJF3DmnOLfa2wUpkZ+/zHr1Y992AxK5t6kst9+TE2P8zMzIZUz1DkzzejIGZmNnoMNeX+1QzRGo+IwwspkZmZld5QLZf/bFopzMxsVBlqyv0bm1kQMxucbxpmZVPPUGQzM7NhcXIxM7OGc3IxM7OGq2dusVcAnwS2rd4+Ig4osFxmZlZi9VxE+WPgO8B3SRNXmpmZDame5NIbEd8uvCRmNiiPFbOyqafP5WpJH5G0laTNK4/CS2ZmZqVVT8vl+Pzzk1XLAti+8cUxM7PRoJ65xbZrRkHMzGz0GGpusQMi4npJbx9ofUT8ZF0HlzQN+AbQDZwbEV+uWT8OuADYE3gMeFdE3CvpGNZsKe0GvDYi5ku6AdgKeC6vOygiHllXWczMrHmGarnsB1wPHDbAugCGTC6SuoGzgQOBJcBcSbMjYmHVZu8DnoiIHSTNAL5CSjAXARfl4+wKXBkR86v2OyYi5g1dNbPRw7O/WNkMNbfYZ/PPE0d47L2BRRFxN4CkS4DpQHVymQ58Lj+/HDhLkmLNiZSOBi4ZYRnMzKwFirxCfyLpFskVS/KyAbeJiF7S7ZS3qNnmXcCPapadL2m+pH+TpIGCS5opaZ6keUuXLh1pHczMbATaevoXSa8Dno2IP1UtPiYidgX2yY9jB9o3ImZFxNSImDphwoQmlNbMzCqKTC73A9tUvZ6Ulw24jaQeYBNSx37FDGpaLRFxf/65HLiYdPrNzMzaSD3XuSDpDcAU1pxb7IJ17DYX2FHSdqQkMgN4d802s0nX0fwOOBK4vtLfIqkLeCepdVIpRw+waUQ8KmkM8FbgunrqYGZmzVPPxJUXAi8H5rN6brEgDSEeVET0SjoJuJY0FPm8iFgg6XRgXkTMBr4HXChpEfA4KQFV7AssrgwIyMYB1+bE0k1KLN9ddzXNys7Dxaxc6mm5TAV2rhnBVZeImAPMqVl2WtXz54GjBtn3BuD1NcueIV0TY2ZmbayePpc/AS8ruiBmZjZ61NNyGQ8slPQHYEVlYUQcXlipzMys1OpJLp8ruhBmZja61DNx5Y3NKIiZDc7Tv1jZDDVx5W8j4o2SlrPmUBUBEREbF146MzMrpaHmFntj/rlR84pjZmajQV0XUQJI2hJYr/I6Iu4rpERmZlZ66xyKLOlwSXcB9wA3AvcCPyu4XGZmVmL1XOfyBdLFjH/Jd6V8M3BToaUyM7NSqye5rIyIx4AuSV0R8SvSVftm1iQeLGZlU09yeVLShsCvgYskfQN4pthimbW/Py5+km/+912MYGYks1GvnuQyHXgW+Gfg58BfGfjWx2Yd5WOX3MrXfvkXli5fse6NzTpMPRdRVlop/cAP8lT4R5PvcW/Wqe597FkAXujrb3FJzNrPoC0XSRtLOkXSWZIOUnIScDfpPitmBvT1+7SYWa2hWi4XAk+QbuT1fuAzpKvzj4iI+U0om1lb6xL0B6zsKz65uFvHymao5LJ9vlc9ks4FHgQm53uwmHW8Lon+CHr7fVrMrNZQHforK08iog9Y4sRitlqXBEBvE1ouZmUzVMvlNZKW5ecC1s+vPXGlGZBzCyvdoW+2lqEmruxuZkHMymZVy8Ud+mZrqec6FzMbQJdbLmaDKjS5SJom6U5JiyR9eoD14yRdmtf/XtKUvHyKpOckzc+P71Tts6ek2/M+35QqJyfMmquZfS7hCWCsZApLLpK6gbOBQ4CdgaMl7Vyz2fuAJyJiB+BM4CtV6/4aEbvnx4eqln8b+ACwY35MK6oOZkOpfK3xdS5mayuy5bI3sCgi7o6IF4BLSFPJVJsO/CA/vxx481AtEUlbARtHxE2RJnS6ADii8UU3W7eufF7Mp8XM1lZkcpkILK56vSQvG3CbiOgFngK2yOu2k3SrpBsl7VO1/ZJ1HBMASTMlzZM0b+nSpS+uJmYDcIe+2eDatUO/csHmHsDHgYslDWvoc0TMioipETF1woQJhRTSOps79M0GV2RyuR/Ypur1pLxswG0k9QCbAI9FxIp8Dxki4mbSTMyvyNtPWscxzZpCzezQd+PISqbI5DIX2FHSdpLGAjOA2TXbzAaOz8+PBK6PiJA0IQ8IQNL2pI77uyPiQWCZpNfnvpnjgKsKrIPZoCotF0//Yra2dU65P1IR0ZtnUb4W6AbOi4gFkk4H5kXEbOB7wIWSFgGPkxIQwL7A6ZJWkqb6/1BEPJ7XfQT4PrA+8LP8MGs6UenQd7PCrFZhyQUgIuYAc2qWnVb1/HngqAH2uwK4YpBjzgN2aWxJzYavy0ORzQbVrh36Zm2v0ufiDn2ztTm5mI1QV/7vacZQZLeNrGycXMxGaPX0L265mNVycjEboS65Q99sME4uZiPUny8+8VBks7U5uZiNUGWUmO9EabY2JxezEerPycWnxczW5uRiNkJ9+bRYXxNOi4Xnf7GScXIxG6HKILGVvojSbC1OLmYjtKpD30ORzdbi5GI2Qu7QNxuck4vZCK3q0PdpMbO1OLmYjVCfT4uZDcrJxWyE+jwU2WxQTi5mI+Qr9M0G5+RiNkKVrhbfz8VsbYXeLMxsNLrvsWc545qFVafF3HIxq+WWi9kw/f6ex/jFwodXvfZQZLO1ObmYDVPtabBmDEX27C9WNk4uZsNUe+dJD0U2W1uhyUXSNEl3Slok6dMDrB8n6dK8/veSpuTlB0q6WdLt+ecBVfvckI85Pz+2LLIOZrX6oza5uFlhVquwDn1J3cDZwIHAEmCupNkRsbBqs/cBT0TEDpJmAF8B3gU8ChwWEQ9I2gW4FphYtd8xETGvqLKbDaU2maz0UGSztRTZctkbWBQRd0fEC8AlwPSabaYDP8jPLwfeLEkRcWtEPJCXLwDWlzSuwLKa1a22z8UtF7O1FZlcJgKLq14vYc3WxxrbREQv8BSwRc027wBuiYgVVcvOz6fE/k3KNzKvIWmmpHmS5i1duvTF1MNsDX01p8WacZ1L4ARm5dLWHfqSXk06VfbBqsXHRMSuwD75cexA+0bErIiYGhFTJ0yYUHxhrWOsNVrMHfpmaykyudwPbFP1elJeNuA2knqATYDH8utJwE+B4yLir5UdIuL+/HM5cDHp9JtZ01SfBuvu0lqjx8ys2OQyF9hR0naSxgIzgNk128wGjs/PjwSuj4iQtClwDfDpiPifysaSeiSNz8/HAG8F/lRgHczWUn1abFxPl1suZgMoLLnkPpSTSCO97gAui4gFkk6XdHje7HvAFpIWAR8HKsOVTwJ2AE6rGXI8DrhW0m3AfFLL57tF1cFsIH1Vo8PG9XQ1pUP/sadfIHwlpZVIoXOLRcQcYE7NstOqnj8PHDXAfmcAZwxy2D0bWUaz4ao+DTa2p4unV/QWHvOMa+6gtz/40H4vLzyWWSO0dYe+WTvqr0kuK/uiKa2Ke5Y+U3gMs0ZxcjEbpuqWy7iebqA5w5Gb0UIyaxQnF7Nh6utfs0Mf1p5vrAjLnl9ZeAyzRnFyMRum6uQyprt5ycUtFysTJxezYVozuaQJIpoxM/Ly551crDycXMyGqbqV0t2VksvKJgxHXu7TYlYiTi5mw1Q9Wqynq3JarPiWy9NuuViJOLmYDVN1y6Wrq3JarPiWyzMv9DVlVJpZIzi5mA1T3xotl8ppseZMAeNOfSsLJxezYeoboM+lWS0K97tYWTi5mA1T74Atl+YkF7dcrCycXMyGqdJ5P66ni9dss+kay4owtqeLQ3fdCvBwZCsPJxezYXqht5+9pmzGwtOn8cqXbQQU03LZdouXcMTuW3PH6dP4wL7bAz4tZuXh5GI2TCt6+xnX0013l1ZfoV9Ah35/BJLo7hIbjksTmLvlYmXh5GI2TCt6+xib5xSr9LkUMf1LBCg/33g9JxcrFycXs2FasbJ/1YSVPd3FDUWOACkdf8OcXNyhb2Xh5GI2TOm0WKXlUjktVkTLJci5hfXHpNNw7nOxsnByMRumF3KfC6xuuRRyWgzIZ92QxEbr9fi0mJWGk4vZMK3o7WPcmPSvs3rK/YI69Ff1usCG43o8v5iVhpOL2TCt6O1nbHdNh34hp8VYdVoMYKP1xrDMycVKotDkImmapDslLZL06QHWj5N0aV7/e0lTqtadkpffKengeo9pVrQVvf1rtVyK6NDvr+rQB9hoXA9Pr3Cfi5VDYclFUjdwNnAIsDNwtKSdazZ7H/BEROwAnAl8Je+7MzADeDUwDfiWpO46j2lWmN6+fvr6oyl9LhA1LRf3uVh59BR47L2BRRFxN4CkS4DpwMKqbaYDn8vPLwfOUvqqNh24JCJWAPdIWpSPRx3HbJhTr7ydP9zzeBGHftH6+oMnnl3Jej1dbDCuZ40PIStO5Ur89WpaLl+99k7O++09Dfs9RMCjT79Ad9UBN1yvh7vuepqDzryxMUGsI5x73F5M3uIlTY9bZHKZCCyuer0EeN1g20REr9v5NhgAABB2SURBVKSngC3y8ptq9p2Yn6/rmABImgnMBJg8efKIKrDVJuvz8gkbjmjfokmw6UvGsmJlP8++4G+zzdLdJfbfaQJv3W1rALbYYCynHPJK7nxoOc+t7GtorKlTNuPdr1v9t/uuqduwsq+f8C1dbBgqF/w2W5HJpaUiYhYwC2Dq1Kkj+nf8hzft0NAy2egjiQ/u9/KmxHrDDuN5ww7jmxLL7MUqMqXdD2xT9XpSXjbgNpJ6gE2Ax4bYt55jmplZixWZXOYCO0raTtJYUgf97JptZgPH5+dHAtdHROTlM/Josu2AHYE/1HlMMzNrscJOi+U+lJOAa4Fu4LyIWCDpdGBeRMwGvgdcmDvsHyclC/J2l5E66nuBf4iIPoCBjllUHczMbGQUHdA7OHXq1Jg3b16ri2FmViqSbo6IqSPZ11fom5lZwzm5mJlZwzm5mJlZwzm5mJlZw3VEh76kpcDfgPHAoy0sSivjd3LdHd9/e44/MttGxISR7NgRyaVC0ryRjnwoe/xOrrvj+2/P8Zsf36fFzMys4ZxczMys4Totuczq4PidXHfH99+e4zdZR/W5mJlZc3Ray8XMzJrAycXMzBrOyWWY8m2YOyp+O9W51WVppVbXvdXxrXVG8rt3cqmDpIMkfQkgWtBJJemlksZX4rfgn3yDmvI0++9mM0ndsKr+TYvfgrrWxt8030ivJb97SeMlbdjC+K+QtF4zY1bF3kPSXq2IneMfIOmDLYx/mKTzYWSfe04u6yDpIOA7wOsk7diC+IcAPwfOkvQdaO4/uaSDgcslnSrptBy/v1kfupIOB64j1f+7lfhNin0A8G5JmzUj3gDxDybdDO/bks6E5n65kTQN+C/gm5JmtSD+ZODPwD80+3eQ634+8HzN8mb93x0O/D9q7rTbxPgHAv8B7CbpLSM6SET4McgDOBi4GXgbcB7wsSbHfx3phmlvAXYAzgG6qtZ3FRx/b+BO4O+BVwG/AC5rYvwdgT8CbwK2Jt0kbjawYdHxgb8D+oFfAu8CNmvy7/4tpA/Ww4A9gR8D725y/IXAIcArgYuBlzTrd59jvDS/B9cB/wxs2qS6H0D6UN8rvx5Xs77ov/txwIXAfvn1hsDmTfzdH5T/76YBJwOnjeQ4hd2Jsszyt4MJwEeBj0fEjZKeAL4l6TcRcUuTivIS4PKIuE7STqR/9H+XtGVEnBipBaHIfxEFxb8sIuYA5CbyNyVdEhEzovgWxJPAXcAdEfEQcLCki4EfAYcVVf/cKtsMeGdeNB3okvTziHgib1PY+y5pfWAf4OSIuDr/Pd5CSrCFy/H3BD6c//Z3Ad4A/LOkzSLiE0342yMiHpZ0Lulv4CRgmaQ7gCci4o4iYkoaA+wB3A48JGlT4OuSngI2iIj3N6HufcAmQLeklwKXAk/m04Mfj4iFRcXP9X0ncFJE/EbSY8B/5c+9Xw3nWD4tNoBIHgGOyf9cY4D/JX1r3gWg0gdQsJXAkZI+R/r29l3gXGCypCsqZS0wvoD3Sdo1v94B+AQwRtK7CgmYm/25n6EPeAJYNS9SRLwbGCfprPy6YfWvxM5J89fANRFxOem05CHA30vavNFxB4j/HOmb69yqD5E7SC3JwtTE/3b+298Y+Awpof8U2L2ov72q331X1WnXbUl/h9OBDwO/AbZpZNzq2BGxktRKvJJ0WmgB6b2/BNhe0qV5u6Lq3h0RvaT3ehfgVOCHEXEEqSX5tSLjR8STpAT2G0k9ETGX9D68UVLPcE6HO7nUkPR6SR+U9Bpg/by4NyJeAP4KfErSBhHRV2D8mZJeExG/BY4Bfgf8KiK+EBF3AW8HnsvfMIuM/yvSB8uPJV1EOl1wAXAj6ZtVEV4GEBG9EfE4Kal+XtIbq7b5CNBbVOwcfxn5fHtEXEQ6PXYwsKekT0n694LjL4qIB6o+RFaQZrdF0rGSPlZkfGB5/vk88MWIOCUiFgLvAZ7OX7gKiZ+Te6Vv4ULSgJItSafJ5gM75KTX8Ng5/n3AnBzrixHxlYj4Hanuzxb0xbJS98rnyiLgzcArSEmFiPg4EJK2Lyp+9kyOV/kfu5P0t79ZpdVW1xGLPn9XpgdwKHAv8G3gB6QOvak125wPfJ48u0GB8X+YY702r5sD7JyfHw/8ltz3UHD8iaREsg3Qk7f7DHA66QOgYe8D6dtpP/CJmuXvB+aRkupk4L2kluRLioxN+vKlqtcHks5FLwb2aPB7P2R8YCfgG/l3NLfyt1Bw/J4Btns/qe+tYe/9EPFFOj23GHiI1X1vVwJbNOHvbpPq9wD4QLPqnpcfRPpy9UlS6/1t+e+vYXUf6r2v2WYWqfXWXfdxG1nIsj+Az5FOhVX+mU8CrgH2rNrm/cCZwHpNiP+POf5k4ETgcdIH/+3Aq5sQ/2M5/tSqbd6b/9lf2eDYk0gJ9GRgCfDJmvXvII2e+QnwB2C3ZsSu+YB/G+kbfaM/2IeMn39OyR8A8xv9u68z/nr5b3B+M+uf138EOLTqdcP+99ZR9+ovFscDt7Wg7vsB/0Q6Y/CzRv7d1/O7r/r97wN8E9io7mM3sqBlfwBfAr5b9Xo8KcGcA0zIyzYGtmpi/I8C5+TX+wL7A9u1oP5bkAY5fAZ4VQGxu4CD8/OdgaUD/KONI3W0b9ns2HndwY3+cBlG3bcmtVZf0aL4k0l9fi2pf17XTeNby/XU/WX5g7WVdV+fAkbLDSP+BsD4YR270YUt8wPYlNSRW9083IXUHGz4P/Uw4v+4xfFX1R8YW0Bc1T4nDX1e9YdOGpY9qYWxX1bQe15v/PVIo5VaFX+jFvzuP1UVv+Ff6IZR942BMS2M3/C/+2bEb3iBy/pgdfPv9aSRGidXrbuIgq8xKEH8yumyhvc1DVCWMfnnq0hN9Z+R+hkKaTG2S+wh4t/c4fFb+buf1wl1LyJ+4QUuw4Oai6KAvYDrSRdOnkIaJVbIqSjHXyOxrZW4SAMHHgV2HW2xHd+/+9EcvyMvosxTGxxI6iD/YUQsyePL+/IFY5sBh5M6EgM4IiLucfxC43dFGua4GymRXZWHg/8f4ICIuL3ssR3fv/uOil9UVmzXB2ko563AvwDfAmZWrXszaSTW/o7f0vh/l193McxOxHaN7fj+3Xda/I66E2W++OlbwBUR8QulGUe3J51bvIvUedUTEZcVMb2C4w8rfnc08ELVVsZ2fP/uOzF+p50WE2nkx4GSHgE+TuqsPIrUz3BkRNxXxAer4w87fqNnQGhlbMf3777j4ndEy0XSy0jT8TwsaVvg66TpQx6IiI/lbb4K3BsRZzv+6InfyXXv9PidXPd2iD/qWy6S3kG6wnWMpKuBn0fE2yQdSZqIsSIoYL4sx29d/E6ue6fH7+S6t0N8GOUtF0lbkObmeS9phuEDSWO4rwZuyutmA38DPgS8JyL+7Pjlj9/Jde/0+J1c93aIXzHaWy7dwDLgnoh4UuneBG8h3YDpIeDdwGmkaU5OLOANdvzWxe/kund6/E6uezvEB0Z5ywVA0jdI8+J8LCKekbQ1cBzQFxFfVbpvSH8UdOMrx29d/E6ue6fH7+S6t0N8GMX3c9Hqm9qcTcriJyvdh+UB0u1y3ypp80j3DSnij9vxWxS/k+ve6fE7ue7tEL/aqEsu0hp3E4Q0dclPSLOKfkfSeNINeHop4IZTjt+6+J1c906P38l1b4f4A5ZptJwWU7r97PMR8WzVsrER8YKkScDmpHsy7JyffzgibnH88sfv5Lp3evxOrns7xB9SNHCKgVY9SHdSu5qUqY9lzZtbvZk0Zf3k/HoTGjx1ueO3Ln4n173T43dy3dsh/jrL18xghVQgNfVuJ2XmfYGvAj8i3TltDGno3Tscf/TF7+S6d3r8Tq57O8Sv5zEahiKPB5ZExEIASZuRZvN9J+mmN9MjXaFa1JQmjt+6+J1c906P38l1b4f46zQaOvT/BCyTdGp+vQfwF+B5YEpEPAxpDgTHH3XxO7nunR6/k+veDvHXqZQd+rmjanlEPKU0XvtNwD+QkmVExHRJRwPTgBMa/QY7fuvid3LdOz1+J9e9HeIPW7TwnNxIHsARwJ9JM3tOqFreBbwM6M6vZwLfcPzRE7+T697p8Tu57u0QfySPUrVcJE0ALgHuI93j+RHgkohYWrPdPwEnkubMaeSd3By/RfE7ue6dHr+T694O8UeqbH0uT5Fm+vwQMJ80u+cMSVvC6guJSHPrHFPAG+z4rYvfyXXv9PidXPd2iD8ipWi5SJpMmnCtJ9a8WOgdwH7AXRHx/yS9Ngq4QMjxWxe/k+ve6fE7ue7tEP/FavuWi6RDgTnAWcD5kl5ZWRcRVwA3AhMkXQncKGmi44+O+J1c906P38l1b4f4DdHqTp/BHoCAbUgXCu0PvBT4BPAg8OqabX8I3Avs6vjlj9/Jde/0+J1c93aI38hHywuwjje6G5gFTGT1KbyPAfcDr8ivtwIWArs7/uiJ38l17/T4nVz3dojfsHq0ugCDvLk7AHsBWwCXAp+qWf8p4PvA+vn1ho4/OuJ3ct07PX4n170d4jf60fICDPAGvxW4jXRO8SzgcFLT75SqbaYA55CzuuOPjvidXPdOj9/JdW+H+EU8Wl6Amjf4DcAdwB759SzgDGBr0hjvU0nZ/QRgHrCZ44+O+J1c906P38l1b4f4RT1aXoAB3uQTql5PAK7Jz7cHzgO+BdxMAZ1Yjt+6+J1c906P38l1b4f4RT1aXoCaN7kb2Ljq+STgVmCrvGxboAfYxPFHV/xOrnunx+/kurdD/KIebXWdS0T0RcSy/FLAk8DjEfGgpPcAnwHGRMRTjj+64ndy3Ts9fifXvR3iF6Xtr9CX9H3SGO+DSE3Hpk5t4Piti9/Jde/0+J1c93aI3whtm1zyfDljSB1dY4A3R8Rdjj/643dy3Ts9fifXvR3iN1LbJpcKSScAcyNigeN3VvxOrnunx+/kurdD/EYoQ3Jp2W06Hb+18Tu57p0ev5Pr3g7xG6Htk4uZmZVPW40WMzOz0cHJxczMGs7JxczMGs7JxczMGs7JxaxBJPVJmi9pgaQ/SvoXSUP+j0maIundzSqjWbM4uZg1znMRsXtEvBo4EDgE+Ow69pkCOLnYqOPkYlaAiHgEmAmcpGSKpN9IuiU/3pA3/TKwT27x/LOkbklflTRX0m2SPgggaStJv87b/UnSPq2qm1k9fJ2LWYNIejoiNqxZ9iSwE7Ac6I+I5yXtCPwoIqZK2h/4RES8NW8/E9gyIs6QNA74H+Ao4O3AehHxRUndwEsiYnnzamc2PD2tLoBZhxgDnCVpd6APeMUg2x0E7CbpyPx6E2BHYC5wnqQxwJURMb/oApu9GE4uZgWRtD0pkTxC6nt5GHgN6XT084PtBvxjRFw7wPH2BQ4Fvi/paxFxQSEFN2sA97mYFUDSBOA7wFl5jqhNgAcjoh84lnRTKEinyzaq2vVa4MO5hYKkV0jaQNK2wMMR8V3gXOC1TaqK2Yi45WLWOOtLmk86BdYLXAh8La/7FnCFpOOAnwPP5OW3AX2S/gh8H/gGaQTZLXn69aXAEcD+wCclrQSeBo5rQn3MRswd+mZm1nA+LWZmZg3n5GJmZg3n5GJmZg3n5GJmZg3n5GJmZg3n5GJmZg3n5GJmZg33/wGLB4EOIhRZlAAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Filter down to last ten days\n",
"days_to_plot = 10\n",
"start_date = datetime.today() - timedelta(days=days_to_plot)\n",
"recent_rain = rain_vals[rain_vals.valid > start_date]\n",
"\n",
"# Plot rainfall amount to make sure it actually rained in that period\n",
"plt.plot(recent_rain.valid, recent_rain.p01i)\n",
"plt.xticks(rotation=45)\n",
"plt.title(\"Rain For Last Ten Days At KDSM\")\n",
"plt.xlabel(\"Dates\")\n",
"plt.ylabel(\"Rain in Inches\")"
]
},
{
"cell_type": "code",
"execution_count": 77,
"metadata": {},
"outputs": [],
"source": [
"# Cool, that worked. Now Let's make a rule\n",
"def moisture_content(previous_value, precip_in, decay_coefficient):\n",
" \"\"\"\n",
" Simple function to do a 'Newton's Law of Cooling' eque calculation for soil moisture\n",
" Args:\n",
" previous_value (float): This should be the previous moisture reading\n",
" precip_in (float): This is the precip for the current time period\n",
" decay_coefficient (float): how much to decay the moisture per reading period. Between 0-1\n",
" \"\"\"\n",
" return previous_value*decay_coefficient + precip_in"
]
},
{
"cell_type": "code",
"execution_count": 78,
"metadata": {},
"outputs": [],
"source": [
"# Strip these out into their own list because I don't know how to do this well in Pandas\n",
"# Also it's Saturday. I don't have time/motivation to figure it out. I'd rather be biking\n",
"dates = recent_rain.valid.values\n",
"precip = recent_rain.p01i.values"
]
},
{
"cell_type": "code",
"execution_count": 79,
"metadata": {},
"outputs": [],
"source": [
"# Actually calculate the moisture values\n",
"moisture = []\n",
"decay_rate = 0.9 # Should be between 0 and 1 unless the group is magically making water\n",
"for index, value in enumerate(precip):\n",
" if index > 0:\n",
" moisture.append(moisture_content(moisture[index-1], value, decay_rate))\n",
" elif index == 0:\n",
" moisture.append(value)"
]
},
{
"cell_type": "code",
"execution_count": 80,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.legend.Legend at 0x7f727a07cef0>"
]
},
"execution_count": 80,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAE8CAYAAADe7fZ4AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd3wc1bn/8c9X3UVylXFvwQ1cwTY2pgVCKKG30EMIkPIDnBB6uEBMkpsAyU0ITggEAiRcU2+AgMFU0wLEBQMugA22sQwY94qxpH1+f8ysvBIqK3lnVys979drbe3s7Dxndnb32XPOzDkyM5xzzrm65GS6AM4555o3TxTOOefq5YnCOedcvTxROOecq5cnCuecc/XyROGcc65enihcrSQ9Jek79Tx+m6T/SnJbMyWdl7rS1Rpji6SBUcZwu0bSYEkbMl0O13ieKFoRScskfSOZdc3sCDO7J3zeOZJerfH4D8zshhSU6XpJJmlyjeWTw+XXJ1ne9mb2UQOxDpJUtgvFbTRJd0vaIWlzeJsv6b8ldUhnOeoiqW+YZOM3k7Q14f7+u7DtzyTtF79vZh+YWcfUlNylkycK1xx8AJxdY9l3wuXNhqS8Jj71RjMrBkqB7wITgNcktUtZ4ZrIzD4Ok2x7M2sfLh6VsOyVjBbQNQueKFqpeC1B0s2S1ktaKumIhMdnSjpP0jDgNmBi+AtzQ/j43ZJ+Ef7dSdITklaH23pCUu9GFGcW0FbSnuH29gSKwuWJZT5f0hJJ6yQ9LqlnwmMmaffw7yMlLQx/wa+UdGn4pfwU0DPh13LPxP0In1ut1hHWwq6Q9A6wVVJe+LxHwv1dKuniZHbSzLab2SzgGKALQdKIxzlX0qLw9ZshqV/CY3tKejbc71WSrg6Xj5f0uqQNkj6VdKukgvCxqZJ+W+P1e1zST5Ipa43ntZH0e0krwlrCHyUVho91l/R0WIa1kl4Ilz8EdAOeCV/riyUNlVSRsN03JF0X/r9J0nRJnRIePy+MuVrS5TVrKC59PFG0bvsA7wNdgRuBOyUpcQUzWwT8AHg9/IVZW9NBDvA3oB/QF/gCuLWRZfk7O2sV3wnvV5F0MPDfwClAD2A5cH8d27oT+H74K3448IKZbQWOAD5J+LX8SZJlOw34FtARiAH/At4GegGHAD+WdFiS28LMNgPPAvuH+3YscDVwAkGt4xVgWvhYMfAc8DTQE9gdeD7cVCXwE4LjNzEsy4/Cx+4BTpOUE26nK/AN4H+TLWeC3wG9gRHAEGAwcGX42BXsfA/1AK4P9/Fk4HPgm+FrfUsd2z4dOCN8bkdgcljeMWHck8PYvcMYLgM8UbRuy83sDjOrJPhi6QHs1tiNmNlaM3vEzLaFX4K/BA5s5Gb+QfDFlg+cGt5PdAZwl5nNNbMvgasIajn9a9lWObCHpBIzW29mcxtZlppuMbMVZvYFMA4oNbMpZrYj7Be5IyxzY3wCdA7//gHw32a2yMwqgF8Bo8NaxVHAZ2b227BGstnM3gQwszlm9oaZVZjZMuAvhK+7mf0H2EiQPAjLN9PMVjWmkGFz2/eAyWa2wcw2Ar9O2N9yggTWN3w9Xm7k63CHmX0YJvKHgdHh8pOBR8L9+xK4Bv++yhh/4Vu3z+J/mNm28M/2daxbJ0ltJf1F0nJJm4CXgY6ScpPdhpl9DCwh+JJcbGYraqzSk6AWEV9/C7CW4Fd9TScCRwLLJb0kaWLj9ugrEsvSj6D5akP8RlAbaGyC7QWsS9jmHxK2tw5QuE4f4MPaNqDgLKInwiaZTQSvXeKv7nuAM8O/z6RGLS1JPYF8YEFC+R4laFaC4EfBJ8CLYbPgJY3c/mcJf29j5/uvJwmvu5ltIkh8LgM8UbhkNDTE8E8JmiT2MbMS4IBwuep+Sq3uDbd1by2PfULwhRpsOOhz6AKs/EphzWaZ2bEEX2aPAg/GH6plu1uBtgn3u9eyTuLzVgBLzaxjwq3YzI6se7eqk9SeoBko3lG8gqCpLHGbbczs3+FjdZ32+2fgPWBQ+LpfTfXX/B/AsZJGAcMIXovG+hSoAL6WULYOZtYFwMw2mtlkM+tHkKCvkTQpfO6uDE39KUFzEwCSSoBmcaZYa+SJwiVjFdA73lFai2KCfokNkjoD1zUxzgPAN9n5xZ5oGvBdSaPDjtRfAW+GTS5VJBVIOkNSBzMrBzYR9CvE96OLqp+aOg84UlJnSd2BHzdQxv8Am8MO7jaSciUNlzSuoZ2TVChpb4Iv7PUE/ToQnCxwlXZ25neQdHL42BNAD0k/Dp9fLGmf8LHicP+2SBoK/DAxnpmVEZwQ8HeCZpwvGipjTeFreBdBjaerAn0kHRqW9RhJA8O+rY0E/SaJr3dTr215EDhR0rjwfTclYbsuzTxRuGS8ACwAPpO0ppbHfw+0AdYAbxB0vDaamX1hZs/V9oVmZs8B/wU8QvBr82vU3S9wFrAsbI75AUH/Bmb2HkHC+ShsRulJ8CX6NrAMeIYgWdVXxkqCfoPRwFKCff4r9f/avVzSZoKmsnuBOcC+Ybs8ZvZP4DfA/WGZ5xN0vMc7vg8FjiZoplkMfD3c7qUEncGbCfpJaiv7PQSd0E1pdor7MUGNbjZBMniaoFMdgprKi2EZXgZuNrPXw8d+CfwyfK0vbExAM3sLuAz4J0Gt8dMw9pe7sB+uieQTFznXckk6gKAJqp9l8Yc9PG12HdDTzD7NdHlaG69RONdChWeQTQb+mo1JImzWahP26fyOoKnRk0QGeKJwrgVScKHkBoJTnn+f4eI01ckEzW1lBGeAnZHZ4rRe3vTknHOuXl6jcM45V6+mDnKWMV27drX+/ftnuhjOOZdV5syZs8bMSpvy3KxLFP3792f27NmZLoZzzmUVScsbXqt23vTknHOuXp4onHPO1csThXPOuXplXR+Fc655KC8vp6ysjO3bt2e6KC5BUVERvXv3Jj8/P2Xb9EThnGuSsrIyiouL6d+/PzXmu3IZYmasXbuWsrIyBgwYkLLtetOTc65Jtm/fTpcuXTxJNCOS6NKlS8preZ4onHNN5kmi+YnimHiicC6Nnl24is83e5u+yy6eKJxLk/LKGOffO5vT73gz00VpMSRx5plnVt2vqKigtLSUo446qt7nPf744/z617+u8/F58+Yxffr0lJUz23micC5NKmPBAJwfr9vWwJouWe3atWP+/Pl88UUw19Wzzz5Lr161TaNe3THHHMOVV15Z5+NNSRRmRizWMifh80ThXJrEB2r2Vv3UOvLII3nyyScBmDZtGqeddlrVY+vWreO4445j5MiRTJgwgXfeeQeAu+++mwsvDCbde+ihhxg+fDijRo3igAMOYMeOHVx77bU88MADjB49mgceeIDrr7+em2++uWq7w4cPZ9myZSxbtowhQ4Zw9tlnM3z4cFasWMEzzzzDxIkT2WuvvTj55JPZsmVLGl+NaPjpsc6lWUvs//35vxaw8JNNKd3mHj1LuO7oPRtc79RTT2XKlCkcddRRvPPOO5x77rm88sorAFx33XWMGTOGRx99lBdeeIGzzz6befPmVXv+lClTmDFjBr169WLDhg0UFBQwZcoUZs+eza233grA9ddfX2f8xYsXc8899zBhwgTWrFnDL37xC5577jnatWvHb37zG373u99x7bXXNv2FaAY8UTiXJobP/RKFkSNHsmzZMqZNm8aRRx5Z7bFXX32VRx55BICDDz6YtWvXsmlT9YQ2adIkzjnnHE455RROOOGERsfv168fEyZMAOCNN95g4cKFTJo0CYAdO3YwceLEpuxWs+KJwrk02dn01PKqFMn88o/SMcccw6WXXsrMmTNZu3Zto55722238eabb/Lkk0+y9957M2fOnK+sk5eXV63/IfE6hXbt2lX9bWYceuihTJs2rQl70Xx5H4VzadYSm54y7dxzz+W6665jxIgR1Zbvv//+3HfffQDMnDmTrl27UlJSUm2dDz/8kH322YcpU6ZQWlrKihUrKC4uZvPmzVXr9O/fn7lz5wIwd+5cli5dWms5JkyYwGuvvcaSJUsA2Lp1Kx988EHK9jNTPFE4lybe8BSd3r17c/HFF39l+fXXX8+cOXMYOXIkV155Jffcc89X1rnssssYMWIEw4cPZ99992XUqFF8/etfZ+HChVWd2SeeeCLr1q1jzz335NZbb2Xw4MG1lqO0tJS7776b0047jZEjRzJx4kTee++9lO9vumXdnNljx441n7jIZaPN28sZcf0ztCvIZcGUwzNdnF22aNEihg0bluliuFrUdmwkzTGzsU3ZntconEuT+E8yH/bCZRtPFM6liV9H4bKVJwrn0qWqSpHRUjjXaJ4onEsTv47CZatIE4WkwyW9L2mJpK8MrCKpr6QXJb0l6R1JR9a2HedaAm96ctkqskQhKReYChwB7AGcJmmPGqtdAzxoZmOAU4E/RVUe5zLNO7NdtoqyRjEeWGJmH5nZDuB+4Nga6xgQv/qlA/BJhOVxLqPip6J7nkid3NxcRo8ezfDhwzn66KPZsGFDg8/Zd999GxXjnHPOYcCAAYwePZpRo0bx/PPPN7W4X9G+fXsAPvnkE0466aSUbTfVokwUvYAVCffLwmWJrgfOlFQGTAcuqm1Dki6QNFvS7NWrV0dRVuciF/Omp5Rr06YN8+bNY/78+XTu3JmpU6c2+Jx///vfjY5z0003MW/ePH7/+9/zgx/8oClFrVfPnj15+OGHU77dVMl0Z/ZpwN1m1hs4Evi7pK+UycxuN7OxZja2tLQ07YV0LhW8MztaEydOZOXKlQBs2bKFQw45hL322osRI0bw2GOPVa0X/xU/c+ZMDjroIE466SSGDh3KGWecQUMXICfGAJgzZw4HHngge++9N4cddhiffvopAHfccQfjxo1j1KhRnHjiiWzbFsxBsnTpUiZOnMiIESO45pprqrazbNkyhg8fDgRDoJ9wwgkcfvjhDBo0iMsvv7xqvTvvvJPBgwczfvx4zj///Kqh0qMW5aCAK4E+Cfd7h8sSfQ84HMDMXpdUBHQFPo+wXM5lRkvOE09dCZ+9m9ptdh8BR9Q9C12iyspKnn/+eb73ve8BUFRUxD//+U9KSkpYs2YNEyZM4JhjjvlK/9Bbb73FggUL6NmzJ5MmTeK1115jv/32qzPO008/zXHHHQdAeXk5F110EY899hilpaU88MAD/OxnP+Ouu+7ihBNO4Pzzzwfgmmuu4c477+Siiy5i8uTJ/PCHP+Tss8+ut/Yzb9483nrrLQoLCxkyZAgXXXQRubm53HDDDcydO5fi4mIOPvhgRo0aldTrs6uiTBSzgEGSBhAkiFOB02us8zFwCHC3pGFAEeBtS65Fshr/u133xRdfMHr0aFauXMmwYcM49NBDgaA/6Oqrr+bll18mJyeHlStXsmrVKrp3717t+ePHj6d3794AjB49mmXLltWaKC677DKuvvpqysrKeP311wF4//33mT9/flXMyspKevToAcD8+fO55ppr2LBhA1u2bOGwww4D4LXXXqsa9vyss87iiiuuqHW/DjnkEDp06ADAHnvswfLly1mzZg0HHnggnTt3BuDkk09O24CDkSUKM6uQdCEwA8gF7jKzBZKmALPN7HHgp8Adkn5C8Pk5x7Jt8CnnkhR/Z8diLfAtnuQv/1SL91Fs27aNww47jKlTp3LxxRdz3333sXr1aubMmUN+fj79+/evNjR4XGFhYdXfubm5VFRU1Brnpptu4qSTTuKPf/wj5557LnPmzMHM2HPPPasSR6JzzjmHRx99lFGjRnH33Xczc+bMqseSOest2XKlS6R9FGY23cwGm9nXzOyX4bJrwySBmS00s0lmNsrMRpvZM1GWx7lMivdRtMA0kXFt27bllltu4be//S0VFRVs3LiRbt26kZ+fz4svvsjy5ctTEufCCy8kFosxY8YMhgwZwurVq6sSRXl5OQsWLABg8+bN9OjRg/Ly8qphziGYJOn+++8HqLY8GePGjeOll15i/fr1VFRUVNVM0iHTndnOtRpeV47WmDFjGDlyJNOmTeOMM85g9uzZjBgxgnvvvZehQ4emJIYkrrnmGm688UYKCgp4+OGHueKKKxg1ahSjR4+uOqPqhhtuYJ999mHSpEnVYv/hD39g6tSpjBgxolqneDJ69erF1Vdfzfjx45k0aRL9+/evap6Kmg8z7lyalK3fxn6/eZHiwjze/flhmS7OLvNhxtNvy5YttG/fnoqKCo4//njOPfdcjj/++K+s58OMO5elqvoosuzHmWs+rr/++qoLDAcMGFB1BlbUfM5s59LM04Rrqptvvjkjcb1G4VyaxCsSLalCkW1N161BFMfEE4VzaRI/66mlND0VFRWxdu1aTxbNiJmxdu1aioqKUrpdb3pyLk2qahSZLUbK9O7dm7KyMnz8tealqKio6iLCVPFE4VyatJQEEZefn8+AAQMyXQyXBt705FyaxFpalcK1Gp4onEsTPz3WZStPFM6ljQ/h4bKTJwrn0mTn6bGeKlx28UThXJr4MOMuW3micC5NWuIFd6518EThXJp4J7bLVp4onEsTzxMuW3micC5NzHsnXJbyROFcmniNwmUrTxTOOefq5YnCuTTxGoXLVp4onEsT76Nw2coThXNpEvM84bKUJwrn0sSH7nDZyhOFc2niacJlK08UzqWJ1yhctvJE4VyaeJ5w2coThXNp4p3ZLlt5onAuTXxQQJetPFE4lyaeJ1y28kThXJp4Z7bLVp4onEsT76Nw2coThXNp4kN4uGzlicK5NPEahctWniicSxPvo3DZyhOFc2niecJlK08UzqWJX0fhslWDiULSpGSW1fHcwyW9L2mJpCvrWOcUSQslLZD0v8ls17ls5HnCZatkahR/THJZNZJyganAEcAewGmS9qixziDgKmCSme0J/DiJ8jiXlbxG4bJVXl0PSJoI7AuUSrok4aESIDeJbY8HlpjZR+H27geOBRYmrHM+MNXM1gOY2eeNK75z2cPPenLZqr4aRQHQniCZFCfcNgEnJbHtXsCKhPtl4bJEg4HBkl6T9Iakw2vbkKQLJM2WNHv16tVJhHauOfJM4bJTnTUKM3sJeEnS3Wa2PML4g4CDgN7Ay5JGmNmGGmW5HbgdYOzYsf5pc1nJaxQuW9WZKBIUSrod6J+4vpkd3MDzVgJ9Eu73DpclKgPeNLNyYKmkDwgSx6wkyuVcVvE+CpetkkkUDwG3AX8FKhux7VnAIEkDCBLEqcDpNdZ5FDgN+JukrgRNUR81IoZzWcPzhMtWySSKCjP7c2M3bGYVki4EZhB0ft9lZgskTQFmm9nj4WPflLSQIAldZmZrGxvLuWzgNQqXrZJJFP+S9CPgn8CX8YVmtq6hJ5rZdGB6jWXXJvxtwCXhzbkWzfOEy1bJJIrvhP9flrDMgIGpL45zLZePHuuyVYOJwswGpKMgzrV0sVimS+Bc0zSYKCSdXdtyM7s39cVxruXyPgqXrZJpehqX8HcRcAgwF/BE4VwjJKYJM0NSxsriXGMk0/R0UeJ9SR2B+yMrkXMtVOJ8FJUxIy/XE4XLDk0ZZnwr4P0WzjVSYstThV+m7bJIMn0U/2JnrTkXGAY8GGWhnGuJEnNDpScKl0WS6aO4OeHvCmC5mZVFVB7nWqzEzmyvUbhs0mDTUzg44HsEI8d2AnZEXSjnWqLE1OA1CpdNkpnh7hTgP8DJwCnAm5KSGWbcOZfAqtUo/KIKlz2SaXr6GTAuPqmQpFLgOeDhKAvmXEsTi1U/68m5bJHMWU85NWaeW5vk85xzCRJTQ0WlJwqXPZKpUTwtaQYwLbz/beCp6IrkXMvkZz25bJXMBXeXSToRmBQuut3M/hltsZxreczPenJZKpkaBWb2iKRn4+tL6pzMMOPOuZ3MaxQuSyVzwd33gZ8D24EYIHyYcecaLeZnPbkslUyN4lJguJmtibowzrVk3kfhslUyZy99CGyLuiDOtXSJExd5H4XLJsnUKK4C/i3pTapPhXpxZKVyrgXyPgqXrZJJFH8BXgDeJeijcM41QeIFd+WV/lFy2SOZRJFvZpdEXhLnWjgf68llq2T6KJ6SdIGkHpI6x2+Rl8y5FsZHj3XZKpkaxWnh/1clLPPTY51rpGp9FD6Eh8siyVyZ/ZXZ7CQVRFMc51ouvzLbZaukB/dT4BBJdwIrIiyTcy2SX0fhslUy81FMkHQLsBx4DHgZGBp1wZxraapfR+FnPbnsUWeikPQrSYuBXwLvAGOA1WZ2j5mtT1cBnWspvEbhslV9fRTnAR8Afwb+ZWZfSvJ3t3NN5Gc9uWxVX9NTD+AXwNHAh5L+DrSRlNSIs865GrxG4bJUnV/6ZlYJPE0wcVEhcBTQBlgp6XkzOz1NZXSuRfAahctWyc5H8SXwCPCIpBLguEhL5VwLVK2PwofwcFmk0c1IZrYJuDeCsjjXoiVecOc1CpdNkr6Owjm3a2Jm5OcK8D4Kl108UTiXJmZGXk7wkfMahcsmyVxw11bSf0m6I7w/SNJR0RfNuZbFwGsULislU6P4G8GERRPD+ysJTpt1zjVCzIy83Bwkn4/CZZdkEsXXzOxGoBzAzLYBirRUzrVAMYMcQX5ODuU+eqzLIskkih2S2hBeLiTpayRMiVofSYdLel/SEklX1rPeiZJM0tikSu1cFjIDSeTnymsULqskc3rsdQQX3vWRdB8wCTinoSdJygWmAocCZcAsSY+b2cIa6xUDk4E3G1d057KLmSEgPy+HCk8ULovUW6OQJOA94ASC5DANGGtmM5PY9nhgiZl9ZGY7gPuBY2tZ7wbgN8D25IvtXPYxgxyJ/NwcdnjTk8si9SYKC2ZamW5ma83sSTN7wszWJLntXlSft6IsXFZF0l5AHzN7sr4NhVOxzpY0e/Xq1UmGd655iZmRIyjIzfGmJ5dVkumjmCtpXKoDS8oBfgf8tKF1zex2MxtrZmNLS0tTXRTn0iLmfRQuSyXTR7EPcIak5cBWgjOezMxGNvC8lUCfhPu9w2VxxcBwYGbQwkV34HFJx5jZ7CTL71zWMAwJ8nJzqPCmJ5dFkkkUhzVx27OAQZIGECSIU4GqEWfNbCPQNX5f0kzgUk8SrqWq3kfhNQqXPZJperI6bvU/yawCuBCYASwCHjSzBZKmSDqm6UV2LjvFLKhRFHjTk8syydQoniRIDAKKgAHA+8CeDT3RzKYD02ssu7aOdQ9KoizOZa3EGoUnCpdNGkwUZjYi8X54ptKPIiuRcy1UZVijCBKF91G47NHo0WPNbC5BB7dzrhFiMSNXIs+bnlyWabBGIemShLs5wF7AJ5GVyLkWqjJm5OaIgtwcdlR4onDZI5kaRXHCrZCgz6K2K6ydc/UILrgT3UqKWLnhC8y8+cllh2Q6sxea2UOJCySdDDxUx/rOuVrEaxS7d2vPhm3lrN26g67tCzNdLOcalEyN4qoklznn6lFpkBMmCoAln2/JcImcS06dNQpJRwBHAr0k3ZLwUAlQEXXBnEunWMy489WlfHt8H0qK8iOLkSsYlJAoJgzsEkks51KpvhrFJ8BsglFd5yTcHqfpV2s71yzNK9vAL6cv4tIH344sRkUsRl5ODj06FNGuINdrFC5r1FmjMLO3gbcl/a+ZlQNI6kQw2uv6dBXQuXQoyssFYO7HGyKLEYtBTk4wMODXurX3ROGyRjJ9FM9KKpHUGZgL3CHpfyIul3NpFQvPQFqzJanJG5uk0oLObIDdSz1RuOyRTKLoYGabCCYvutfM9gEOibZYzqVXZWznqapRnbZaGQtOjwXYfbf2fLZpO5u2l0cSy7lUSiZR5EnqAZwCPBFxeZzLiIqERLE6olpFLKFGMbhbMQCLV3mtwjV/ySSKKQQjwC4xs1mSBgKLoy2Wc+mVWKP44LNovrwrwyE8AIb2CBLFe59tiiSWc6nUYKIws4fMbKSZ/Si8/5GZnRh90ZxLn4rYziE1Pli1OZIYlTEjJ6xR9OrYhuLCPN77NJpYzqVSfddRXG5mN0r6I7XMP2FmF0daMufSKCFP8P5n0SWKvDBRSGJoj2KvUbisUN8QHovC/33GOdfixWsUbfJzI/vyrrSdNQqAod1LePStlZgZ4XTAzjVL9V1H8a/w/3sAJLUP73vvm2tx4n0UI3p34O0VG6iojJGX2+hR+OsVS+ijgKCfYvMbFazc8AW9O7VNaSznUqnBT4Kk4ZLeAhYACyXNkdTg7HbOZZP4WU8je3Xgy4oYH63ZmvIYiddRQFCjgOiaupxLlWR+Mt0OXGJm/cysL/BT4I5oi+VcesUSahQACz9JffNTZWX1RDGke3FksZxLpWQSRTszezF+x8xmAu0iK5FzGRCvUQzpXkxBXg4LPtmY8hjlMSM/d2eiaF+Yx8Cu7Xh3ZepjOZdKySSKjyT9l6T+4e0a4KOoC+ZcOsX7KArzchnavZiFn6b+V35FZTAoYKI9e3VgvicK18wlkyjOBUqB/wtvpeEy51qMeI0iV2LPniUs+GRTyofyKK808mt0kI/oVcInG7ezNsIxppzbVclccLfezC42s73C22QfPda1NJXh6bG5uWLPnh3YsK2csvVfpDRGeWWsWtMTwPBeQZ+INz+55qy+C+4er++JZnZM6ovjXGZUhhfc5eWI0X06AvB22Qb6dE7daasVMSOvjkQxf+VGDhrSLWWxnEul+i64mwisAKYBbwJ+RZBrsapqFDmq6tB+e8UGjhrZMyXbj8UsvDK7eiW+pCif/l3aeo3CNWv1JYruwKHAacDpwJPANDNbkI6COZdO8T6KvByRn5vD8J4lvL0idV/e5WEiKsj7amvvyN4dmbVsXcpiOZdqdfZRmFmlmT1tZt8BJgBLgJmSLkxb6ZxLk/hZT/EhNkb27si7KzdSURmr72lJq6jcmYhqGtO3I59u3M6nG1PbJ+JcqtTbmS2pUNIJwD+A/wfcAvwzHQVzLp0qY9W/yEf36cgX5ZUsTtEsdFWJopZhQfbq2wmAucujm4bVuV1RZ6KQdC/wOrAX8HMzG2dmN5jZyrSVzrk0qTo9NiFRAMz9ODUn+MWbnmqe9QQwrEcJhXk5vJWiWM6lWn01ijOBQcBk4N+SNoW3zZJ8zAHXouysUQQfiX5d2tK1fQFzlqUoUfjJlyIAAB1QSURBVFTGE8VXP3IFeTmM6NUhZUnJuVSrr48ix8yKw1tJwq3YzErSWUjnohavUcS7ECQxtl9nZi1PTSfzo299AsDKOq7N2KtfJ+av3MSXFZUpiedcKqV2HGXnslRlLEZujqrNCzG2fydWrPuCVZu27/L2H5sXtNgu/rz2kWL36tuJHZUx3i3z02Rd8+OJwjmCC+5ya5yRNK5/ZwBmp6D5qaJG01ZN4wcEsd5c6qfJuubHE4VzBDWKmqeu7tGzhDb5uSm5xiHeB1LXRHad2xUwZLdi3vho7S7Hci7VPFE4R/CLv2aNIj83h737dUrJl/dBQ0oBmLR71zrX2WdgZ+YsX1/V8e1cc+GJwjmC6xxqJgqAiV/rwnufbWbNLo7uOjKcEGnCwC51rrPPgC5s21HJAp/IyDUzkSYKSYdLel/SEklX1vL4JZIWSnpH0vOS+kVZHufqsqMiRmEtw2vEawC7WquIj1he34Bp8X6K1z/05ifXvESWKCTlAlOBI4A9gNMk7VFjtbeAsWY2EngYuDGq8jhXny8rKinKz/3K8uE9SyguzOO1JSlKFPVkitLiQoZ2L+bVJat3KZZzqRZljWI8sMTMPjKzHcD9wLGJK5jZi2a2Lbz7BtA7wvI4V6ft5bXXKPJyc9hnYGf+/eGaXdp+fAokNTAI8wGDS5m1dD1f7PDrKVzzEWWi6EUwTHlcWbisLt8DnqrtAUkXSJotafbq1f5ry6XWjooYy9ZurbVGAbDf7l1ZvnYby9dubXKM+Gx59dUoAPYf1JUdlTHeXOrNT675aBad2ZLOBMYCN9X2uJndbmZjzWxsaWlpegvnWryf/fNd3vtsM3XNfBqfUGjm+03/kZLspKrj+nemIC+HVxbvWg3GuVSKMlGsBPok3O8dLqtG0jeAnwHHmJlPHOzS7qUPggSwdUdFrY/379qOAV3b8eL7nzc9SBJ9FABF+bnsM6DzrsVyLsWiTBSzgEGSBkgqAE4Fqk2vKmkM8BeCJOGfDJcRHdrkA7Dty7r7BQ4aUsrrH65tct+BEW96aniiyG8M242PVm/lo9WpGeLcuV0VWaIwswrgQmAGsAh40MwWSJoiKT7f9k1Ae+AhSfMamqfbuSiUhIli65e11ygAvj6kG19WxJrcqZ3M6bFxhwwLmrqeX+S/nVzzUN9UqLvMzKYD02ssuzbh729EGd+5ZMRrFJvrSRQTBnahuDCPGQs+45BhuzU6RtVZT0lkit6d2jK0ezHPLlrF+QcMbHQs51KtWXRmO5dJ7Qob/r1UkJfDIcO68ezCVU2aHjWoURjt5t4B65c1uP43hu3G7GXrWLd1R6NjOZdqnihcq9cmP/gYXNDAr/fDh3dn/bZy/tOEQQINowfrKHnpv2D6ZQ2uf/jw7sQMZiz4rNGxnEs1TxSu1auMwW4lhVx1xNB61ztgcClF+Tk89W7jv7zNoEhh7WDthw2uv2fPEgZ0bceT73za6FjOpZonCtfqlVfGaFuQ1+AZSW0L8jhk6G5Mf/fTRo/wakA+YR9IbkGD60viWyN68O8P1+zygITO7SpPFK7VK6+MkZ+bzPlIcOzonqzduoNXlzTu7CczS0gUyZ1D8q2RPYgZTH/XaxUuszxRuFYvSBTJfRQOGtKNDm3yeeytr1w7Wi8zKKQ8uJNEjQJgaPdihnYv5pG5jYvlXKp5onCt3o5KSzpRFOTlcOSI7jy7cFWjmp/MjAIl3/QEQfPTSXv35u0VG1i8qva5tp1LB08UrtUrr4hRkGSiADhwcClbd1TyTtmGpJ8Ts8Q+ivykn3fcmF7k5YiH55Ql/RznUs0ThWv1yitj5Ocl10cBwUx0Evy7EXNUGFDQyKYngK7tC/n60G48PKeMLyt86HGXGZ4oXKvXmD4KgE7tChjWvYTXGzHrnZlR0IiznhKdNaEfa7fu4On5fk2FywxPFK7Va0wfRdyoPh1Y9OmmqnkmkrGzRpF80xME82EM6NqOe19f3qjnOZcqnihcq1de2bg+CoDBuxWzfls5q5O8xiHWhM7suJwccdaEfsxZvp63Pl7fqOc6lwqeKFyr15jrKOIG71YMwOJVyQ0FXu302JzG1SgAThnXh5KiPG5/+aNGP9e5XeWJwrV65RWN66MAGLRbewA+SPK01aAzu/FnPcW1L8zjrIn9eHrBZz5PhUs7TxSu1dtRaeTnNe6jUNq+kE5t85NOFDGznX0UeYWNLSIA5+w7gMK8HG59YUmTnu9cU3micK1eU/ooJDGkezFvr9iY1PpmjRvrqTalxYWcPbE/j85byZLP/QI8lz6eKFyr15Q+CoCDh3Zj4aeb+HjttqTWr+rMTmb2ojp8/4CBFOXn8ttnPmjyNpxrLE8UrtVr7HUUcUeO6AHA9PkND9oXiyU0PTXilNqaurQv5PsHfI2n5n/Gf5Y2fl4M55rCE4Vr1cyM8iZcRwHBlKW9OrZJqp/CLMaYnNT0LVxwwEB6dChiyhMLqIw1Pek4lyxPFK5VK68MvmgLGtmZHde1fQGrNzd8LcXYj+9ifM77wZ1dqFEAtCnI5aojhzF/5SbufX3ZLm3LuWR4onCtWnwE2Kb0UQDs0bMDby5d12A/xW5bFiXc2/VawNEje3Dg4FJumvE+ZeuT6yNxrqk8UbhWbWeiaNpHYfIhg8jLETc8ubDe9WQJA/rtYo0CgrOufnHccARMvn9eo2fcc64xPFG4Vm3HLiaK7h2KuPDg3Xl24Spe+mB1nevJEr/IU9Ov0KdzW351wgjmLF/P/zzrZ0G56HiicK1aVR9FExMFwPf2G0D/Lm35xRN11ypEamsUcceO7sVp4/vwp5kf8nI9icq5XeGJwrVq5RVhjaIR81HUVJiXy5kT+rH48y2s27qj1nWiqFHEXXvUngzZrZgfPzCPpWu2pnTbzoEnCtfKxdv228S2wht/hu2bmrSdYT1KALjz1doH7Ut1H0WiNgW5/PnMvQA46843WbVpe0q375wnCteqxfso+q14FJ6+Et76e5O2M3FgF04d14epL37I1Be/er1ElDUKgIGl7bn7u+NYv3UHZ9/5HzZuK095DNd6eaJwrVpVH0UsvBZic9NmkcvJEb88fgTHje7JTTPe5/aXP6z2eJQ1iriRvTty+9ljWbpmK2fe+WadzWDONZYnCteqxZue8uKdzeVNvyYhN0fcfPIovjWiB7+a/h7XPjafivhpqxElh5om7d6V287aiw9Wbeacv/0nLTFdy+eJwrVq8c7sgopwjodtyc+DXZu83BxuOW0MFxwwkHtfX853754VNAMl1igiaHpKdPDQ3bjw67vzTtlGtnxZEWks1zp4onCtWtV1FPFEsXnVLm8zN0dcfeQwbjxxJG98tJajb32Vyh1f7FwhDbWLgaXBxEor1vlV227XeaJwrVq8jyK/PBzYb0vT+ihqc8q4Ptx/wQR2VMTILU+clS76RNG3c1uAei8CdC5Znihcq1ZeUc51effQYemTwYIU1CgS7d2vM09evB8dcxMGDkxDjWKPniXsP6grv37qPU8Wbpd5onCtWuHGZXw3b8bOBeVb4cvUzh7XpX0h7UloekpDjSI3R/z1O2Pp1bENtzy/GEtTZ7prmTxRuFZNX9Qy+c+Wz1MfKJbQqfzWffDZu6mPUUNhXi7fP3Agc5av55bnl7DeT5d1TeSJwrVquds3fHVhE6+lSJpVwrTT09IEdeq4vuy3e1f+57kPOOjmmUz7z8d+fYVrNE8UrlXL/bKWRJHCDu06bfwYPl/U8Hq7qCAvh3+ctw9PXLQf/bu246r/e5e9bniW7/99NsvX+rhQLjl5UW5c0uHAH4Bc4K9m9usajxcC9wJ7A2uBb5vZsijL5Fyi3O0J100M+Ra8/2TKO7TrtPgZ2G2PtIQa3qsDj/5oX+at2MBzi1Zx16vLmLFgJr06tuGMCX0Z2asj3UoK6VZcSIc2+UhNHyTRtTyRJQpJucBU4FCgDJgl6XEzSxyL+XvAejPbXdKpwG+Ab0dVJudqylm/lA20p+Npd0Kf8fDbZ6OpUYw4Bd59EM57AYpK4KHvBolivx+nPlYdJDGmbyfG9O3EmRP68eQ7n/L8os+58en3q61XkJdDaftCurYvoKRNPsVFeZQUVf+/bUEehfk5FOblUJifG/yfl0tRfvB/YV4OBXk55OaIXImcHJGXo+B+wjKXHaKsUYwHlpjZRwCS7geOBRITxbHA9eHfDwO3SpJFcIrGrP/7A6Xz70j1Zon2rR5NG7ai2m4WnlkzPLaWFQVfo+OQw4MF7XeD2X+D954EBFW/rBP/JqF/wWrcr2PZllXQZRD03ju4P+hQeO33MHVCancoST2A88Jbea8Y5ZVGRWWMiphRGbPg/41GbL1RaUYsZsTMiCVxiCuBbeGtIQr/Cf5XrZ+nWj9jauDxhEfU8IpZY93ek9n7W+elPW6UiaIXsCLhfhmwT13rmFmFpI1AF2BN4kqSLgAuAOjbt2+TCpPXvgvr2g5o0nMbFt070KLadqRNC9FsO4rXYhWQN/KknQsOvBw+egmw8Iu+xv/VXrfw72SXDTxo56K9vwMbPobKzHcs54e3ZMQMKmIxKmPxBEK1RFIZI/zfMDOM4GWr828s/H/nciB4uWvEthoPfOVxq7m+Jd5Jw0nJVaEiU9C+c4Rbr5uiOr9a0knA4WZ2Xnj/LGAfM7swYZ354Tpl4f0Pw3XW1LZNgLFjx9rs2bMjKbNzzrVUkuaY2dimPDfKs55WAn0S7vcOl9W6jqQ8oANBp7ZzzrlmIspEMQsYJGmApALgVODxGus8Dnwn/Psk4IUo+iecc841XWR9FGGfw4XADILTY+8yswWSpgCzzexx4E7g75KWAOsIkolzzrlmJNLrKMxsOjC9xrJrE/7eDpwcZRmcc87tGr8y2znnXL08UTjnnKuXJwrnnHP18kThnHOuXpFdcBcVSauB5UBXalzBnWatOX5r3neP7++9bI3fz8xKm/LErEsUcZJmN/UqQ4+fvbE9fuuO35r3PZPxvenJOedcvTxROOecq1c2J4rbPX6rjO3xW3f81rzvGYuftX0Uzjnn0iObaxTOOefSwBOFc865erX6RKEMziKfidiZ3N+a8TNdlkzK9L5nOr7LnKYc+1aZKCR9U9KvANI9/4Wk3SR1jcfOwAe2XY3ypPs90ElSLlTtf1rjZ2B/E2N3DCfoysixl9RVUvsMxh8sqSidMWvEHyNpXIZiHyzp+5mIHcY/WtLfoGnfea0uUUj6JnAbsI+kQWmOfQTwNHCrpNsgvR9YSYcBD0u6RtK1YfxYur48JR0DPEew/3fE46cjdhj/YOB0SZ3SFTMh9mEEE3X9WdL/QHp/pEg6HHgCuEXS7RmI3xd4D/h/GXr9Dwf+BmyvsTzyz174vv8jNWb4TOPn/lDgRmCkpG80aSNm1mpuwGHAHOB44C5gchpj7wMsBL4B7A78BchJeDwn4vjjgfeBI4FhwDPAg2mMPwh4G/g60JNgQqvHgfZpij8JiAHPAt8GOqXx2H+D4EvyaGBv4CHg9DTHXwgcAQwF/hdom65jH8bYLXwNngN+AnRM4/4fTPAlPS68X1jj8cj2HygE/g4cGN5vD3RO475/M/zcHQ5cAVzblO1EOnFRcxFm7lLgYuASM3tJ0nrgT5JeMbO5aShGW+BhM3tO0hCCD+1/S+pmZt+14Je9LDy6EcV/0ILJpAirobdIut/MTrXof9lvABYDi8zsM+AwSf8LTAOOjnL/wxpTJ+CUcNGxQI6kp81sfbhOVLHbAPsDV5jZv8L34lyCZBm5MP7ewA/D9/1wYF/gJ5I6mdmlaXjvYWarJP2V4D1wIbBJ0iJgvZktiiqupHxgDPAu8JmkjsDvJW0E2pnZeRHvfyXQAciVtBvwALAhbIK7xMwWRvje60jwnr/QzF6RtBZ4IvzOe7Ex22oVTU8W+Bw4I/yw5AP/JvhFOxwg3m4eoXLgJEnXE/yqugP4K9BX0iPxckYYX8D3JI0I7+8OXArkS/p2JAHDqnXYLl8JrAeqxqkxs9OBQkm3hvdTuv/x+GESfBl40sweJmj+OwI4UlLniGN/QfCLclbCF8IighpeZGrE/3P4vi8BriZIzv8ERkf13ks49jkJTZv9CN6HxwI/BF4B+qQybs34ZlZOUIN7lKD5ZQHB638/MFDSA+F6Kdv/hH3PNbMKgtd6OHAN8A8zO46ghve7VMdOjG9mGwiS0SuS8sxsFsFrsJ+kvMY0Obf4RCFpgqTvSxoFtAkXV5jZDuBD4HJJ7cysMqLYF0gaZWavAmcArwMvmtkNZrYYOAH4IvzlF2X8Fwm+JB6SdB9Bdfxe4CWCXzxR6A7B/Olmto4gQf5c0n4J6/wIqIgyfliGTYTt02Z2H0ET1GHA3pIul/TfEcZeYmafJHwhfEkwCiiSzpI0OcWxq8UHNof/bwd+aWZXmdlC4ExgS/jDKZL4YZKOt8X/neBkim4ETVHzgN3DBBZJ/LAMHxNMyTyPYP9/Y2avE+z/tgh+JMb3Pf6dsgQ4BBhMkCAws0sAkzQwxbGr4oe2hvHin7H3Cd73neI1qaS2GHUbWSZvwLeAZcCfgXsIOrPG1ljnb8DPCa9Sjyj2P8I4e4WPTQf2CP/+DvAqYVt9xPF7ESSFPkBeuN7VwBSCD3PKXgOCX40x4NIay88DZhMkyL7AuQS1u7apil1XfIIfRkq4fyhB++0KYEy6YgNDgD+Ex2hW/L0Qcfy8WtY7j6CvKh2vvQiawFYAn7Gzr+pRoEvU8cPlHRJfB+D8VO9/PbG/SfBD6TKCWvXx4Xsv8n2v+bkmGAbkfiA36e2mspDN7QZcT9DcFP9wXgg8CeydsM55wP8ARRHHviiM3Rf4LrCO4Ev8XWDPNOz75DD+2IR1zg0/uENTHLs3QTK8AigDLqvx+IkEZ4H8H/AfYGS64tf4wj6e4Nd2yr6oG4od/t8//DDPS/WxTzJ+UfgenJfKfU/y2P8I+FbC/VR/7urb/8QfCd8B3knXsQ8fPxD4MUFN/qkMvO/jx39/4BagOOltp7Kgze0G/Aq4I+F+V4Jk8RegNFxWAvRIU+yLgb+E9w8ADgIGZGDfuxB07l8NDIsgdg5wWPj3HsDqWj40hQQdzN0yET987LAIviiT2feeBLXIwRl67fsS9JGldN8b+drnkuJabCP2v3v4RZn2Yx8+1oYIzvpqRPx2QNdGbTvVhW1ON6AjQSdmYjVsOEG1K+Uf0iRjPxR17GT3HSiIIK5q/k1wOm7Vm5bgVOHeEe13svG7ZzB2EcEZN5mKX5yBY395QvyU/zBr5P6XAPkZip3p932T4qe8wM3lxs5q1gSCsw6uSHjsPiI8jz2TsZOMH2+SSumvuTrKkh/+P4ygOvwUQbt8JF8WzSl+HbHnZHjfMx0/08d+dgaPfab3vcnxIy9wJm7UuIAGGAe8QHCR3VUEZztF1eSTsdjNKX5tSYig03wNMKIlxm/N+97a47f0fW8RF9yFl6gfStBB/A8zKwvPYa4MLzDqBBxD0JFmwHFmtjTbYzfj+DkWnHo3kiApPRaenjwRONjM3m0J8Vvzvrf2+K1u36PKcOm6EZxi+BbwU+BPwAUJjx1CcFbRQS0tdhbFnxTez6GRHWjNOX5r3vfWHr817ntWz3AXXijzJ+ARM3tGweiMAwna4xYTdN7kmdmDqb5MPpOxszB+rqX4gsZMxm/N+97a47fWfc/2picRnMFwqKTPgUsIOutOJmibP8nMPo7iizLDsbMtfsqves9w/Na87609fqvc96ysUUjqTjBEyipJ/YDfEwwD8YmZTQ7XuQlYZmZTW0psj+/H3uP7sc9E/KyrUUg6keDqxnxJ/wKeNrPjJZ1EMNBdnJHiMYwyGdvj+7H3+H7sMxEfsqxGIakLwXgp5xKMxnoowXnC/wLeCB97HFgO/AA408zey/bYHt+Pvcf3Y5+J+HHZVqPIBTYBS81sg4Lx1b9BMCHMZ8DpwLUEw1V8N8UvWCZje3w/9h7fj30m4gNZVqMAkPQHgrFKJpvZVkk9gbOBSjO7ScHcBzGLYCKeTMb2+H7sPb4f+0zEhyyaj0I7J9mYSpBhr1Awj8QnBNNqHiWpswVzH6T0BctkbI/vx97jZyZ+a973mpp9opCqzVIGwRAU/0cwAuNtkroSTAhSQYonwMlkbI/vx97j+7HPRPxay9Rcm54UTFG53cy2JSwrMLMdknoDnQnGlN8j/PuHlqK5rzMZ2+P7sff4fuwzEb9elsJLy1N1I5il6V8EWfQsqk+2cwjBUN19w/sdSOGQzZmM7fH92Ht8P/aZiN9g+dIZLMkXbDDBWCV7EEzucxPBZPD7A/kEp4Sd2NJie3w/9h7fj30m4idza46nx3YFyiyY/B1JnQhGPj2FYBKOYy24OjGKoSkyGdvj+7H3+H7sMxG/Qc2xM3s+sEnSNeH9McAHwHagv5mtguBa9hYW2+P7sff4fuwzEb9BzaIzO+yo2WxmGxWcE/x14P8RJDIzs2MlnQYcDpyTyhcsk7E9vh97j+/HPhPxG80y2O4V7vtxwHsEoyCWJizPIZgEPTe8fwHwh5YS2+P7sff4fuwzEb8pt4zWKCSVAvcDHxPM6/o5cL+Zra6x3o+B7xKMY5KqGaoyFtvj+7H3+H7sMxG/qTLdR7GRYFTEHwDzCEZCPFVSN9h54QnBeCdnpPgFy2Rsj+/H3uP7sc9E/CbJSI1CUl+CAa3yrPrFJScCBwKLzeyPkvayFF9QksnYHt+Pvcf3Y5+J+Lsq7TUKSd8CpgO3An+TNDT+mJk9ArwElEp6FHhJUq+WENvj+7H3+H7sMxE/JdLVGQII6ENwYclBwG7ApcCnwJ411v0HsAwYke2xPb4fe4/vxz4T8VN5S2+woN3tdqAXO5u9JgMrgcHh/R7AQmB0S4nt8f3Ye3w/9pmIn7L9SEuQoMNmHNAFeAC4vMbjlwN3A23C++1bQmyP78fe4/uxz0T8VN+iDwBHAe8QtMPdChxDUMW6KmGd/sBfCDNuS4jt8f3Ye3w/9pmIH8Ut2o3DvsAiYEx4/3bgF0BPgvOIryHIvOcAs4FOLSG2x/dj7/H92GciflS3aDcevGjnJNwvBZ4M/x4I3AX8CZhDijtxMhnb4/ux9/h+7DMRP6pbtBsPOnJKEv7uDbwF9AiX9QPygA4tKbbH92Pv8f3YZyJ+VLdIr6Mws0oz2xTeFbABWGdmn0o6E7gayDezjS0ptsf3Y+/x/dhnIn5U0n5ltqS7Cc4j/iZBFS1tl6hnMrbH92Pv8f3YZyJ+KqQtUYRjmOQTdPTkA4eY2eKWHtvj+7H3+H7sMxE/lTJRozgHmGVmC9IaOMOxPb4fe4/vxz4T8VMhE4kiY9P5ZTK2x/dj7/H92GerZjHDnXPOueYr0/NROOeca+Y8UTjnnKuXJwrnnHP18kThnHOuXp4onKuFpEpJ8yQtkPS2pJ9KqvfzIqm/pNPTVUbn0sUThXO1+8LMRpvZnsChwBHAdQ08pz/gicK1OJ4onGuAmX0OXABcqEB/Sa9Imhve9g1X/TWwf1gT+YmkXEk3SZol6R1J3weQ1EPSy+F68yXtn6l9cy4Zfh2Fc7WQtMXM2tdYtgEYAmwGYma2XdIgYJqZjZV0EHCpmR0Vrn8B0M3MfiGpEHgNOBk4ASgys19KygXamtnm9O2dc42Tl+kCOJeF8oFbJY0GKoHBdaz3TWCkpJPC+x2AQcAs4C5J+cCjZjYv6gI7tys8UTiXBEkDCZLC5wR9FauAUQTNt9vrehpwkZnNqGV7BwDfAu6W9DszuzeSgjuXAt5H4VwDJJUCtwG3hmP2dAA+NbMYcBbBBDUQNEkVJzx1BvDDsOaApMGS2knqB6wyszuAvwJ7pWlXnGsSr1E4V7s2kuYRNDNVAH8Hfhc+9ifgEUlnA08DW8Pl7wCVkt4G7gb+QHAm1NxwyOnVwHHAQcBlksqBLcDZadgf55rMO7Odc87Vy5uenHPO1csThXPOuXp5onDOOVcvTxTOOefq5YnCOedcvTxROOecq5cnCuecc/X6/4QAiwf3uhmkAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.plot(dates, moisture, label='Moisture')\n",
"plt.plot(dates, precip, label='Rain Reading')\n",
"plt.xticks(rotation=45)\n",
"plt.title(\"Initial Moisture Decay Testing\")\n",
"plt.xlabel(\"Dates\")\n",
"plt.ylabel(\"Moisture Amount\")\n",
"plt.legend()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Woo Hoo!\n",
"\n",
"Okay, that worked pretty well. The rain is decaying as expected. If it rains, the rain gets piled on top of what is already there. Very high amounts of rain decay quickly and smaller values stick around for a while. This is roughly what I've seen. It'll be refined further later, but it's a start. There is one issue with that data input - see if you can figure out what it is. The graph above has the issue in it. Have you identified the tiny giant problem? \n",
"\n",
"\n",
"## The Tiny Giant Problem\n",
"\n",
"The issue is that this is time series data. Which means the time steps have to be very exact. If we have one reading now, one reading in an hour and another reading in a day, this equation becomes completely meaningless. There are some points in the graph above that are possibly displaying this type of artifact. The decay value around '2019-08-13' drops off very quickly, suggesting more data points were around this time period resulting in inconsistent decay. I'll fix this in the next post. \n",
"\n",
"### THANKS FOR READING"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.8"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment