diff --git a/tutorials/site-monitoring-foundations.ipynb b/tutorials/site-monitoring-foundations.ipynb new file mode 100644 index 0000000..0f65ec3 --- /dev/null +++ b/tutorials/site-monitoring-foundations.ipynb @@ -0,0 +1,2447 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "623fe031", + "metadata": {}, + "source": [ + "### Getting Started with Time-Series Monitoring in the Planetary Computer\n", + "\n", + "
\n", + "📘 Outline\n", + "\n", + "1. Introduction \n", + "2. Learning Objectives \n", + "3. Core Concepts\n", + "4. Environment configuration\n", + "4. Load/Define Your Area of Interest \n", + "5. Define Monitoring Conditions \n", + "6. Query the Planetary Computer's STAC Catalog \n", + "7. Apply NDSI Index \n", + "8. Visualize Time Series \n", + "9. Explore and Adjust Parameters \n", + "10. Export Results \n", + "11. What's Next \n", + "\n", + "
\n", + "\n", + "Monitoring landscapes over time can help detect crop stress, measure snow cover changes, or identify wildfire recovery. This notebook is your starting point for learning how to track those changes using Microsoft’s Planetary Computer platform. In this foundational notebook, you'll define an area of interest, query satellite imagery from the Open Planetary Computer (OPC), apply indices like NDSI, and identify changes across a time series.\n", + "The goal is to help you understand key remote sensing concepts and build hands-on skills for real-world geospatial monitoring.\n", + "\n", + "➡️ Next notebooks in this series:\n", + "- [Vegetation Monitoring (Intermediate)](/tutorials/site-monitoring-hls.ipynb)\n", + "- [Monitoring Extreme Weather Events (Advanced)](/tutorials/climate-risk.ipynb)\n", + "\n", + "### Learning Objectives\n", + "\n", + "By the end of this notebook, you should be able to:\n", + "\n", + "- 🧭 **Understand the purpose** of satellite-based site monitoring.\n", + "- 🛰 **Define and import an Area of Interest (AOI)** for analysis.\n", + "- 📦 **Retrieve** Earth observation data using the STAC API.\n", + "- 🧮 **Calculate** remote sensing indices (like NDSI) to track change.\n", + "- 📈 **Visualize and interpret** changes in in your data over time.\n", + "- 🔍 **Detect** significant changes using thresholds or time-series trends\n", + "- 🛠 **Export your results** to resuse in tools like Fabric or Power BI.\n", + "\n", + "### Core Concepts\n", + "\n", + "Before we begin, here are a few key concepts we'll use in this notebook:\n", + "\n", + "- **Raster data**: Satellite imagery is made up of pixels (grids), each with numeric values representing the Earth's surface.\n", + "- **Bands**: Satellite images contain multiple spectral bands—each measuring light reflected at specific wavelengths (e.g., Red, Near-Infrared).\n", + "- **NDSI**: A snow index calculated using the green and shortwave infrared bands. Higher values typically indicate the presence of snow or ice.\n", + "- **Area of Interest (AOI)**: A polygon or boundary that defines where you're monitoring.\n", + "- **Timeseries**: A sequence of observations (images) over time for your AOI.\n", + "- **[STAC](https://stacspec.org/en)**: A standard for searching and describing satellite imagery and geospatial data.\n", + "- **Planetary Computer**: A catalog of open datasets we’ll query using STAC.\n", + "\n", + "Want to explore the catalog? Check out the [Planetary Computer](https://planetarycomputer.microsoft.com/)." + ] + }, + { + "cell_type": "markdown", + "id": "1dfce91c", + "metadata": {}, + "source": [ + "### Configuring your Notebook Environment\n", + "\n", + "You'll need the following dependencies installed to follow along the tutorial (this may take a minute!):" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "db6baeb1", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Requirement already satisfied: rioxarray>=0.19.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (0.19.0)\n", + "Requirement already satisfied: ipyleaflet>=0.18 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (0.19.2)\n", + "Requirement already satisfied: ipywidgets>=8.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (8.1.7)\n", + "Requirement already satisfied: jupyterlab-widgets>=3.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (3.0.15)\n", + "Requirement already satisfied: jupyterlab>=4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (4.4.2)\n", + "Requirement already satisfied: geopandas>=1.0.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (1.0.1)\n", + "Requirement already satisfied: ipykernel>=6.29.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (6.29.5)\n", + "Requirement already satisfied: pystac>=1.13.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (1.13.0)\n", + "Requirement already satisfied: pystac-client>=0.8.6 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (0.8.6)\n", + "Requirement already satisfied: planetary-computer>=1.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (1.0.0)\n", + "Requirement already satisfied: odc>=0.1.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (0.1.2)\n", + "Requirement already satisfied: dask>=2025.4.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (2025.4.1)\n", + "Requirement already satisfied: odc-stac>=0.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (0.4.0)\n", + "Requirement already satisfied: packaging in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rioxarray>=0.19.0) (25.0)\n", + "Requirement already satisfied: rasterio>=1.4.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rioxarray>=0.19.0) (1.4.3)\n", + "Requirement already satisfied: xarray>=2024.7.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rioxarray>=0.19.0) (2025.4.0)\n", + "Requirement already satisfied: pyproj>=3.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rioxarray>=0.19.0) (3.7.1)\n", + "Requirement already satisfied: numpy>=1.23 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rioxarray>=0.19.0) (2.2.5)\n", + "Requirement already satisfied: branca>=0.5.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipyleaflet>=0.18) (0.8.1)\n", + "Requirement already satisfied: jupyter-leaflet<0.20,>=0.19 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipyleaflet>=0.18) (0.19.2)\n", + "Requirement already satisfied: traittypes<3,>=0.2.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipyleaflet>=0.18) (0.2.1)\n", + "Requirement already satisfied: xyzservices>=2021.8.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipyleaflet>=0.18) (2025.4.0)\n", + "Requirement already satisfied: comm>=0.1.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipywidgets>=8.0) (0.2.2)\n", + "Requirement already satisfied: ipython>=6.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipywidgets>=8.0) (8.36.0)\n", + "Requirement already satisfied: traitlets>=4.3.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipywidgets>=8.0) (5.14.3)\n", + "Requirement already satisfied: widgetsnbextension~=4.0.14 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipywidgets>=8.0) (4.0.14)\n", + "Requirement already satisfied: async-lru>=1.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (2.0.5)\n", + "Requirement already satisfied: httpx>=0.25.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (0.28.1)\n", + "Requirement already satisfied: jinja2>=3.0.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (3.1.6)\n", + "Requirement already satisfied: jupyter-core in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (5.7.2)\n", + "Requirement already satisfied: jupyter-lsp>=2.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (2.2.5)\n", + "Requirement already satisfied: jupyter-server<3,>=2.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (2.15.0)\n", + "Requirement already satisfied: jupyterlab-server<3,>=2.27.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (2.27.3)\n", + "Requirement already satisfied: notebook-shim>=0.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (0.2.4)\n", + "Requirement already satisfied: setuptools>=41.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (65.4.1)\n", + "Requirement already satisfied: tomli>=1.2.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (2.2.1)\n", + "Requirement already satisfied: tornado>=6.2.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab>=4.0) (6.4.2)\n", + "Requirement already satisfied: anyio>=3.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (4.9.0)\n", + "Requirement already satisfied: argon2-cffi>=21.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (23.1.0)\n", + "Requirement already satisfied: jupyter-client>=7.4.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (8.6.3)\n", + "Requirement already satisfied: jupyter-events>=0.11.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.12.0)\n", + "Requirement already satisfied: jupyter-server-terminals>=0.4.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.5.3)\n", + "Requirement already satisfied: nbconvert>=6.4.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (7.16.6)\n", + "Requirement already satisfied: nbformat>=5.3.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (5.10.4)\n", + "Requirement already satisfied: overrides>=5.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (7.7.0)\n", + "Requirement already satisfied: prometheus-client>=0.9 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.21.1)\n", + "Requirement already satisfied: pyzmq>=24 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (26.4.0)\n", + "Requirement already satisfied: send2trash>=1.8.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.8.3)\n", + "Requirement already satisfied: terminado>=0.8.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.18.1)\n", + "Requirement already satisfied: websocket-client>=1.7 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.8.0)\n", + "Requirement already satisfied: babel>=2.10 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (2.17.0)\n", + "Requirement already satisfied: json5>=0.9.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (0.12.0)\n", + "Requirement already satisfied: jsonschema>=4.18.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (4.23.0)\n", + "Requirement already satisfied: requests>=2.31 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (2.32.3)\n", + "Requirement already satisfied: pyogrio>=0.7.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from geopandas>=1.0.1) (0.11.0)\n", + "Requirement already satisfied: pandas>=1.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from geopandas>=1.0.1) (2.2.3)\n", + "Requirement already satisfied: shapely>=2.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from geopandas>=1.0.1) (2.1.0)\n", + "Requirement already satisfied: appnope in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipykernel>=6.29.5) (0.1.4)\n", + "Requirement already satisfied: debugpy>=1.6.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipykernel>=6.29.5) (1.8.14)\n", + "Requirement already satisfied: matplotlib-inline>=0.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipykernel>=6.29.5) (0.1.7)\n", + "Requirement already satisfied: nest-asyncio in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipykernel>=6.29.5) (1.6.0)\n", + "Requirement already satisfied: psutil in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipykernel>=6.29.5) (7.0.0)\n", + "Requirement already satisfied: python-dateutil>=2.7.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pystac>=1.13.0) (2.9.0.post0)\n", + "Requirement already satisfied: click>=7.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from planetary-computer>=1.0.0) (8.1.8)\n", + "Requirement already satisfied: pydantic>=1.7.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from planetary-computer>=1.0.0) (2.11.4)\n", + "Requirement already satisfied: pytz>=2020.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from planetary-computer>=1.0.0) (2025.2)\n", + "Requirement already satisfied: python-dotenv in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from planetary-computer>=1.0.0) (1.1.0)\n", + "Requirement already satisfied: GeoAlchemy2<0.16.0,>=0.15.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (0.15.2)\n", + "Requirement already satisfied: SQLAlchemy<3.0.0,>=2.0.35 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (2.0.40)\n", + "Requirement already satisfied: func-timeout<5.0.0,>=4.3.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (4.3.5)\n", + "Requirement already satisfied: h3<5.0.0,>=4.1.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (4.2.2)\n", + "Requirement already satisfied: igraph<0.12.0,>=0.11.6 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (0.11.8)\n", + "Requirement already satisfied: matplotlib<4.0.0,>=3.9.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (3.10.3)\n", + "Requirement already satisfied: networkx<4.0,>=3.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (3.4.2)\n", + "Requirement already satisfied: osmnx<3.0.0,>=2.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (2.0.3)\n", + "Requirement already satisfied: psycopg2-binary<3.0.0,>=2.9.9 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (2.9.10)\n", + "Requirement already satisfied: pymannkendall<2.0.0,>=1.4.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (1.4.3)\n", + "Requirement already satisfied: scikit-learn<2.0.0,>=1.5.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (1.6.1)\n", + "Requirement already satisfied: scipy<2.0.0,>=1.14.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (1.15.3)\n", + "Requirement already satisfied: tqdm<5.0.0,>=4.66.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc>=0.1.2) (4.67.1)\n", + "Requirement already satisfied: texttable>=1.6.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from igraph<0.12.0,>=0.11.6->odc>=0.1.2) (1.7.0)\n", + "Requirement already satisfied: contourpy>=1.0.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (1.3.2)\n", + "Requirement already satisfied: cycler>=0.10 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (0.12.1)\n", + "Requirement already satisfied: fonttools>=4.22.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (4.57.0)\n", + "Requirement already satisfied: kiwisolver>=1.3.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (1.4.8)\n", + "Requirement already satisfied: pillow>=8 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (11.2.1)\n", + "Requirement already satisfied: pyparsing>=2.3.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from matplotlib<4.0.0,>=3.9.2->odc>=0.1.2) (3.2.3)\n", + "Requirement already satisfied: tzdata>=2022.7 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pandas>=1.4.0->geopandas>=1.0.1) (2025.2)\n", + "Requirement already satisfied: certifi in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pyproj>=3.3->rioxarray>=0.19.0) (2025.4.26)\n", + "Requirement already satisfied: six>=1.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from python-dateutil>=2.7.0->pystac>=1.13.0) (1.17.0)\n", + "Requirement already satisfied: affine in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rasterio>=1.4.3->rioxarray>=0.19.0) (2.4.0)\n", + "Requirement already satisfied: attrs in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rasterio>=1.4.3->rioxarray>=0.19.0) (25.3.0)\n", + "Requirement already satisfied: cligj>=0.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rasterio>=1.4.3->rioxarray>=0.19.0) (0.7.2)\n", + "Requirement already satisfied: click-plugins in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from rasterio>=1.4.3->rioxarray>=0.19.0) (1.1.1)\n", + "Requirement already satisfied: joblib>=1.2.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from scikit-learn<2.0.0,>=1.5.2->odc>=0.1.2) (1.5.0)\n", + "Requirement already satisfied: threadpoolctl>=3.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from scikit-learn<2.0.0,>=1.5.2->odc>=0.1.2) (3.6.0)\n", + "Requirement already satisfied: typing-extensions>=4.6.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from SQLAlchemy<3.0.0,>=2.0.35->odc>=0.1.2) (4.13.2)\n", + "Requirement already satisfied: cloudpickle>=3.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (3.1.1)\n", + "Requirement already satisfied: fsspec>=2021.09.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (2025.3.2)\n", + "Requirement already satisfied: partd>=1.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (1.4.2)\n", + "Requirement already satisfied: pyyaml>=5.3.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (6.0.2)\n", + "Requirement already satisfied: toolz>=0.10.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (1.0.0)\n", + "Requirement already satisfied: importlib_metadata>=4.13.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from dask>=2025.4.1) (8.7.0)\n", + "Requirement already satisfied: odc-geo>=0.4.7 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc-stac>=0.4.0) (0.4.10)\n", + "Requirement already satisfied: odc-loader>=0.5.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc-stac>=0.4.0) (0.5.1)\n", + "Requirement already satisfied: exceptiongroup>=1.0.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from anyio>=3.1.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.2.2)\n", + "Requirement already satisfied: idna>=2.8 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from anyio>=3.1.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (3.10)\n", + "Requirement already satisfied: sniffio>=1.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from anyio>=3.1.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.3.1)\n", + "Requirement already satisfied: argon2-cffi-bindings in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from argon2-cffi>=21.1->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (21.2.0)\n", + "Requirement already satisfied: httpcore==1.* in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from httpx>=0.25.0->jupyterlab>=4.0) (1.0.9)\n", + "Requirement already satisfied: h11>=0.16 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from httpcore==1.*->httpx>=0.25.0->jupyterlab>=4.0) (0.16.0)\n", + "Requirement already satisfied: zipp>=3.20 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from importlib_metadata>=4.13.0->dask>=2025.4.1) (3.21.0)\n", + "Requirement already satisfied: decorator in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (5.2.1)\n", + "Requirement already satisfied: jedi>=0.16 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (0.19.2)\n", + "Requirement already satisfied: pexpect>4.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (4.9.0)\n", + "Requirement already satisfied: prompt_toolkit<3.1.0,>=3.0.41 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (3.0.51)\n", + "Requirement already satisfied: pygments>=2.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (2.19.1)\n", + "Requirement already satisfied: stack_data in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from ipython>=6.1.0->ipywidgets>=8.0) (0.6.3)\n", + "Requirement already satisfied: wcwidth in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from prompt_toolkit<3.1.0,>=3.0.41->ipython>=6.1.0->ipywidgets>=8.0) (0.2.13)\n", + "Requirement already satisfied: parso<0.9.0,>=0.8.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jedi>=0.16->ipython>=6.1.0->ipywidgets>=8.0) (0.8.4)\n", + "Requirement already satisfied: MarkupSafe>=2.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jinja2>=3.0.3->jupyterlab>=4.0) (3.0.2)\n", + "Requirement already satisfied: jsonschema-specifications>=2023.03.6 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema>=4.18.0->jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (2025.4.1)\n", + "Requirement already satisfied: referencing>=0.28.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema>=4.18.0->jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (0.36.2)\n", + "Requirement already satisfied: rpds-py>=0.7.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema>=4.18.0->jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (0.24.0)\n", + "Requirement already satisfied: platformdirs>=2.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-core->jupyterlab>=4.0) (4.3.8)\n", + "Requirement already satisfied: python-json-logger>=2.0.4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (3.3.0)\n", + "Requirement already satisfied: rfc3339-validator in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.1.4)\n", + "Requirement already satisfied: rfc3986-validator>=0.1.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.1.1)\n", + "Requirement already satisfied: fqdn in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.5.1)\n", + "Requirement already satisfied: isoduration in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (20.11.0)\n", + "Requirement already satisfied: jsonpointer>1.13 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (3.0.0)\n", + "Requirement already satisfied: uri-template in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.3.0)\n", + "Requirement already satisfied: webcolors>=24.6.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (24.11.1)\n", + "Requirement already satisfied: beautifulsoup4 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (4.13.4)\n", + "Requirement already satisfied: bleach!=5.0.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from bleach[css]!=5.0.0->nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (6.2.0)\n", + "Requirement already satisfied: defusedxml in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.7.1)\n", + "Requirement already satisfied: jupyterlab-pygments in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.3.0)\n", + "Requirement already satisfied: mistune<4,>=2.0.3 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (3.1.3)\n", + "Requirement already satisfied: nbclient>=0.5.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.10.2)\n", + "Requirement already satisfied: pandocfilters>=1.4.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.5.1)\n", + "Requirement already satisfied: webencodings in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from bleach!=5.0.0->bleach[css]!=5.0.0->nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (0.5.1)\n", + "Requirement already satisfied: tinycss2<1.5,>=1.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from bleach[css]!=5.0.0->nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.4.0)\n", + "Requirement already satisfied: fastjsonschema>=2.15 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from nbformat>=5.3.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (2.21.1)\n", + "Requirement already satisfied: cachetools in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from odc-geo>=0.4.7->odc-stac>=0.4.0) (5.5.2)\n", + "Requirement already satisfied: locket in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from partd>=1.4.0->dask>=2025.4.1) (1.0.0)\n", + "Requirement already satisfied: ptyprocess>=0.5 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pexpect>4.3->ipython>=6.1.0->ipywidgets>=8.0) (0.7.0)\n", + "Requirement already satisfied: annotated-types>=0.6.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pydantic>=1.7.3->planetary-computer>=1.0.0) (0.7.0)\n", + "Requirement already satisfied: pydantic-core==2.33.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pydantic>=1.7.3->planetary-computer>=1.0.0) (2.33.2)\n", + "Requirement already satisfied: typing-inspection>=0.4.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from pydantic>=1.7.3->planetary-computer>=1.0.0) (0.4.0)\n", + "Requirement already satisfied: charset-normalizer<4,>=2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from requests>=2.31->jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (3.4.2)\n", + "Requirement already satisfied: urllib3<3,>=1.21.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from requests>=2.31->jupyterlab-server<3,>=2.27.1->jupyterlab>=4.0) (2.4.0)\n", + "Requirement already satisfied: cffi>=1.0.1 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from argon2-cffi-bindings->argon2-cffi>=21.1->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.17.1)\n", + "Requirement already satisfied: pycparser in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from cffi>=1.0.1->argon2-cffi-bindings->argon2-cffi>=21.1->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (2.22)\n", + "Requirement already satisfied: soupsieve>1.2 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from beautifulsoup4->nbconvert>=6.4.4->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (2.7)\n", + "Requirement already satisfied: arrow>=0.15.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from isoduration->jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (1.3.0)\n", + "Requirement already satisfied: types-python-dateutil>=2.8.10 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from arrow>=0.15.0->isoduration->jsonschema[format-nongpl]>=4.18.0->jupyter-events>=0.11.0->jupyter-server<3,>=2.4.0->jupyterlab>=4.0) (2.9.0.20241206)\n", + "Requirement already satisfied: executing>=1.2.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from stack_data->ipython>=6.1.0->ipywidgets>=8.0) (2.2.0)\n", + "Requirement already satisfied: asttokens>=2.1.0 in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from stack_data->ipython>=6.1.0->ipywidgets>=8.0) (3.0.0)\n", + "Requirement already satisfied: pure-eval in /Users/kcarini/Documents/GitHub/PlanetaryComputerExamples/.venv/lib/python3.10/site-packages (from stack_data->ipython>=6.1.0->ipywidgets>=8.0) (0.2.3)\n", + "Note: you may need to restart the kernel to use updated packages.\n" + ] + } + ], + "source": [ + "%pip install \\\n", + " \"rioxarray>=0.19.0\" \\\n", + " \"ipyleaflet>=0.18\" \\\n", + " \"ipywidgets>=8.0\" \\\n", + " \"jupyterlab-widgets>=3.0\" \\\n", + " \"jupyterlab>=4.0\" \\\n", + " \"geopandas>=1.0.1\" \\\n", + " \"ipykernel>=6.29.5\" \\\n", + " \"pystac>=1.13.0\" \\\n", + " \"pystac-client>=0.8.6\" \\\n", + " \"planetary-computer>=1.0.0\" \\\n", + " \"rioxarray>=0.19.0\" \\\n", + " \"odc>=0.1.2\" \\\n", + " \"dask>=2025.4.1\" \\\n", + " \"odc-stac>=0.4.0\"" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "d1e63d70-3611-4b9b-a42e-3a3728865e34", + "metadata": {}, + "outputs": [], + "source": [ + "import os\n", + "\n", + "import pyproj\n", + "import rioxarray\n", + "import pandas as pd\n", + "import geopandas as gpd\n", + "import planetary_computer\n", + "import pystac_client\n", + "from odc.stac import load\n", + "from ipyleaflet import Map, GeomanDrawControl\n", + "from ipywidgets import widgets\n", + "from IPython.display import display, clear_output, Image\n", + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "id": "738033fc", + "metadata": {}, + "source": [ + "### 🛠️ Fixing EPSG/PROJ errors with `pyproj`\n", + "\n", + "Some environments (especially with software like QGIS installed) may cause conflicts in how spatial reference systems are resolved. \n", + "If you see an error like:\n", + "\n", + "```bash\n", + "CRSError: The EPSG code is unknown. PROJ: internal_proj_create_from_database: \n", + "… contains DATABASE.LAYOUT.VERSION.MINOR = 1 whereas a number >= 3 is expected.\n", + "```\n", + "\n", + "It likely means `pyproj` is picking up an outdated `proj.db` file from another installa tion.\n", + "\n", + "To fix this, we explicitly tell `pyproj` to use the correct internal path:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e10a7ef4", + "metadata": {}, + "outputs": [], + "source": [ + "os.environ[\"PROJ_DATA\"] = pyproj.datadir.get_data_dir()" + ] + }, + { + "cell_type": "markdown", + "id": "60df8bcc", + "metadata": {}, + "source": [ + "### 🗺️ Define Your Own Area of Interest\n", + "\n", + "We start by initializing an empty `GeoDataFrame` to store any drawn geometries. \n", + "An interactive map widget (`ipyleaflet`) is displayed so you can **draw your own** area of interest (AOI) using the polygon tool.\n", + "\n", + "Once a shape is drawn, it is automatically stored in the global dataframe for use in downstream analysis.\n", + "\n", + "👉 ACTION: Draw an AOI\n", + "- Make sure your AOI is **large enough** to contain multiple pixels at Sentinel-2’s 10–30 m resolution.\n", + "- Avoid drawing a region that's too large, as it may slow down processing in this tutorial." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "58d795c0-b006-4284-aa2a-b93f565a47e1", + "metadata": {}, + "outputs": [], + "source": [ + "# global GeoDataFrame\n", + "gdf = gpd.GeoDataFrame(columns=[\"geometry\"], geometry=\"geometry\", crs=\"EPSG:4326\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "bcf607b8", + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "c7c8c9bb437c485792ab6c3cac9d041c", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "Map(center=[45.48510207245395, -73.65836430652902], controls=(ZoomControl(options=['position', 'zoom_in_text',…" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Set up map\n", + "m = Map(center=(45.48510207245395, -73.65836430652902), zoom=10, scroll_wheel_zoom=True)\n", + "\n", + "# Set up draw control\n", + "draw_control = GeomanDrawControl()\n", + "draw_control.polygon = {\n", + " \"pathOptions\": {\n", + " \"fillColor\": \"#6be5c3\",\n", + " \"color\": \"#6be5c3\",\n", + " \"fillOpacity\": 1.0\n", + " }\n", + "}\n", + "m.add_control(draw_control)\n", + "\n", + "# Define callback\n", + "def handle_draw(target, action, geo_json):\n", + " global gdf\n", + " if action == \"create\":\n", + " # Build a GeoDataFrame from the list of features\n", + " new_gdf = gpd.GeoDataFrame.from_features(geo_json, crs=\"EPSG:4326\")\n", + " gdf = pd.concat([gdf, new_gdf], ignore_index=True)\n", + "\n", + "# Attach callback\n", + "draw_control.on_draw(handle_draw)\n", + "\n", + "m" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "54f31066-6ffa-47a1-a927-69591dab00da", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
geometrystyletype
0POLYGON ((-73.91602 45.36951, -73.25409 45.388...{'fillColor': '#6be5c3', 'color': '#6be5c3', '...polygon
\n", + "
" + ], + "text/plain": [ + " geometry \\\n", + "0 POLYGON ((-73.91602 45.36951, -73.25409 45.388... \n", + "\n", + " style type \n", + "0 {'fillColor': '#6be5c3', 'color': '#6be5c3', '... polygon " + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "gdf" + ] + }, + { + "cell_type": "markdown", + "id": "2056334b", + "metadata": {}, + "source": [ + "### 📅 Select a Time Range\n", + "\n", + "Use the interactive date picker to define the time period for analysis.\n", + "\n", + "Make sure your selected range includes months where snow is likely present (e.g., winter or early spring). \n", + "To keep computations efficient, avoid selecting a very long date range — a few weeks to a couple of months is usually sufficient for detecting snow events." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "306d5a5d", + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "ac805aeb2d284604af8ba2d5acbf2b49", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "DatePicker(value=Timestamp('2023-01-01 00:00:00'), description='Date', step=1)" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "01e5f1710af545d0af492955e007fb60", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "DatePicker(value=Timestamp('2024-01-02 00:00:00'), description='Date', step=1)" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "start_date = widgets.DatePicker(\n", + " description='Date',\n", + " disabled=False,\n", + " value=pd.to_datetime(\"2023-01-01\")\n", + ")\n", + "end_date = widgets.DatePicker(\n", + " description='Date',\n", + " disabled=False,\n", + " value=pd.to_datetime(\"2024-01-02\")\n", + ")\n", + "display(start_date, end_date)" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "a9bcb0b5", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "'2023-01-01/2024-01-02'" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "datetime_range = f\"{start_date.value.strftime('%Y-%m-%d')}/{end_date.value.strftime('%Y-%m-%d')}\"\n", + "datetime_range" + ] + }, + { + "cell_type": "markdown", + "id": "4b6224c0", + "metadata": {}, + "source": [ + "### 🌍 Searching the STAC Catalog\n", + "\n", + "We use the Planetary Computer's STAC API to search for Sentinel-2 Level 2A imagery. \n", + "The `pystac-client` library allows us to query the catalog for scenes that intersect our area of interest and fall within a selected date range.\n", + "\n", + "The `get_items` function wraps this search logic, returning all matching items as a collection we can use for loading and analysis.\n", + "\n", + "We then use the first geometry drawn on the map (i.e., the first row of our `GeoDataFrame`) along with the selected `datetime_range` to search for relevant satellite imagery.\n", + "This ensures we query only the area and time window that we interactively defined." + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "id": "1289cada-971e-460a-96bf-51b5d0bbb59f", + "metadata": {}, + "outputs": [], + "source": [ + "catalog = pystac_client.Client.open(\n", + " \"https://planetarycomputer.microsoft.com/api/stac/v1\",\n", + " modifier=planetary_computer.sign_inplace,\n", + ")\n", + "collection_id = \"sentinel-2-l2a\"\n", + "collection = catalog.get_collection(collection_id)" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "id": "5a300ae4", + "metadata": {}, + "outputs": [], + "source": [ + "def get_items(geom, collection_id, datetime_range, max_cloud_cover=100):\n", + " # Get the items from the collection\n", + " items = catalog.search(\n", + " collections=[collection_id],\n", + " intersects=geom,\n", + " datetime=datetime_range,\n", + " max_items=10000,\n", + " query={\"eo:cloud_cover\": {\"lte\": max_cloud_cover}}, \n", + " )\n", + " return items.item_collection()" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "id": "87064c32", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'type': 'Polygon',\n", + " 'coordinates': (((-73.916016, 45.369514),\n", + " (-73.254089, 45.388806),\n", + " (-73.248596, 45.211069),\n", + " (-73.924255, 45.224612),\n", + " (-73.916016, 45.369514)),)}" + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# first geometry of first row\n", + "aoi = gdf.geometry.iloc[0].__geo_interface__\n", + "aoi" + ] + }, + { + "cell_type": "markdown", + "id": "3dad7331", + "metadata": {}, + "source": [ + "> 📝 **Note on Temporal Resolution and Weather Constraints**\n", + ">\n", + "> Sentinel-2 provides high-resolution optical imagery with a revisit time of **5 days** at the equator (with both Sentinel-2A and 2B satellites active). However, the *usable* temporal resolution can be much lower in practice due to weather conditions.\n", + ">\n", + "> Persistent **cloud cover**, seasonal storms, or regional weather patterns can significantly reduce the number of clear observations—especially in mountainous or tropical regions. This makes detecting changes over time more difficult despite the nominal 5-day revisit rate.\n", + ">\n", + "> Filtering scenes by cloud cover improves data quality, but may result in sparse coverage depending on the location and time of year.\n", + ">\n", + ">Modify the cloud cover slider below to visualize the number of items changing!" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "id": "a9d52e1b", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Number of items without filtering clouds: 360\n", + "Number of with filtering clouds: 46\n" + ] + } + ], + "source": [ + "items = get_items(aoi, collection_id, datetime_range)\n", + "# uncomment to see the items, output to large when default rendering not enabled such as on github\n", + "print(f\"Number of items without filtering clouds: {len(items)}\")\n", + "items = get_items(aoi, collection_id, datetime_range, max_cloud_cover=10)\n", + "print(f\"Number of with filtering clouds: {len(items)}\")" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "id": "afd780c5", + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "9d7506896af7434f88ac22eda23e8438", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "IntSlider(value=10, continuous_update=False, description='Cloud cover')" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "cloud_slider = widgets.IntSlider(\n", + " value=10,\n", + " min=0,\n", + " max=100,\n", + " step=1,\n", + " description='Cloud cover',\n", + " continuous_update=False,\n", + " orientation='horizontal',\n", + " readout=True,\n", + ")\n", + "\n", + "def on_slider_change(change):\n", + " clear_output(wait=True)\n", + " display(cloud_slider)\n", + " \n", + " max_cloud = change['new']\n", + " items = get_items(aoi, collection_id, datetime_range, max_cloud_cover=max_cloud)\n", + " print(f\"{len(items)} items found with ≤ {max_cloud}% cloud cover\")\n", + " \n", + " # Optional: Plot items per month again\n", + " dates = [item.datetime for item in items if item.datetime is not None]\n", + " df = pd.DataFrame({\"datetime\": dates})\n", + " df[\"month\"] = df[\"datetime\"].dt.to_period(\"M\")\n", + " monthly_counts = df.groupby(\"month\").size()\n", + " \n", + " monthly_counts.plot(kind=\"bar\", figsize=(10, 4))\n", + " plt.title(f\"Items per Month (Cloud ≤ {max_cloud}%)\")\n", + " plt.xlabel(\"Month\")\n", + " plt.ylabel(\"Item Count\")\n", + " plt.xticks(rotation=45)\n", + " plt.tight_layout()\n", + " plt.show()\n", + "\n", + "cloud_slider.observe(on_slider_change, names='value')\n", + "display(cloud_slider)" + ] + }, + { + "cell_type": "markdown", + "id": "c2421228", + "metadata": {}, + "source": [ + "> 🧭 **Note on Coordinate Reference Systems (CRS)**\n", + ">\n", + "> Your area of interest (AOI) is defined in **WGS 84 (EPSG:4326)**, a geographic coordinate system using latitude and longitude.\n", + "> Sentinel-2 imagery is typically stored in **UTM zones** (e.g., EPSG:32618), a projected CRS in meters.\n", + ">\n", + "> Coordinates in the item bbox and geometry fields of a STAC item are always stored in the EPSG:4326 coordinate system (lat/lon). The underlying files may use a different CRS but the STAC metadata will always be EPSG:4326. When filtering STAC items with the bbox or intersects parameters be sure to provide coordinates in lat/lon!\n", + ">\n", + ">The odc-stac.load() function then uses those matched STAC items to load the imagery, preserving the original raster CRS." + ] + }, + { + "cell_type": "markdown", + "id": "48645db9", + "metadata": {}, + "source": [ + "### 📁 Exploring the STAC Items\n", + "\n", + "The returned STAC items contain metadata and asset links for each matching satellite scene. \n", + "Each item includes information like the acquisition time (`datetime`), asset URLs for individual bands (e.g., `\"B03\"`, `\"B11\"`, `\"SCL\"`), cloud cover estimates (`\"eo:cloud_cover\"`), and geometry of the image footprint.\n", + "\n", + "In this tutorial, we specifically use:\n", + "- `\"B03\"` (green band) and `\"B11\"` (SWIR band) to compute NDSI,\n", + "- `\"SCL\"` (Scene Classification Layer) to mask clouds and shadows,\n", + "- and the item geometries to align data spatially with our area of interest.\n", + "\n", + "We can take a quick look at one of the returned items to inspect its properties and preview the image using the rendered_preview asset provided by the Planetary Computer:" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "id": "dd05d682", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Datetime: 2023-12-21 15:56:51.024000+00:00\n", + "Cloud cover: 8.673067\n", + "Available assets: ['AOT', 'B01', 'B02', 'B03', 'B04', 'B05', 'B06', 'B07', 'B08', 'B09', 'B11', 'B12', 'B8A', 'SCL', 'WVP', 'visual', 'preview', 'safe-manifest', 'granule-metadata', 'inspire-metadata', 'product-metadata', 'datastrip-metadata', 'tilejson', 'rendered_preview']\n", + "Geometry type: Polygon\n" + ] + }, + { + "data": { + "text/html": [ + "" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# Show basic info from the first STAC item\n", + "item = items[0]\n", + "print(\"Datetime:\", item.datetime)\n", + "print(\"Cloud cover:\", item.properties.get(\"eo:cloud_cover\"))\n", + "print(\"Available assets:\", list(item.assets.keys()))\n", + "print(\"Geometry type:\", item.geometry[\"type\"])\n", + "\n", + "display(Image(url=item.assets[\"rendered_preview\"].href))" + ] + }, + { + "cell_type": "markdown", + "id": "2d6dd9a2", + "metadata": {}, + "source": [ + "### 📦 Loading Data with ODC-STAC and Xarray\n", + "\n", + "We use `odc-stac` to load the Sentinel-2 imagery directly from the STAC items into an `xarray.Dataset`. \n", + "This library simplifies working with cloud-optimized geospatial data and integrates cleanly with the Planetary Computer.\n", + "\n", + "We configure band aliases (e.g., `\"green\"` for B03 and `\"swir\"` for B11) to make the code more readable. \n", + "The result is an `xarray.Dataset` with time, latitude, and longitude dimensions that supports efficient analysis and lazy loading with Dask." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "d4be725a", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
<xarray.Dataset> Size: 61GB\n",
+       "Dimensions:      (y: 10980, x: 20982, time: 22)\n",
+       "Coordinates:\n",
+       "  * y            (y) float64 88kB 5.1e+06 5.1e+06 5.1e+06 ... 4.99e+06 4.99e+06\n",
+       "  * x            (x) float64 168kB 5e+05 5e+05 5e+05 ... 7.098e+05 7.098e+05\n",
+       "    spatial_ref  int32 4B 32618\n",
+       "  * time         (time) datetime64[ns] 176B 2023-01-17T15:45:59.024000 ... 20...\n",
+       "Data variables:\n",
+       "    green        (time, y, x) float32 20GB dask.array<chunksize=(1, 10980, 20982), meta=np.ndarray>\n",
+       "    swir         (time, y, x) float32 20GB dask.array<chunksize=(1, 10980, 20982), meta=np.ndarray>\n",
+       "    SCL          (time, y, x) float32 20GB dask.array<chunksize=(1, 10980, 20982), meta=np.ndarray>
" + ], + "text/plain": [ + " Size: 61GB\n", + "Dimensions: (y: 10980, x: 20982, time: 22)\n", + "Coordinates:\n", + " * y (y) float64 88kB 5.1e+06 5.1e+06 5.1e+06 ... 4.99e+06 4.99e+06\n", + " * x (x) float64 168kB 5e+05 5e+05 5e+05 ... 7.098e+05 7.098e+05\n", + " spatial_ref int32 4B 32618\n", + " * time (time) datetime64[ns] 176B 2023-01-17T15:45:59.024000 ... 20...\n", + "Data variables:\n", + " green (time, y, x) float32 20GB dask.array\n", + " swir (time, y, x) float32 20GB dask.array\n", + " SCL (time, y, x) float32 20GB dask.array" + ] + }, + "execution_count": 18, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "cfg = {\n", + " \"sentinel-2-l2a\": {\n", + " \"aliases\": {\n", + " \"green\": \"B03\",\n", + " \"swir\": \"B11\"\n", + " }\n", + " }\n", + "}\n", + "\n", + "ds = load(\n", + " items,\n", + " bands=[\"green\", \"swir\", \"SCL\"],\n", + " stac_cfg=cfg,\n", + " chunks={},\n", + " dtype=\"uint16\",\n", + ")\n", + "ds" + ] + }, + { + "cell_type": "markdown", + "id": "30f9ae17", + "metadata": {}, + "source": [ + "### 🧾 Understanding the Xarray Output\n", + "\n", + "The output from `odc.stac.load()` is an `xarray.Dataset` containing multi-dimensional arrays (DataArrays) for each selected band. \n", + "The key dimensions are:\n", + "\n", + "- `time`: one entry per satellite scene that matched your query\n", + "- `y` and `x`: the spatial grid in projected coordinates\n", + "- `spatial_ref`: defines the coordinate reference system (CRS)\n", + "\n", + "Each band (e.g., `\"green\"`, `\"swir\"`, `\"SCL\"`) is stored as a Dask-backed array, enabling efficient, lazy computation on large datasets.\n", + "\n", + "This structure allows you to analyze changes over time and space while only loading data into memory when needed." + ] + }, + { + "cell_type": "markdown", + "id": "e3515fe3", + "metadata": {}, + "source": [ + "### ☁️ Masking Clouds and Shadows with the Scene Classification Layer (SCL)\n", + "\n", + "To improve the accuracy of snow detection, we use the Scene Classification Layer (SCL) to filter out cloudy or shadowed pixels. \n", + "Each pixel in the SCL band is assigned a class (e.g., cloud shadow, medium/high probability cloud, cirrus, snow).\n", + "\n", + "We define a list of classes to exclude, then create a `valid_mask` that flags only the clear pixels. \n", + "This mask is applied to the computed NDSI values so that we only analyze regions with reliable surface observations.\n", + "\n", + "The result is `ndsi_masked`, an `xarray.DataArray` that contains NDSI values for valid pixels and `NaN` elsewhere." + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "id": "fd2f4f3a", + "metadata": {}, + "outputs": [], + "source": [ + "# SCL values to mask out\n", + "cloud_codes = [3, 8, 9, 10] # shadow, medium/high clouds, cirrus\n", + "\n", + "# Valid where SCL is not in cloud codes\n", + "valid_mask = ~ds.SCL.isin(cloud_codes)" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "id": "be8def29", + "metadata": {}, + "outputs": [], + "source": [ + "ndsi = (ds.green - ds.swir) / (ds.green + ds.swir + 1e-6)\n", + "ndsi_masked = ndsi.where(valid_mask)" + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "id": "c1aece3f", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
<xarray.DataArray (time: 22, y: 10980, x: 20982)> Size: 20GB\n",
+       "dask.array<where, shape=(22, 10980, 20982), dtype=float32, chunksize=(1, 10980, 20982), chunktype=numpy.ndarray>\n",
+       "Coordinates:\n",
+       "  * y            (y) float64 88kB 5.1e+06 5.1e+06 5.1e+06 ... 4.99e+06 4.99e+06\n",
+       "  * x            (x) float64 168kB 5e+05 5e+05 5e+05 ... 7.098e+05 7.098e+05\n",
+       "    spatial_ref  int32 4B 32618\n",
+       "  * time         (time) datetime64[ns] 176B 2023-01-17T15:45:59.024000 ... 20...
" + ], + "text/plain": [ + " Size: 20GB\n", + "dask.array\n", + "Coordinates:\n", + " * y (y) float64 88kB 5.1e+06 5.1e+06 5.1e+06 ... 4.99e+06 4.99e+06\n", + " * x (x) float64 168kB 5e+05 5e+05 5e+05 ... 7.098e+05 7.098e+05\n", + " spatial_ref int32 4B 32618\n", + " * time (time) datetime64[ns] 176B 2023-01-17T15:45:59.024000 ... 20..." + ] + }, + "execution_count": 21, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ndsi_masked" + ] + }, + { + "cell_type": "markdown", + "id": "33286a3e", + "metadata": {}, + "source": [ + "### 🔍 Comparing Raw and Masked NDSI\n", + "\n", + "The unmasked NDSI array includes values for all pixels, regardless of cloud or shadow contamination. \n", + "In contrast, `ndsi_masked` applies a filter using the SCL layer to remove unreliable observations, replacing them with `NaN`.\n", + "\n", + "By comparing the two arrays:\n", + "- The shape and dimensions remain the same.\n", + "- The masked version is **sparser**, containing only valid surface pixels.\n", + "- The difference in data volume highlights how much of the imagery was affected by cloud or shadow cover.\n", + "\n", + "This step ensures that subsequent analyses or visualizations reflect only high-quality observations." + ] + }, + { + "cell_type": "markdown", + "id": "40e07e67", + "metadata": {}, + "source": [ + "### 📈 Generating an NDSI Time Series\n", + "\n", + "To understand how snow presence changes over time, we compute the mean NDSI for each scene by averaging across all valid (non-masked) pixels. \n", + "This produces a 1D time series showing the evolution of surface conditions in your selected area.\n", + "\n", + "Plotting this time series provides a clear view of when snow was present, and how it varied across the selected date range.\n", + "It also helps identify scenes that were too cloudy to yield valid NDSI values (which will appear as gaps or drops in the plot)." + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "id": "6ef79505", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/Users/zacdez/Documents/github/PlanetaryComputerExamples/.venv/lib/python3.11/site-packages/rasterio/warp.py:387: NotGeoreferencedWarning: Dataset has no geotransform, gcps, or rpcs. The identity matrix will be returned.\n", + " dest = _reproject(\n" + ] + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA/UAAAGJCAYAAADLxBNqAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAmJ5JREFUeJzt3Qd408UbB/C3u7TQMsoehbL3BtlDEBQQnKhsFRXEhQsQQeSvOBFUEAUZCiqKojIE2TJlD9mjbOgAuunO//le+YWkTdKmNM36fp4nTzMuySW5/Jr37r07D51OpxMiIiIiIiIicjqe9q4AEREREREREeUPg3oiIiIiIiIiJ8WgnoiIiIiIiMhJMagnIiIiIiIiclIM6omIiIiIiIicFIN6IiIiIiIiIifFoJ6IiIiIiIjISTGoJyIiIiIiInJSDOqJiIiIiIiInBSDeiIiIqICVrVqVRk6dKi9q0FERG6AQT0REbml+fPni4eHhzpt2bIlx+06nU4qV66sbu/du7c4egCJer7wwgs5btu4caO6bcmSJSZfO07+/v5SoUIF6dGjh3z++ecSHx9v8nnwPt17771SsWJFdZ8qVapInz595IcffjAqh8ccNWpUnuqelpamnrNly5ZSrFgxKVq0qDqP63Cbo9Dex7yciIiICpN3oT4bERGRg0FwiqC0ffv2Rtdv2rRJLl68KH5+fuIsZs+eLWPHjlUBel68++67Uq1aNRU8X716VQWuL7/8skydOlX+/PNPadSokb7sL7/8Iv3795cmTZrISy+9JCVKlJDw8HD5559/1PM+8cQTVtc3MTFRevXqpd5rdJxgZNvT01NWrVqlnuO3336TFStWSGBgoNhb3bp15fvvvze6Du81OiHeeuutHOWPHz+uXgsREZGtMagnIiK3dt9996mAFSPD3t63/y0i0G/evLlER0eLM6hfv74KJD/44AP1WvICo+4tWrQwClLXr1+vAuz7779fjh49KkWKFFG3vfPOO1KvXj3ZsWOH+Pr6Gj1OZGRkvuo8evRoFdB/8cUXRiP7I0aMkBkzZqjrXnvtNfnqq6+ksCBDIzk5Wf+6NWXLlpWBAwcaXYf3OiQkJMf14EydQURE5NzYhUxERG7t8ccfl2vXrsmaNWv016Wmpqp0dXOjz5mZmTJt2jQVSGOkHwHfs88+Kzdu3DAq98cff6iRaIycI8irXr26TJ48WTIyMozKde7cWRo0aCBHjhyRLl26SEBAgEpx/+ijj6xKwR88eLAaNb98+bLkV9euXeXtt9+Wc+fOycKFC/XXnz59WqXFZw/ooUyZMlY/D7Igvv32W/V8plL1n3/+efVezJkzR5UFvEe4ztTngffr4YcftvozwvuGTozVq1erDg4E819//bUU9Jx6bcoDpjC8+OKLUrp0aSlevLiqE9pbTEyM+vyQAYHTG2+8oToYsr/OvLwmIiJyLwzqiYjIrSH4atOmjfz444/66/766y+JjY2Vxx57zOR9EEi9/vrr0q5dO5k+fboMGzZMFi1apOakG84DRyCH9GyMSKMcRv4nTJggY8aMyfGYCMx69uwpjRs3lk8//VTq1Kkjb775pqpLXiENPD09XY0g34lBgwapv3///bf+utDQUFm3bp0+wL5TeF3o3EAgaw5uw+tBOj4g/R/p/pgqYAiBMjoyDD+vvH5GgAwHdO50795dlcUUA1vBugcnT56USZMmqWyIb775RnWiYG0CvB/vv/++mgry8ccf50j3t+Y1ERGRG9ERERG5oXnz5mEYVLdr1y7dl19+qStWrJguKSlJ3fbII4/ounTpos6HhobqevXqpb/f5s2b1f0WLVpk9HirVq3Kcb32eIaeffZZXUBAgC45OVl/XadOndR9v/vuO/11KSkpunLlyukeeuihXF+LYR2HDRum8/f3112+fFld3rBhg3rsX375xeRrNyc4OFjXtGlT/eVvv/1W3cfX11e9N2+//bZ6LzIyMnLcF+Wef/55i3V++eWXVbl9+/aZLbN3715VZvTo0ery8ePH1eUvvvjCqNzIkSN1RYsW1b/f1nxGeO9wHW6zVv369dVnZwoed8iQITne8x49eugyMzP117dp00bn4eGhe+655/TXpaen6ypVqmT02Na8JiIici8cqSciIrf36KOPys2bN2X58uVq5Xf8NZd6j/n3wcHBalQX8+21E0bhMSq/YcMGfVnDedl4XJTr0KGDJCUlybFjx4weF/c1nJuNNPdWrVrJmTNnrHot48ePL5DRetTHcBX8J598Uo2YY6oARsYxjQCvpWbNmrJt2zarH197bKx4b452W1xcnPpbq1YtNYq+ePFifRmMbmOqBEa6tffbms8IsFggRrsLw1NPPWW0Qn7r1q1Vmj2u13h5eampAIafvbWviYiI3AcXyiMiIreH+c3dunVTi+Mh4EagaDg/2xBSp5Gab24eueGicYcPH1ZBNhaf0wJTDR7DUKVKlXJsh4a51QcPHrTqtYSFhan0eaR1m0rzz6uEhIQcrxGBL054j/bs2aOC61mzZqk56eiksGZuvRawm9s+z1zgjxT8cePGyaVLl9Q8eqzYj/cc1+fnM9KC+sKCbQANIVAHbJ+Y/XrDufLWviYiInIfDOqJiIhE1Mj88OHD1XxtrAqPRcxMwWJlCKwwl9lcBwFg4bNOnTpJUFCQ2joOi+RhcbO9e/equfJ4HEMYnTUl+2JpeZ1bj/nYH374ofTr18/q+2PePALIGjVqmLwdC/lhlB4nrP6O+eGYIz9kyBCrtogDdFqYm8OudWhg1X0Ngnes0o+Ra2y/9/PPP6sAGOsRWPsZabKvdG9L5j5nU9cbfvbWviYiInIfDOqJiIhE5IEHHlALkWHLNsP07uwQnK9du1YtVmYpGMQIMlbVx17rHTt21F+Pvd1tDXVEKj9WcUd6t7W0BdrykpKubYl35coVq54DHScIZPFc5hbL++6779Q2g4YBO0bVMS0BnxFWzcf7i44Lwy3k8voZORNXfE1ERFQwOKeeiIjo1hxy7IeO/dgxP9vS/Huk52NOeXaYy44ResORV8PRVmxdNnPmTCkMSPvHiujWbIsHmCqA14bgecCAAfrrsfK9KStXrlR/a9eubdXzIN0cq7cjUDW1Dz3S+lEXzDXH1ARDGK1H58vcuXPVvHLD1HtrPiNn4oqviYiICgZH6omIiG7JS/o4Uuoxoj9lyhTZv3+/3HPPPeLj46PmPCMlHFuNYT5+27Zt1Zx4PCb2Jcd8eYxK5yed/k5G6xcsWGC2DFLmMRceQWFERIQKotesWaO2r/vzzz/VdAFN3759VaCPDg88dmJiogrIly1bpvavt9QRYs5nn32mnn/kyJFqET5tRB57xv/xxx/qvcb2fqYC3Ndee02dSpYsqdZDyM9n5Exc8TUREVHBYFBPRERkJYwiY9VxpLdj0TakiGO/ewTRSI+GUqVKqVX0X331VTVqjgAft999992FttI6nnfhwoVqhNeUCRMm6FfaR3DcsGFDmTZtmhpBz74q/Zw5c1SgjTns2BMenRNYlA/z97FGAN6D/GRHIAMA2QuoJ/Zgx+PWqVNH1QPBPgLX7DByj06TrVu3ytNPP22yTF4+I2fjiq+JiIjunAf2tSuAxyEiIiIiIiKiQsY59UREREREREROikE9ERERERERkZNiUE9ERERERETkpBjUExERERERETkpBvVERERERERETopBPREREREREZGT4j71ucjMzFT78WK/Xg8PD3tXh4iIiIiIiFycTqeT+Ph4qVChgnh6Wh6LZ1CfCwT0lStXtnc1iIiIiIiIyM1cuHBBKlWqZLEMg/pcYIQe5syZI/369RMfHx97V4mcWFpamvz9999yzz33sC2RTbCNkS2xfVFhYVsjW2MbI0dvX3FxcWpwWYtHLWFQnwst5T4gIECCgoL4pac7/oKzLZEtsY2RLbF9UWFhWyNbYxsjZ2lfeZkCzoXyiIiIiIiIiJyU0wX1M2bMkKpVq4q/v7+0bt1adu7cabbs/PnzVc+G4Qn3IyIiIiIiInIFThXUL168WEaPHi0TJ06UvXv3SuPGjaVHjx4SGRlp9j5Iebhy5Yr+dO7cuUKtMxEREREREZGtOFVQP3XqVBk+fLgMGzZM6tWrJ7NmzVJzFebOnWv2PhidL1eunP5UtmzZQq0zERERERERka04zUJ5qampsmfPHhk7dqz+OuzX161bN9m+fbvZ+yUkJEhoaKjab75Zs2by/vvvS/369c2WT0lJUSfDVQcNFzwguhNaG2JbIlthGyNbYvuiwsK2RrbGNkaO3r6sua/TBPXR0dGSkZGRY6Qdl48dO2byPrVr11aj+I0aNZLY2Fj55JNPpG3btnL48GGze/1NmTJFJk2aZPK2NWvWFMArIWJbIttjGyNbYvuiwsK2RrbGNkaO2r6SkpJcL6jPjzZt2qiTBgF93bp15euvv5bJkyebvA8yATBvP/v+gNC9e3eH3PIiI1Mnu8/dkMj4FClTzE9ahJYQL8/ctz6gwoceN3y5HbUtkfNjGyNbYvuiwsK2RrbGNkaO3r4MM8ZdJqgPCQkRLy8viYiIMLoelzFXPi/whjZt2lROnTpltoyfn586mbu/o33pV/13RSYtOyJXYpP115UP9peJfepJzwbl7Vo3Ms8R2xK5FrYxsiW2LyosbGtka2xj5Kjty5r7Oc1Ceb6+vtK8eXNZt26d/jrMk8dlw9F4S5C+f+jQISlf3jWCXQT0IxbuNQro4WpssroetxMREREREZHrcpqgHpAWP3v2bFmwYIEcPXpURowYIYmJiWo1fBg8eLDRQnrvvvuu/P3333LmzBm1Bd7AgQPVlnZPP/20ODuk3GOEXmfiNu063I5yRERERERE5JqcJv0e+vfvL1FRUTJhwgS5evWqNGnSRFatWqVfPO/8+fNqRXzNjRs31BZ4KFuiRAk10r9t2za1HZ6z2xl+PccIvSGE8rj9zSUHpXOd0lKjTFGpWipQ/H28CrWeREREREREZDtOFdTDqFGj1MmUjRs3Gl3+7LPP1MkVRcabD+gNLdl7UZ0Aa+dVLhkg1UsXVUF+9dKB+vPFA3xtXGMiIiIiIiISdw/qKUuZYv55KtepVojEJafLqcgEiU9Ol3PXktRp/bFIo3KlAn2lugr0jQP+isWLiCdX0iciIiIiInJIDOqdVKtqJdUq91gUz9SseYTh5YL9Ze7QVmp7O51OJ1EJKXI6MlFORSXI6cgEOX3r7+XYZLmWmCrXwq+rtH5D/j6eEhaiBfq3/pYJZCo/ERERERGRA2BQ76QQqGPbOqxyjwDeMLDXxtVxu7ZfvYeHhxrdx6lN9VJGj5WYki5nohJVkI8Rfe3v2WuJkpyWKUeuxKmTIS2Vv0ZpBPlZI/ta4M9UfiIiIiIiosLBoN6JYR/6rwY2y7FPfTkr96kP9POWhpWC1clQekamnL+eJKdNBPyGqfzrsqXyhxT1lTB9Gv/tdP4KwUzlJyIiIiIiKkgM6p0cAvfu9cqptHksnoeReKTmayP0d8Lby1MF5zh1l6wdBkBL5c8K8hNzpPJHJ6RKdELOVP4iPl4SZrA4n/a3akiA+HkzlZ+IiIiIiMhaDOpdAAL47Cn1tmSYyt+2eojRbQkp6RIehXn78Vnz928F/Ejlv5mWIYcvx6mTIfQ/VLm1Kj9S+bNS+gOlRuliEhzgk+96ZmTqbNLZQURERERE5CgY1FOBKmomlT8tI1MuXE+6PbqvpfMjlT8lXc5eS1InU6n8xsF+1uh++SB/i6n8q/67kmNaQnkrpyUQERERERE5Ogb1VCh8DFL5DalU/vgUgxX5b4/uXzFI5f/XTCq/8bz9rFT+Dcci1QKC2XcFwE4BuB7rEDCwJyIiIiIiV8CgnuxKpfIH+auTqVT+Mwj29aP6WdvxnY02n8qPsXuM4Jva5k9363aM4GMdAqbiExERERGRs2NQTw6dyt+oUnF1yp7Kr1blj0y4NcKfFeyfuZXKj7n05uAWZABgrn1hrkNARERERERkCwzqySlT+dU8+9JF5Z5sqfwLd5yTt/84nOtjYPE8IiIiIiIiZ+dp7woQFWQqf40yxfJUFqvhExEREREROTsG9eRSsG0dVrm3NFu+bJCfKkdEREREROTsGNSTS8Hid9i2DswF9lg5PzU9s1DrRUREREREZAsM6snlYLs6bFtXLtg4xb50MT8J8PWSs9eS5OXF+ywuqEdEREREROQMuFAeuWxgj23rsMo9FsXDHHqk3O85d0MGzvlXVh+OkCkrj8r43lmj+kRERERERM6II/Xk0qn42Laub5OK6i8uI7D/+JFG6vY5W8Jlwbaz9q4mERERERFRvjGoJ7eDIP/1HrXV+UnLDsvaIxH2rhIREREREVG+MKgntzSyc3V5rGVlwbT6F37cJ4cuxtq7SkRERERERFZjUE9uu6f95H4NpEPNELmZliFPLtglF28k2btaREREREREVmFQT27Lx8tTZg5oJnXKFZOo+BR5cv4uib2ZZu9qERERERER5RmDenJrxfx9ZO7QllI2yE9ORCTIyEV7uIc9ERERERE5DQb15PYqFC+iAvtAXy/ZeuqavLX0kOh03MOeiIiIiIgcH4N6IhGpXyFYvhzQTG1798uei/Ll+lP2rhIREREREZHrBfUzZsyQqlWrir+/v7Ru3Vp27tyZp/v99NNPanG0fv362byO5Jy61C4j7/atr85/uuaELN130d5VIiIiIiIicp2gfvHixTJ69GiZOHGi7N27Vxo3biw9evSQyMhIi/c7e/asvPbaa9KhQ4dCqys5pwGtQ+XZTmHq/BtLDsqOM9fsXSUiIiIiIiKzvMWJTJ06VYYPHy7Dhg1Tl2fNmiUrVqyQuXPnypgxY0zeJyMjQwYMGCCTJk2SzZs3S0xMTL6eOznZSxITRXx8jK/38hLx9799GWXM8fQUKVIkf2WTkkTMTfP28BAJCMhf2Zs3RTItrAsXGJi/ssnJeO8Lpizqi3pDSopIenrBlMX7i/cZUlNF0m4tfD+qfR0Jv5wiqw5flafn7JMfhreWBqHFTJY1Be0B7cJUWZw3bEuGZXEbypvj5yfi7W19WbwHeC/M8fW93a6tKYvPDJ+dOSiH8taWRRtDWyuIsngP8F4AvhP4bhREWWu+94V9jMjexoDHiII9RtxpWUvHCEtlHeUYYe7/YfayPEY45jHCFEc9RlhqazxGOO4xwtl+R5hrY8BjhGMfIxz5d0RamvFryc8xwtLnl4POSaSkpOi8vLx0S5cuNbp+8ODBuvvvv9/s/SZMmKDr16+fOj9kyBBd3759LT5PcnKyLjY2Vn+6cOECmqxOJFaX1XyNT/fem6FLTU3VnwICMk2Ww6ljR+OyISHmyzZvblw2NNR82bp1M43K4rK5sngcw7J4HnNlUT/Dsqi/ubJ43YZl8b6YK4uTYdkHH7Rc9saN22UHDbJc9tKl22Wfey7dYtkTJ26XHT3actkNWxP1ZcePt1x227Y0fdkpUyyXXbPmdtnp0y2X/f3322XnzEmzWPaHH26XxXlLZfFYWlk8h6WyqKNWFnW3VBavXSuL98RSWbynWtl9+1ItlsVnpZXFZ2ipLNqAVhZtw1JZtC2tLNqcpbJos4Zt2FJZHiPc4xiBdusOx4hff71psSyPEVknHiPu/BjRo4fl9s5jhGMeI5zpd0RiYqLu669XWyzLY4TjHiOc4XfE11+vVu0s/8cIxJ+iYtLcOM1IfXR0tBp1L1u2rNH1uHzs2DGT99myZYt8++23sn///jw/z5QpU9Sofl4h9X/lyn/1lzMyeplNgLh27bqsXLlVfzk1tSf6Qk2WjY2NlZUr/9FfTkrqjr4hk2UTEuJl5coNBpe7iEiQybJJSTdl5co1Bs/TUURKmCybmpoqK1euMqh/OxEJMVkWn83KlSv1lyMjW4tIOTHHsOzVqy1EpKLZsqtXrxZ//6zutYsXm4pIFbNl165dK8HBWV3P5841EpFqZstu2LBBypbN6so9c6aeiNQ0W/a1H3fKm1djxddL5OTJ2iJSx2zZrVu3SmRkVkbIsWM1sAyf2bI7duyQxMSsFP/Dh1FX1Nm03bt3i0iEOn/gQGURaWa27L59+yQg4PKt8xVEpKXZsgcOHJSVKy/ceg58v+4yW/bw4cOycmW4On/oUCkRaW+2LL6XK1dmLTh48mRxEelktuzJkydl5crj6vz588VEpKvZsmfOnJGVK4+o8xER6GK+x2zZc+fOy8qVB9X52Fh0y99rtuzFixdl5cp96jxGh0R6my179eoVWbkSn4emr9myPEa4xzECmWDnzsWr8658jNi3j8cI4DHC9seI6GgeI5zxGOFsvyNEDIaqTeAxwnGPEc7wOwLWrFlT4McIUzwQ2YsTuHz5slSsWFG2bdsmbdq00V//xhtvyKZNm+Tff283doiPj5dGjRrJzJkz5d57s76AQ4cOVen3v//+u9nnSUlJUSdNXFycVK5cWebOXSy9e/cWn2z5OUyJcc6UmLymwp29liiD5u6RuPRkuadeGfniscaSke5xB+n3abJ+/Xrp2rWraktMm8tZlqm1d5p+b9zGgMeI/JVlam3O7/3Nm2ny11/G7ctcWR4jHPMYYYojHiPi49NkzRrzbY3HCMc8RjjT7wj8v1y9eo20b9/dZBsDHiMc9xjh6L8j0tLSZMuWNdKjR1b7ys8xAnFoaGiI6oAJCjLdiaJxmpH6kJAQ8fLykoiIrB5GDS6XK5ezl+b06dNqgbw+ffror8u81Uq8vb3l+PHjUr169Rz38/PzU6fs0HtTvLiP2S+9pjg6EvPImrLBwbYpm8vLcfuyTYoXl7nDm8qA2f/KmqOR8vGaU/J273r5flx8Qc21JVw0PAjm9rjWlDU8aBdkWcN/MgVVFkx8BQukrPYPv6DL2up7n59jhKU2lr2ss38/3bGsvY8R1vw/5DHCMY8RBV3WVu29WDHr2pot6uCMZe19jHC23xEIovLSxjQ8Rlhf1lG+G4VdFr/H0L7QtrJO1j+up6eP661+7+vrK82bN5d169YZBem4bDhyr6lTp44cOnRIpd5rp/vvv1+6dOmizmP0nSgvWlYtKZ882lid/3ZLuMzfmpU2RkREREREZG9OM1IP2M5uyJAh0qJFC2nVqpVMmzZNEhMT9avhDx48WKXoY1489rFv0KCB0f2L3+quyn49UW7ub1xBLt5Iko9WHZd3lx+RiiUCpHs94/UdiIiIiIiICptTBfX9+/eXqKgomTBhgly9elWaNGkiq1at0i+ed/78efHUJiMQFbARnarLhetJ8uPOC/Lij/tk8bN3SaNKVuQ1ERERERERuXNQD6NGjVInUzZu3GjxvvPnz7dRrcgdeHh4yLt9G8ilmGT550SUPDl/t/z+fFupVCKPk9KIiIiIiIgKGIe1iazg4+UpM55oKnXKFZPohBQZNm+XxN60sJQlERERERGRDTGoJ7JSMX8fmTespZQN8pOTkQkyYuEeSU23sP8GERERERGRjTCoJ8qH8sFFZO7QlhLo6yXbTl+Tsb8dEp25DTuJiIiIiIhshEE9UT7VrxAsXw5oJl6eHvLr3ovy+bpT9q4SERERERG5GQb1RHegS+0yMrlv1haJn609IUv3XbR3lYiIiIiIyI0wqCe6Q0+0riLPdaquzr+x5KBsP33N3lUiIiIiIiI3waCeqAC80aO29GpUXtIydPLs97vlVGS8vatERERERERugEE9UQHw9PSQTx9pLM1DS0hccroMnbdLouJT7F0tIiIiIiJycQzqiQqIv4+XzB7cQqqWCpCLN27K0wt2yc3UDHtXi4iIiIiIXBiDeqICVDLQV+YNayUlAnzkwMVYeemnfZKRya3uiIiIiIjINhjUExWwaiGBasTe19tT/j4SIe+tOGrvKhERERERkYtiUE9kAy2qllRz7GHu1nCZtzXc3lUiIiIiIiIXxKCeyEb6NK4gb/aso86/u/yI/H34qr2rRERERERELoZBPZENPdcpTB5vVUV0OpEXf9onBy/G2rtKRERERETkQhjUE9mQh4eHTO5bXzrVKi3JaZny7KJ9ci3Z3rUiIiIiIiJXwaCeyMa8vTxlxoBmUrd8kEQnpMrXx7wk7maavatFREREREQugEE9USEo6uctc4e2kLJBfhJx00Oe/3G/pKZn2rtaRERERETk5BjUExWS8sFFZPbAZuLnqZMd4TdkzG8HRYfJ9kRERERERPnEoJ6oENUtX0yG1coUL08P+W3vJZm+7qS9q0RERERERE6MQT1RIatbQifv9K6rzk9be1KW7Llo7yoREREREZGTYlBPZAePtawkIzpXV+fH/nZQtp2OtneViIiIiIjICTGoJ7KT1++pLb0blZe0DJ08+/0eORkRb+8qERERERGRk2FQT2Qnnp4e8skjjaVFaAmJT06XofN2SWQ8N7EnIiIiIqK8Y1BPZEf+Pl7yzeAWUi0kUC7F3JSnF+yWpNR0e1eLiIiIiIicBIN6IjsrGegr84a2lBIBPnLwYqy89NN+ycjkVndEREREROSCQf2MGTOkatWq4u/vL61bt5adO3eaLfvbb79JixYtpHjx4hIYGChNmjSR77//vlDrS5QXVUMCZc6QFuLr7SlrjkTI/1YcsXeViIiIiIjICThVUL948WIZPXq0TJw4Ufbu3SuNGzeWHj16SGRkpMnyJUuWlLfeeku2b98uBw8elGHDhqnT6tWrC73uRLlpHlpSPnu0iTo/b+tZmbsl3N5VIiIiIiIiB+ctTmTq1KkyfPhwFZjDrFmzZMWKFTJ37lwZM2ZMjvKdO3c2uvzSSy/JggULZMuWLaozwJSUlBR10sTFxenPp6WlFeCrIXektSFzbemeuiHyRo+a8tHqkzJ5xREpH+Qr3eqWKeRakiu3MaI7wfZFhYVtjWyNbYwcvX1Zc18PnU7nFJN3U1NTJSAgQJYsWSL9+vXTXz9kyBCJiYmRP/74w+L98TLXr18v999/v/z+++/SvXt3k+XeeecdmTRpUo7rf/jhB/X8RLaGb+TP4Z6yLcJTfDx18kL9DAktau9aERERERFRYUlKSpInnnhCYmNjJSgoyDVG6qOjoyUjI0PKli1rdD0uHzt2zOz98CZUrFhRjb57eXnJzJkzzQb0MHbsWJXibzhSX7lyZXUe9/Px8SmQ10PuCT1ua9asybUt9cjIlGcX7ZN/Tl6TBWcC5JdnW0nlEuxUooJrY0T5wfZFhYVtjWyNbYwcvX0ZZoznxmmC+vwqVqyY7N+/XxISEmTdunUqYA8LC8uRmq/x8/NTJ1PwgfBLTwUht7aEm2YObCGPzNouR6/EyTML98uvz7WV4AC2P8obHq/Ilti+qLCwrZGtsY2Ro7Yva+7nNAvlhYSEqJH2iIgIo+txuVy5cmbv5+npKTVq1FAr37/66qvy8MMPy5QpUwqhxkR3pqift9rqrlyQv5yKTJBnF+6W1PRMe1eLiIiIiIgciNME9b6+vtK8eXM12q7JzMxUl9u0aZPnx8F9DBfCI3Jk5YL9Zd6wlirA33Hmuoz59aBaH4KIiIiIiMipgnpA6vzs2bPVCvZHjx6VESNGSGJion41/MGDB6s58RqMyGMuw5kzZ1T5Tz/9VO1TP3DgQDu+CiLr1C0fJDMGNBMvTw/5bd8lmbb2pL2rREREREREDsKp5tT3799foqKiZMKECXL16lWVUr9q1Sr94nnnz59X6fYaBPwjR46UixcvSpEiRaROnTqycOFC9ThEzqRTrdLyv34NZOxvh2T6upNSuWSAPNy8kr2rRUREREREduZUQT2MGjVKnUzZuHGj0eX//e9/6kTkCh5vVUUuXE+SmRtPqzT88sH+0q5GiL2rRUREREREduRU6fdE7u61e2pLn8YVJD1TJ88t3CMnIuLtXSUiIiIiIrIjBvVETsTT00M+friRtKxaQuKT02XYvF0SGZ9s72oREREREZGdMKgncjL+Pl7yzaAWUi0kUC7F3JSn5u+WpNR0e1eLiIiIiIjsgEE9kRMqEeir9rAvGegrhy7Fyos/7peMTG51R0RERETkbhjUEzmpqiGBMntwC/H19pS1RyNk8vIj9q4SEREREREVMgb1RE6seWgJmda/iTo/f9tZmbsl3N5VIiIiIiKiQsSgnsjJ3dewvIy7r446P3nFEVl9+Kq9q0RERERERIWEQT2RCxjeIUwG3lVFdDqRl37aJ/svxNi7SkREREREVAgY1BO5AA8PD3mnT33pUru0JKdlytMLdsmF60n2rhYREREREdkYg3oiF+Ht5SlfPtFM6lcIkuiEVBk6b6fEJqXZu1pERERERGRDDOqJXEign7fMHdpSygf7y+moRHl24W5JSc+wd7WIiIiIiMhGGNQTuZiyQf4qsC/q5y07zlyXMb8eEh0m2xMRERERkcthUE/kguqWD5KZA5qJl6eHLN13ST5be9LeVSIiIiIiIhtgUE/kojrWKi3v9Wugzn++7qT8svuCvatEREREREQFjEE9kQt7rFUVeb5LdXV+7G+HZOupaHtXiYiIiIiIChCDeiIX92r32nJ/4wqSnqmT577fI8evxtu7SkREREREVEAY1BO5OE9PD/n4kUbSqmpJiU9Jlyfn75LIuGR7V4uIiIiIiAoAg3oiN+Dn7SVfD2ouYSGBcinmpjy5YJckpqTbu1pERERERHSHGNQTuYkSgb4yb1hLKRnoK/9dipMXf9wnGZnc6o6IiIiIyJkxqCdyI6GlAmX24Bbi5+0p645FyrvLDnMPeyIiIiIiJ8agnsjNNA8tIdP6NxEPD5EF28/J3K1n7V0lIiIiIiLKJwb1RG7o3oblZdy9ddX5/604Iqv+u2rvKhERERERUT4wqCdyU093qCaD7goVZN+/vHif7Dt/w95VIiIiIiKymYxMnWw/fU3+2H9J/XWV9aW87V0BIrIPDw8PmdinnloNf/2xSHl6wW5ZOrKdVCkVYO+qEREREREVqFX/XZFJy47IldjbWzuXD/ZXv4d7NigvzszpRupnzJghVatWFX9/f2ndurXs3LnTbNnZs2dLhw4dpESJEurUrVs3i+WJ3I23l6d88XhTqV8hSK4lpsrQ+TslJinV3tUiIiIiIirQgH7Ewr1GAT1cjU1W1+N2Z+ZUQf3ixYtl9OjRMnHiRNm7d680btxYevToIZGRkSbLb9y4UR5//HHZsGGDbN++XSpXriz33HOPXLp0qdDrTuSoAv28Ze7QllIh2F/ORCXKs9/vkZT0DHtXi4iIiIjojmVk6tQIvalEe+063O7MqfhOlX4/depUGT58uAwbNkxdnjVrlqxYsULmzp0rY8aMyVF+0aJFRpfnzJkjv/76q6xbt04GDx5s1XNnZGRIampqju2/PD09xdv79tuIMpbSnX18fPJVNi0tzezWY7YqC76+vvkqm56eLpmZmQVSFvVFvW1ZFp8vTgVRFu0B7cJUWbyHhm3JUllrHtdSWbwHeC/M8fLykrJB/jJ3WEt59Kttsic8St78ea98+FAj/Ws2LItTXh83P2XxvuB9Koiyht9PW5XN7btc2MeI7G3MUllrHtccHiOsL1sY3/uCPkYYfpfN/T80VZbHCMc7RhRm2Ts9RlhqazxGFF5ZW32XHeUYYa6NAY8Rti1bWL8jtp+KlKjYRLOBb7p4qhH8neHXpWVocIEcI7K/lvwcIyx9fk4b1ONF7dmzR8aOHWvUyJFSj1H4vEhKSlJvcMmSJc2WSUlJUSdNXFyc+nv48GE5ffp0jvLVq1eX/v376y9/8sknZr/EVapUkYEDB+ovT5s2TW7evGmybPny5fWdF9q0g9jYWJNlQ0JC5JlnntFf/uabbyQ6Otpk2eDgYHn++ef1l+fNmydXrphONylSpIi88sor+ssLFy6U8+fPm22Ir7/+ulFWhan3SzNu3Dj9+d9++02OHTtmtuxrr72m/2IuW7ZMDh06ZLbsSy+9JIGBger8qlWrVEaHOSNHjpTixYur8+jo+ffff82WRWdS6dKl1fl//vlHtmzZYrbs0KFDpUKFCur8jh07ZP369TnKaK9hwIABEhoaqs7v3r1b/v77b7OP++ijj0qNGjXU+YMHD8ry5cvNln3ggQekbt2s1e2PHj0qS5cuNVu2d+/e0qhRI6leqohM7FRC/tu8SuTkPvngg5U5yiLTpUWLFur8uXPncnScGeratavcdddd6vzly5dl/vz5Zsu2b99eOnbsqM5HRUWpqTPmYNrN3Xffrc7HxMTIzJkzzZZt1qyZ9OzZU51PTEyU6dOnmy3bsGFD6dOnj/54g++yOXXq1JEHH3xQf3nKlClmy9rrGGH4PeExwvmOERpHOkbAiRMn1Ptr7j3mMcJ5jhGGHPEYsWTJEjlz5ozZtsZjhGMeI06dOiU///yz2bKOdIzA9yy37zKPEa7xO2LXxlUyqIj5LZy/v9lU0sVLrsQkyp8HNhXYMQLfIe2zys8xIjnZeKqASwT1aDjosShbtqzR9bhs6UBu6M0331QHSXQEmIMv1aRJk/JcLxw4Vq68HfxY6oG5du2aUVlLvS/4IWJYFh0S5iQkJBiVxWVz8DiGZfE85qB+hmVRf3Pwug3L4n2xxLCsuS+6ZvXq1foe1YsXL1osiy+M1puZW1lMy/Dz81Pnc5uSgS8XDjx5qe+2bdskICBrsTlzU0M0+GeNDqO8vGe7du1SP6hz+yxg3759Eh4enutnDAcOHNC/V4lmDuYa1FV7TfHx8RbL4nt5/fr1XNsvnDx5Ut9uzf3z0eBHntbxZtgBZwr+MWhtzVIPP+A90Mpa+h5rbcCwDVvCY0QWHiNc4xhh7gefhseILDxG3PkxwlzAoOExIguPEXd2jMgNjxGu8TsiMdZyG9acObxffGMK7hgBa9asKfBjhCkeOks5Dg4EPXQVK1ZUB7o2bdror3/jjTdk06ZNFns+4IMPPpCPPvpIzbPXehPzOlKPufjff/+96ok0TCcBpsSYLsvUWvPp9+hxR+8zHtNR0+a+2HBa5mw5K96eHjLjiSZyV7WSDps2x9TanOn3hm3MUllrHtccHiPcK7UW/x8xEmjYvsyV5THCMY8RhVn2To4RCMzWrl1rtq3xGFF4ZR0hpd4WxwiUwfGsc+fOJtsY8Bhh27KF9TsiJTVNun+2WaISTHfmZIinlAv2lw2jO4ouM6PA0u8RdyI7BffJzzECcSgGpNFZFhQUJC4xUo+0D3wBIyIijK7H5XLlylm8L9JUENTjn4OlgB7Q46r1uhrCcyPVwtyXXpPb7Szr3mXxBTfXlmxVBzDVpi2VfbNXQ7mSkCF/7L8sL/9yRJaMaCu1yxW748fNK8MDsTOUdaS2ZqmNFVYdWNa5ylr7/czr/0NrH5fHCJbNjr+9HKOsLb/L9j5GIIjK6/HMmse1tqwjfM6uXrZepZKy7ljOzBdt9aiJfeqLv1/eP7Pc6oDfY1oHh3ay9nFzy/hwytXv8cVo3ry5SnnQoGcElw1H7rPD6PzkyZPVnAdtDg8RWYaD0EcPN5JWVUtKfEq6DJu3UyLi8j6vh4iIiIjIEZyJSpBNJ7LS9UsEGAfuGKH/amAzp9+n3mlG6gHb2Q0ZMkQF561atVKLP2DRCm2RB6xojxR9bbGJDz/8UCZMmCA//PCD2tv+6tWr6vqiRYuqExGZ5+ftJd8Mbi4PfrVNbXX31IJdsviZNmoLPCIiIiIiZ/D+ymOSnqmTLrVLy5whLdUq95HxyVKmmL+0qlZSvDyNd3tyRk716xwrP2JRBATqCNCbNGmiRuC1xfOwmIU2/we++uorNZfk4YcfNnoc7HP/zjvvFHr9iZxN8QBfmT+0lTwwc6v8dylOXvhxn3wzqLl4ezlNkg8RERERuamtp6Jl7dEIFbi/1auu+tumeilxNU4V1MOoUaPUyRQsRmDo7FnzWxcQUd5UKRUgs4e0kMe/2SHrj0XKpGVH5N2+9XPsYU9ERERE5CgyMnUyefkRdX7QXaFSo0zO9aFcBYfbiChXzaqUkGn9mwji+O93nJNvt2Rtc0NERERE5IgW77ogx67GS3ARH3np7priyhjUE1Ge3NuwvLx1X111/r2VR+WvQ9bvoUlEREREZGtxyWny6d/H1XkE9CUCrVvZ3tkwqCeiPHuqfTUZ3CZUsIXoy4v3y77zN+xdJSIiIiIiIzPWn5JriakSVjpQBrUJFVfHoJ6I8gzz6Cf0rid31ykjKemZ8vSC3XL+WpK9q0VEREREpJy7lijztmatrTa+V13xcYMFnl3/FRJRgcLK958/3lQaVAxSPaBD5++UmKRUe1eLiIiIiEimrDwmqRmZ0qFmiHSpXUbcAYN6IrIa9qqfO6SlVAj2V3vYP/P9HklJz7B3tYiIiIjIjW0/fU1WHb4q2Hp+fK96brNbE4N6IsqXMkH+Mm9YKynm5y07w6/LG0sOig6T7YmIiIiI7LCF3f9WZG1h90TrKlK7nOtuYZcdg3oiyjccLL8a2Fy8PT3kj/2XZeqaE/auEhERERG5oV/3XJTDl+OkmL+3vNKtlrgTBvVEdEfa1wyR9x9sqM5/sf6U/Lzrgr2rRERERERuJCElXT5anbWF3Ytda0qpon7iThjUE9Ede7RFZXmxaw11ftzSQ7L5ZJS9q0REREREbuKrjackOiFFqpYKkCFtq4q7YVBPRAXile615IGmFSU9UycjF+6VY1fj7F0lIiIiInJxF64nyezN4er82Pvqiq+3+4W47veKicgmsLroBw81lNbVSkp8Sro8OW+XRMQl27taREREROTCPlh1TFLTM6VNWCm5p15ZcUcM6omowPh5e8k3g1pIWOlAuRybLE/O3yWJKen2rhYRERERuaDdZ6/LioNXBDvXvd3bfbawy45BPREVqOAAH5k/tJWUCvRVK5C+8OM+Sc/ItHe1iIiIiMiFZGbq5N3lWVvY9W9RWepVCBJ3xaCeiApclVIBMmdIC/Hz9pT1xyLlnWWHuYc9ERERERWYpfsuycGLsVLUz1tevae2uDMG9URkE02rlJDpjzVR6VALd5yXObcWMCEiIiIiuhNJqdjC7pg6/3yXGlK6mHttYZcdg3oispmeDcrLW/fVVeffW3lUVh66Yu8qEREREZGTm7XpjETEpUilEkVkWDv328IuOwb1RGRTT7WvJkPahKrzryzeL3vO3bB3lYiIiIjISV2OuSnf/HNanR93X13x9/ESd8egnohsCquQTuhTX7rVLSMp6Zky/Lvdcu5aor2rRURERERO6KNVxyQ5LVNaVS0p9zYoZ+/qOAQG9URkc16eHvL5402lYcVguZ6YKsPm7ZKYpFR7V4uIiIiInMje8zfk9/2X3X4Lu+y8xQqjR4/OU7mpU6da87BE5AYCfL3l2yEt5IGZ2+RMdKI8890e+f7pVmpveyIiIiIiS7CT0uRbW9g91KySNKwUbO8qOWdQv2/fvlzLsLeEiMwpE+Qv84a1lIdmbpOdZ6/L678clGn9m4inJ48bRERERGTenwcuy77zMRLg6yWv93DvLezuKKjfsGGDNcWJiHKoVbaYzBrUXIbM3akOzlVKBshrPDATERERkRk3UzPkw7+ytrAb2bm6lA3yt3eVXG9OfXp6uiQkJBTEQxGRG2hXI0SmPNhQnf9ywylZvOu8vatERERERA5q9uYzcjk2WSoWLyJPdwizd3WcO6hftmyZzJ8/3+i69957T4oWLSrFixeXe+65R27csO12VTNmzJCqVauKv7+/tG7dWnbu3Gm27OHDh+Whhx5S5TEtYNq0aTatGxHl3SMtKsuLd9dU58ct/U/+ORFl7yoRERERkYOJiEuWrzZmbWH35r11uIXdnQb1WAAvMfH2VlTbtm2TCRMmyNtvvy0///yzXLhwQSZPniy2snjxYrVY38SJE2Xv3r3SuHFj6dGjh0RGRposn5SUJGFhYfLBBx9IuXLc7oDI0bzSraY82LSiZGTqZOSivXLsapy9q0REREREDuSjVcflZlqGNKtSXPo0Km/v6jh/UI+R77Zt2+ovL1myRLp37y5vvfWWPPjgg/Lpp5+q0XxbQafC8OHDZdiwYVKvXj2ZNWuWBAQEyNy5c02Wb9mypXz88cfy2GOPiZ+fn83qRUT5gwyaDx5qJHeFlZSElHS11R16Y4mIiIiIDl6MkV/3XlTnJ/Spz0XZC2KhvPj4eClVqpT+8pYtW+SRRx7RX65fv75cvnxZbCE1NVX27NkjY8eO1V/n6ekp3bp1k+3btxfY86SkpKiTJi7u9shhWlpagT0PuSetDbEt3YZD85ePNZZHv9mptrobOnen/PB0SynqZ9XhiW5hGyNbYvuiwsK2RrbGNuYcW9i9u+ywOt+3cXmpXy7QaT6vtAJoX9bc16pfzRUrVpSjR49KlSpV1MJ4Bw4ckM8++0x/+7Vr19TIuS1ER0dLRkaGlC1b1uh6XD52LGslxIIwZcoUmTRpksnb1qxZU2DPQ+6NbSmngZVFpsZ6ydGr8TLgy7XydJ1M8WJnbL6xjZEtsX1RYWFbI1tjG3Nc+655yO5zXuLjqZNmXhdk5coL4k7tKykpyTZBPUblX375ZRk3bpysXLlSzVO/66679Lfv3r1batd27q2pkAmAefuGI/WVK1dW5zHVwMfHx461I2eHHjd8udmWTGvUKlYGzt0lR2JEdmVUkUl96jLNykpsY2RLbF9UWNjWyNbYxhxbSlqGfPz5VhFJluc6VZcnutYQd2tfcQYZ4wUa1GNRvEuXLsmLL76oAvqFCxeKl9ft1Qd//PFH6dOnj9hCSEiIeq6IiAij63G5IBfBw9x7c/Pv8YHwS08FgW3JtBbVQmT6Y03luYV75MddF6Va6aLyTMfq9q6WU2IbI1ti+6LCwrZGtsY25phmbz0nF2OSpVyQv4zoUlN8fLzdrn35WHE/q96dIkWKyHfffWf29g0bNoit+Pr6SvPmzWXdunXSr18/dV1mZqa6PGrUKJs9LxEVrh71y8n4XvVk8vIj8v7KY1KxeID04kqnRERERG4hMj5ZZqw/pc6/eW9tCfB1zoC+MHnfyRz3s2fPqtRY7ANvuICerSAtfsiQIdKiRQtp1aqV2nceW+xhNXwYPHiwmvePefHa4npHjhzRn0eWwf79+6Vo0aJSo4ZzpXAQuZMn21WVC9eTZP62s/LKz/ulXLCfNA8tae9qEREREZGNfbr6hCSmZkjjSsHSt3FFe1fH9ba007a169ixo1qgrnXr1iq4LlOmjHTt2rVAF6wzpX///vLJJ5+oaQBNmjRRAfqqVav0i+edP39erly5oi+PlfibNm2qTrge98X5p59+2qb1JKI7g87Ct3vXk251y0hqeqYM/26PnI1OtHe1iIiIiMiGDl+OlZ/3ZC2IN6FPPfH05NpKBT5Sf/XqVenUqZOULl1a7Rlfp04dtdUARsNnz56tgv3//vtPBfm2glR7c+n2GzduNLqMDALUj4icj5enh3z+eFPp//UOOXQpVobN3yW/jWgrQUV8ZGf4dZWaVaaYv7SqVlKVJSIiIiJn38LuiCB869O4ArM0bRXUY/u60NBQ2bp1q/j7++uv79mzp4wYMULat2+vymjp70REdwJzqL4d2kIemLFNwqMT5eFZ2yQxJV2uxqXoy5QP9peJfepJzwacd09ERETkrFYfjpB/w6+Ln7envNnTuXdUc+j0eyzL/+abbxoF9IaL6L3++uuyevXqgqwfEbk5jMbPG9ZS/H085XRUolFAD1djk2XEwr2y6r/bU2+IiIiIyHmkpGfIlL+OqvPDO4RJpRIB9q6S6wb1Z86ckWbNmpm9HQvYoQwRUUGqXrqo2ZVPtQk2k5YdkYxMTrchIiIicjYLtp2Vc9eSpHQxPxnRmdsZ2zSoj4+Pl6CgILO3FytWTBISEqyuBBGRJZhDfz0x1eztCOWvxCarckRERETkPKITUuSLdVlb2L3eo7YE+nELO2tZ/Y4hsDeVfg9xcXFcmI6IChwWxSvIckRERETkGD5bc0LiU9KlfoUgebhZJXtXx/WDegTstWrVsng7tqIiIiroefUFWY6IiIiI7O/Y1Tj5ced5dX5Cb25hVyhB/YYNG2xXEyIiM7BtHVa5x6J45nKBMAcL5YiIiIjI8WFA+H/LjwqWRLq3QTlpHVbK3lVyj6Aee9QTERU27EOPbeuwyj36b00F9vE302TTiUjpWqesHWpIRERERNZYfyxStpyKFl8vTxl7b117V8d9Fsrz9PQULy8viydvby5sQEQFD/vQfzWwmZQLNk6xLxvkJzXKFJXk9Ex5asFu+Wrjaa7tQUREROTAUtMz5b0VWVvYPdm+mlQpxS3s7oRVEfjSpUvN3rZ9+3b5/PPPJTMz844qRERkKbDvXq+cWuUei+JhDj1S7rGV3aRlh2XRv+flw1XH1PysDx9qJP4+XvauMhERERFl8/2Oc3ImOlFCivrK8124hV2hBvV9+/bNcd3x48dlzJgxsmzZMhkwYIC8++67d1wpIiJLqfhtqpfKcd17DzSUOuWDZNKfh+WP/ZclPDpRvh7UXMoHF7FbXYmIiIjI2I3EVJm+9oQ6/+o9taWYv4+9q+Re6feGLl++LMOHD5eGDRtKenq67N+/XxYsWCChoaEFW0MiojwadFeofP9UaykR4CMHL8bK/V9ulb3nb9i7WkRERER0y7S1JyQuOV3qlCsmj7aobO/quGdQHxsbK2+++abUqFFDDh8+LOvWrVOj9A0aNLBNDYmIrIBR/D9HtVf/KKLiU+Sxr3fIkj0X7V0tIiIiIrd3MiJeFv57ews7ZFtSIQf1H330kYSFhcny5cvlxx9/lG3btkmHDh0KoBpERAWncskA+XVEW+lRv6ykZmTKa78ckP8tPyLpGVzzg4iIiMhe3lt5VK2F1L1eWWlbI8Te1XHPOfWYO1+kSBE1So9Ue5xM+e233wqqfkRE+RLo5y1fDWgu09edVKc5W8LleES8fPl4MwkO4NwtIiIiosK04XikbDweJT5eHjLuPm5hZ7egfvDgweLhwRQJInIOnp4e8kr3WlK7XDF59ecDsvlktPSbuVVmD26htsEjIiIiIttLy7i9hd2QNlWlWkigvavkvkH9/PnzbVcTIiIbua9heQktFSDPfLdHrYr/wIyt8vnjTaVLnTIqBSz7Fnmc30VERERUcH7ceV5ORSaoxYxfuLumvavj3kE9EZGzql8hWP4Y1U5GLtwrO89elycX7JJ+TSrK9jPX5Gpssr5c+WB/mdinnvRsUN6u9SUiIiJyBbFJaTJ1TdYWdqPvqS3BRTgN0mG2tCMicjYhRf1k4dOt5fFWVUSnE1m675JRQA+4PGLhXln13xW71ZOIiIjIVWBto5ikNKlVtqg83pJb2NkCg3oiciu+3p4yuW99CfI3naiku/V30rIjKjWfiIiIiPLnTFSCfLf9rDo/vlc98fZi+GkLfFeJyO3sOntD4pLTzd6OUP5KbLKaa09ERERE+fP+yqOSnqmTrnXKSMdape1dHZfFoJ6I3A4WxSvIckRERERkbMvJaFl7NFK8PbmFna0xqCcit4NV7guyHBERERHdlp6RKZOXH1HnB94Vyq2EbYxBPRG5HWxbh1XuzW1ch+txO8oRERERkXUW774gxyPi1Ur3L3fjFna25nRB/YwZM6Rq1ari7+8vrVu3lp07d1os/8svv0idOnVU+YYNG8rKlSsLra5E5JiwDz22rYPsgb12Gbdzv3oiIiIi68Qlp8nUv7O2sHulW00pHuBr7yq5PKcK6hcvXiyjR4+WiRMnyt69e6Vx48bSo0cPiYyMNFl+27Zt8vjjj8tTTz0l+/btk379+qnTf//9V+h1JyLHgn3ovxrYTMoFG6fY4zKu5z71RERERNabsf6UXEtMleqlA2XAXaH2ro5bML2nk4OaOnWqDB8+XIYNG6Yuz5o1S1asWCFz586VMWPG5Cg/ffp06dmzp7z++uvq8uTJk2XNmjXy5ZdfqvtaIzk5WRITE8XHx8foei8vL5UFoEEZczw9PaVIkSL5KpuUlCQ6bKxtgoeHhwQEBOSr7M2bNyUzM9NsPQIDA/NVFu9XRkZGgZRFfVFvSElJkfT09AIpi/cX7zOkpqZKWlpagZRFe0C7MFUW5w3bkmFZ3Iby5vj5+Ym3t7fVZfEe4L0wx9fXV9+urSmLzwyvxRyUQ3lry6KNoa0VRFm8B3gvAN8JfDcMdagWJKtHtZaPVx+V73ZckmbVSsuSEW0FA/SWvp/WfO8L+xiRvY0BjxHOc4ywVNZRjhHm/h+64jEiv2Ud+RhhiqMeIyy1NR4jHPcY4Wy/I8y1MeAxwvpjxPlrSTJ3a7h+C7v01BRJdcPfEWlpaUavJT/HCEufXw46J5GSkqLz8vLSLV261Oj6wYMH6+6//36T96lcubLus88+M7puwoQJukaNGpl9nuTkZF1sbKz+dOHCBbRYs6d7771Xl5qaqj8FBASYLduxY0ejsiEhIWbLNm/e3KhsaGio2bJ169Y1KovL5sricQzL4nnMlUX9DMui/ubK4nUblsX7Yul9Myz74IMPWix748YNfdlBgwZZLHvp0iV92eeee85i2RMnTujLjh492mLZffv26cuOHz/eYtlt27bpy06ZMsVi2TVr1ujLTp8+3WLZ33//XV92zpw5Fsv+8MMP+rI4b6ksHksri+ewVBZ11Mqi7pbK4rVrZfGeWCqL91Qri/faUll8VlpZfIaWyqINaGXRNiyVDWxwt67jh+tVWbQ5S2XRZg3bsKWyPEZknXiMcI1jxK+//uq2xwi0La0sjxG2P0b06NHD4vvGY4RjHiOc6XdEYmKi7uuvv7ZYlscI648RRUuV14W+uVw3cPZ2Fb+58++Ir7/+WrWzOz1GICbNjdOM1EdHR6selrJlyxpdj8vHjh0zeZ+rV6+aLI/rzZkyZYpMmjQpz/VC6r/hPH1LvUDXrl0zKmuphzQ2NtaorKXevoSEBKOyuGwOHsewLJ7HHNTPsCzqbw5et2FZc1MiNIZlLX0esHr1an0P5cWLFy2WXbt2rQQHB6vz586ds1h2w4YN+vZx5swZi2U3b96sf7yTJ09aLLt161b96zfXNjU7duzQ98IdPnzYYtndu3frzx84cMBiWUw30XpUcd4SPJb2eRg+hymoo1b20KFDFsvitWtlc3vPcLtW9vz58xbL4rPSykZERFgsi89MK2uprWsuXE+UP5evlPRUy1vZoc3mdX0OHiOy8BjhGseI3Mq68jECbUsra2nEEHiMuPNjBH73WcJjhGMeI5ztd0RueIyw/hhxMy1DPEUn7QIj5K+//nLr3xGALPGCPkaY4oHIXpzA5cuXpWLFimqefJs2bfTXv/HGG7Jp0yb5999/c9wHqTULFixQ8+o1M2fOVEG7uX/0SKMwTBuKi4uTypUrqxT/3r17M/3eRVJi7Jl+v379eunatSvT7x0otRa3tf14s9zM9JJVL7aTsJAAp02by97GgMeI/JVlam3O7z0+Y/xIM2xf7nCMsLasIx8jTHHEY0R8fLz6MWyurfEY4ZjHCGf6HYH3AMFc+/btmX5fAMeIjEydPDZ7pxyPTJBB7WrKO33quvXviLS0NNmyZYta/w3tKz/HCMShoaGhqmMkKChILHGakfqQkBDVqLMH47hcrlw5k/fB9daUB3wZtS+kIXyZihcvbvZLr0GZvLKmrNYjVNBlc3s9LFuwZfEFNdeWcNnwgJnb41pT1vCgXZBlDf/JFFRZMPUdLIiy2j98U6qXLyH/XYqT8zeSpU6F4hbLFtb3Pj/HCEttLHtZZ/0euXNZex8jrPl/6GrHiDsp60jHiIIua6v2XqxYMavami3q4Ixl7X2McLbfEQii8tLGNDxGmC+7eNd5OXkjXYKLFZVXe9TRv6eO8t0o7LL4PYb2hTLaydrH1YJ8l1r9Hl+M5s2by7p16/TXoScHlw1H7g3hesPygF5fc+WJyH1VCymq/oZHW7EoCREREZGbS0hJl49XZ21h99LdNaVkILewK2xOM1IP2M5uyJAh0qJFC2nVqpVMmzZNpZVoq+EPHjxYpehjXjy89NJL0qlTJ/n000+lV69e8tNPP6m5Pt98842dXwkROZqwkKw0rTNRDOqJiIiI8mrmhlMSnZAiVUsFyOA2Ve1dHbfkVEF9//79JSoqSiZMmKAWPGjSpImsWrVKv0gJFsYwTFNo27at/PDDDzJ+/HgZN26c1KxZU37//Xdp0KCBHV8FETmisNJZQT1H6omIiIjy5sL1JJmzJWsLu7d61RNfb6dJBHcpThXUw6hRo9TJlI0bN+a47pFHHlEnIiJLwm6l35+JNr+iKxERERHd9sFfxyQ1PVPaVi8l3eqWsXd13Ba7UoiIMKf+1kh9dEKqxN40vzopEREREYnsOntdVhy6Ip4eIm/3rqdfFZ4KH4N6IiIRKernLWWKZa2WyxR8IiIiIvMyM3Xy7rIj6nz/llWkbnnLW66RbTGoJyLKNq/+TBRT8ImIiIjM+W3fJTl0KVYNiozuXsve1XF7DOqJiG4JK81t7YiIiIgsSVRb2B1T50d1rSGlb2U6kv0wqCciuoXb2hERERFZ9vWm0xIRlyJVSgbIsHbcws4RMKgnIsqefs+ReiIiIqIcLsXclK//OaPOj723jvh5e9m7SsSgnojotmq3trULj05QC8AQERER0W0frTomKemZ0qpaSenZoJy9q0O3MKgnIrqlcoki4u3pIclpmXIlLtne1SEiIiJyGHvO3ZA/9l8W7Fw3gVvYORQG9UREt3h7eUqVUgHqfDjn1RMREREpyGCcvDxrC7uHm1WSBhWD7V0lMsCgnojIQNitFPwz0dzWjoiIiAiWHbws+y/ESICvl7zeo7a9q0PZMKgnIjJQXb9XPUfqiYiIiG6mZsgHf2VtYfd8lxpSJsjf3lWibBjUExEZqKZta8cV8ImIiIjkm3/OyJXYZKlYvIg81b6avatDJjCoJyIyEFb6Vvp9FNPviYiIyL1djU2WWZtOq/Nj7q0j/j7cws4RMagnIjIxUo99WJPTMuxdHSIiIiK7+Wj1MbmZliHNQ0tI70bl7V0dMoNBPRGRgZCivlLM31t0OpFz15LsXR0iIiIiuzhwIUZ+23tJnecWdo6NQT0RkQH8wwq7NVofzhXwiYiIyA3pdLe3sHuwaUVpXLm4vatEFjCoJyIyM6/+NFfAJyIiIje04tAV2X3uhhTx8ZLXe3ILO0fHoJ6IKJvbI/UM6omIiMi9YE2hKSuztrB7tlOYlA8uYu8qUS4Y1BMRZVNNv1c90++JiIjIvXy7JVwtGFw+2F+e7Vjd3tWhPGBQT0SUTVjIrW3tOFJPREREbiQyLllmbjilzr/Zs44U8eUWds6AQT0RUTZVQwLU35ikNLmRmGrv6hAREREVik/+Pi6JqRnSpHJxub9xBXtXh/KIQT0RUTYBvt5SIdhfnT/DFfCJiIjIDfx3KVZ+2XNRnX+7dz3x9OQWds6CQT0RkYUV8M9wBXwiIiJyky3sdDpRI/TNQ0vYu0pkBQb1REQmVLu1Aj7n1RMREZGrW334qvwbfl38vD3lzXvr2Ls6ZCUG9UREJoTdWgE/nCP1RERE5MJS0jPk/Vtb2D3TMUwqFucWds7GaYL669evy4ABAyQoKEiKFy8uTz31lCQkWJ7r+s0330jnzp3VfTw8PCQmJqbQ6ktErjJSzzn1RERE5Lrmbz0r568nSZlifvJcJ25h54ycJqhHQH/48GFZs2aNLF++XP755x955plnLN4nKSlJevbsKePGjSu0ehKRa6h+a0792WtJkpGps3d1iIiIiApcdEKKfLE+awu7N3rWkUA/b3tXifLBKT61o0ePyqpVq2TXrl3SokULdd0XX3wh9913n3zyySdSoYLp7RZefvll9Xfjxo15fq6UlBR10sTFxenPp6Wl3cGrILrdhtiWHF/pQG/x9faU1PRMORcdJ5VLZG1z5+jYxsiW2L6osLCtka2xjWX5ZPUxSUhJlwYVgqRPgzJu/344Uvuy5r4eOix16ODmzp0rr776qty4cUN/XXp6uvj7+8svv/wiDzzwgMX7I6jv0qWLuj9S9y155513ZNKkSTmu/+GHHyQgwDl+1BNRwZiy30uu3vSQ5+pkSN0SDn+oJCIiIsqzS4kiHx/0Ep14yIv106V6kL1rRNmzzp944gmJjY1V08mdfqT+6tWrUqZMGaPrvL29pWTJkuq2gjR27FgZPXq00Uh95cqV1fnu3buLj49PgT4fuRf0uGEKCduSc1gRu1+uHomUkLB6cl+bUHEGbGNkS2xfVFjY1sjW3L2NYVx36Pw9opPrcm/9svLCY43tXSWXklYA7cswYzw3dg3qx4wZIx9++GGuqfeFyc/PT51MwQfijl96KnhsS86hepliIkci5dz1ZKf7vNjGyJbYvqiwsK2RrblrG1t7JEK2nbkuvl6eMq5XPbd8Dxy9fVlzP7sG9UipHzp0qMUyYWFhUq5cOYmMjDS6Hun3WBEftxER2UIYV8AnIiIiF4P1gt5bmTVw+lSHalK5JKcYOzu7BvWlS5dWp9y0adNGbUe3Z88ead68ubpu/fr1kpmZKa1bty6EmhKRO+Je9URERORqvtt+VsKjEyWkqJ+M7Mwt7FyBU2xpV7duXbU13fDhw2Xnzp2ydetWGTVqlDz22GP6le8vXbokderUUbdrMN9+//79cupU1jYNhw4dUpcxwk9ElJuwkKxt7S7HJktSarq9q0NERER0R64npsrn606q86/dU0uK+TPt3hU4RVAPixYtUkH73Xffrbaya9++vXzzzTdGixEcP35crRKomTVrljRt2lR1BkDHjh3V5T///NMur4GInEuJQF8pEZD1z+5s9O1jCxEREZEzmrb2hMQlp0vd8kHySIusxcDJ+TnF6veAle6xrZw5VatWVas4Zt+eDiciovyqFhIoN87HqHn19SpwrxciIiJyTicj4mXRv+fV+bd71xUvTw97V4ncbaSeiMgewkpnpeBzXj0RERE5s/+tOCoZmTq5p15ZaVs9xN7VoQLEoJ6IKJeRejgTzaCeiIiInNOG45Gy6USU+Hh5yLj76tq7OlTAGNQTEVlQ/dYK+GeiuK0dEREROZ+0jEz53/Ij6vywdtWk6q0BC3IdDOqJiCyodmsFfIzUZ1+3g4iIiMjR/fDveTkdlSglA31lVNca9q4O2QCDeiIiC0JLBYiHh0h8crpEJ6TauzpEREREeRaTlCqfrT2hzo/uXkuCuIWdS2JQT0Rkgb+Pl1QqUUSdD+e8eiIiInIi09edlJikNKldtpg81pJb2Lkqp9nSztFlZGRIWlqavatBDg5txNvbW5KTk1WbcXY+Pj7i5eUl7pCCf+H6TTWvvlW1kvauDhEREVGuTkclyPfbz6nz43vXFW8vjue6Kgb1dwhzbK9evSoxMTH2rgo5SXspV66cXLhwQTyQ0+0Cihcvrl6Tq7weU8JCAuWfE1EcqSciIiKn8f6Ko5KeqZO765SRDjVL27s6ZEMM6u+QFtCXKVNGAgICXDqwoTuXmZkpCQkJUrRoUfH09HT6DoqkpCSJjIxUl8uXLy+uKuzWCvhYZIaIiIjI0W0+GSXrjkWKt6eHjOvFLexcHYP6O4D0aS2gL1WqlL2rQ04S1Kempoq/v7/TB/VQpEjWXHME9vgeuGoqfph+BXxua0dERESOLV1tYXdUnR/UJlSql876HUOuy/mjCjvS5tBjhJ7IXWnt35XXlNBG6s9fS1L/KImIiIgc1U+7LsjxiHgpHuAjL91d097VoULAoL4AMOWe3Jk7tP9yQf7i7+Op5qVduHHT3tUhIiIiMin2ZppMXZO1hd0r3WpJ8QBfe1eJCgGDeiKiXHh6eqgV8CGcKfhERETkoGZsOCXXE1OleulAeaJ1FXtXhwoJg3pyKu+88440adKk0J936NCh0q9fvzt6jPnz50vJkiUL5PW9/fbb8swzz4it3XXXXfLrr7/a/HmcZQV8OMPF8oiIiMgBnY1OlHlbw9X58b3riQ+3sHMb/KTdUFRUlIwYMUKqVKkifn5+ajuyHj16yNatW8XR0rp///13e1fDIXdcmD59urz11ltG18+YMUOqVq2qFuFr3bq17Ny5M8+P+dNPP6n3O3vHxfjx42XMmDFqgT93xxXwiYiIyJG9v/KopGXopFOt0tKldhl7V4cKEYN6N/TQQw/Jvn37ZMGCBXLixAn5888/pXPnznLt2jVxRVht3pXMmTNH2rZtK6GhofrrFi9eLKNHj5aJEyfK3r17pXHjxqqjRttuzpKzZ8/Ka6+9Jh06dMhx27333ivx8fHy119/iburdmuknun3RERE5Gi2nY6Wv49EiJenh4znFnZuh0G9DSQmmj8lJ+e97M1s63GZK2cNbMG3efNm+fDDD6VLly4qMGzVqpWMHTtW7r//fn05jNoieHzggQfU6uY1a9ZUwb+hTZs2qftitB97lGNENz09Xd22fPlyKV68uNr2D/bv368eE2U0Tz/9tAwcONBkPTHiDHh+3E+7rPn+++/VdcHBwfLYY4+pwFODDopRo0bJyy+/LCEhISq4hf/++08FqdgjvmzZsjJo0CCJjo7W32/JkiXSsGFDtU0btijs1q2bJGZ7gz/55BP1WnH7888/b7Ti+40bN2Tw4MFSokQJ9Z7huU6ePGnx8/jggw9UXYoVKyZPPfWUJGdvIGZG1fv06WN03dSpU2X48OEybNgwqVevnsyaNUvVYe7cuRYfC5/PgAEDZNKkSRIWFpbjdmxRd99996nndHdht7aDYfo9EREROZKMTJ1MvrWF3YDWVaRm2WL2rhIVMgb1NlC0qPnTQw8Zly1TxnzZe+81Lou41lQ56+pWVJ2Q1p6SkmKxLAK9Rx99VA4ePKgCOwR/169fV7ddunRJXdeyZUs5cOCAfPXVV/Ltt9/K//73P3U7Rn0RaCMjQOsAQIC9ceNG/ePjOgTgpuzatUv9nTdvnly5ckV/GU6fPq3qj44DnPA4CI4NIQvB19dXTSlAgIvOjK5du0rTpk1l9+7dsmrVKomIiFCvD/Acjz/+uDz55JNy9OhRVc8HH3xQdDqd/jE3bNignht/8fiYI4+T4bx7PDY6P7Zv367ui/fI3FZvP//8s5pD//7776v7obNg5syZFj8TvP9HjhyRFi1aGGUi7NmzR3VCaDw9PdVl1MOSd999V+0vjw4Fc9Bxg44gd6eN1EfGp0hCSlbnFREREZG9/bL7ghy9EidB/t7ycrda9q4O2QGDejfj7e2tAlEEpRhJb9eunYwbN04F7tkhSEWgW6NGDRV4JiQk6OdpI/isXLmyfPnll1KnTh01FxudAJ9++qmaf40RdCz4pgXx+PvKK6+oIB+Pg06BU6dOSadOnUzWs3Tp0uov6og5/9plwOPjNTRo0EB1HmDEfd26dUb3R2bBRx99JLVr11Yn1BMBPV4H6ovzGMVGgI4pCAjqkWWAQB4ZABixHzlypOoA0WAEXnu9vXv3ll69eumfFyPyCOaR3YA6If190aJF6nWaWxdg2rRpKpjGCXVEhwhG2S05f/686iyoUKGC/jpkG2DEHSP+hnAZ8+/N2bJli+qImT17tsXnxHNduHDB7efVBxfxkZCiWdvChHO0noiIiBxAfHKafPL3cXX+pW61pGQgt7BzRwzqbSAhwfwp+0LimPJsrmz2acxnz5oul5859ZcvX1ZBaM+ePVXA3axZM6NRZ2jUqJH+fGBgoAQFBennaGM0u02bNkZ7lKODAAH7xYsX1WUE7HhsBKEY6UXAXLduXRVMYnQdwSKCb2sh6Ea6ugYj3Nnnjjdv3tzoMrIJEMBrmQo4ITgHjL4jCL/77rtVMP/II4+oQBfp9Ibq16+v0tFNPS/eD3SYYIE6DVL0EazjNlNwvWF5wHtqyc1bczKwGN6dQBYFOkPwOpFBYQmmIyCgzy2zwx2E3drW7gzn1RMREZEDmLnxtEQnpKqMwkF33V5vidyLt70r4IoCA+1fNjcICrt3765O2B4N89uxyBpG5zU+Pj5G90EAb81oLVLrMRqOgBqPhSAa1yHQR8BsbpQ+N3mpFzohDKGzAfPQsZZAdgjOEayvWbNGtm3bJn///bd88cUXanX5f//9V6pVq5bn57U1LQDH+6dlL+A61B/TCQzhMrIcTEFHBhbIM5ybr70WdE4cP35cqlevrk/5x/uJ4N7d4R/mzrPXOa+eiIiI7O7C9ST5dnPWFnZv3VdXfL05Xuuu+MmTgrTv7IvCWYIRd23euAbz1zGCXqlSJaN59Z999pk+gNeCepzMzafXIIjWFtq7U8hEOHz4sBrlx3QCw5PWAYAgHdkGmEaAaQKYk7906dI8vx9I30cngAa7CSA4NpdSj/sYlocdO3ZYfB4E2siYwLx6DeqJzATDKQgI0HHZ3Mg/OlgOHTqkFjDUTlgoEYsn4jymVmiwwCCmK9Dtbe3ORDOoJyIiIvua8tdRSc3IlPY1QuTuutzCzp0xqHczCDSxYNzChQvVPPrw8HD55Zdf1Pzzvn375vlxMN8c86xfeOEFOXbsmPzxxx9qpB/bqmGRNm0OOlL4MbdcC+A7duyotlzDPPbcRuoRgCMwxbzw7Knw1sJK9RhxxhoBWHQPI9WrV69Wq8Wj4wDBtbZgHeat//bbbxIVFaUC77zANAK8f1iBHtMLkJ2Alf0rVqxo9n196aWXVCYDFgPE+4H3Dx0PlmgL4OE5DOF9Ryo91kpAWv+IESNUJw1enwYr82OXAy1TA2sSGJ6wfgE6ZXAeHQUaTJ2455578vQ+uMsK+NzWjoiIiOxpZ/h1WXnoqnh6iIzvXddoSiy5Hwb1bgZzyTGPG6PnCLARwCH9HsEoFoHLKwSrK1euVAvnYT76c889pxZ8Gz9+vFE5BO4ImrWgvmTJkmrkGmnhmG9uCRbdQ0o8Ro3vdKQY8/eRSYC6IEDF3HlseYdAFoEyRr//+ecftVp9rVq11OvA82NburxCcI4RcyyihxFyZDHgPcqetq/p37+/eu/feOMNdb9z586pYDw3mCqBLeYMU//xWNhub8KECWqBQoy2Y4V/w8Xz0FmBBQGtgYX+MCXBsHPAnen3qo9KNMpSISIiIiosmZk6eXd51kDQY62qSJ1yQfauEtmZh85JfplilBWjwsuWLVNBGBZ7mz59utHq5NnLY+QT86MRzGD+MVZonzx5slqZPa/i4uJU+R9++EEefvhhowANe4pjpBtzru904TJyDwjE0abQiaBlNFgLX1l0zGA3AWQe2NKbb76psiS++eYbs2Xc6XuQmp4pdSesUvvB/jvubikb5HivF1soojMJHVTmOpSI8ovtiwoL2xrZmjO3MWxh9/qSg1LMz1s2vN5ZQor62btKZIP2pcWhsbGxKnZwiZF67JGO1GSM3GJvcoyqPvPMM2bLY3V3nDB6iTnBWNkdI5eW9uMmcgZIr0KQjTn8toY97NERRlmwAE3lElkLBp6OYgo+ERERFa7ElHT5eHXWFnajutZgQE/Os/o95ggjIMdc6BYtWqjrsDo5ej4QtBvu2a1BWvmvBvvHYYGx9957T81zRjCEFb6JnBVS7HGytVdffdXmz+GM8+rPXkuS8OhEaVvd8naARERERAVp1qbTEhmfIlVKBsjQdlXtXR1yEE4R2WKVdcx91gJ6wGJhSF/GAmcPPPBAnh5HS12wFNBjL27D/biR9mCYRmEIl5EKjZTqwt7ajJyTNttFazeuAK8DrwffB2yt5+pCS2aN1J+KiM9xTHAEWp0csW7k/Ni+qLCwrZGtOWMbuxxzU77554w6/2aPmuKpy5S0NNf4Pelq0gqgfVlzX6cI6rH6OdKADSEwx6JruC0voqOjVRqxpZR9mDJlitrSzBSk/mevAxZ8wx7oqampeaoHEWCrP1eBtn/z5k01JaYwpgTYW/xVrC7rJX/vD5fAmDNSPUinVp51NNmPV0QFie2LCgvbGtmaM7WxBSc8JSXdU2oE6SQtfI+sPGvvGpEt21dSUpJzBPVjxoyRDz/8MNfU+zuF0fZevXqpVdffeecdi2Wx5Re2BzO8r7Znd/fu3XMslIdt3bBYn6svEEYFAyPaCOixdZyrbD2C70GRIkXUbgqu/j1YfThC1h84gr5TuZTkKV8eESkX5Cfj76sjPerf3mnAntCri38g2Y9XRAWB7YsKC9sa2ZqztbG952Nk7/adgp+PnwxoI/UrcMV7V29fcQYZ4w4d1GO+7tChQy2WCQsLU6PhkZGRRtdjRBAr3OM2SxBA9ezZUwVRS5cuzfVN9fPzUydTcF/D+2N7NARmmAaQ35XMyb1oKfdau3EFeB14Pdm/H65m1X9X5IWfDkj27UIi4lLU9V8NbCY9G5QXR+HqnwfZF9sXFRa2NbI1Z2hj2MLu/VUn1PlHm1eWJqGl7F0lKoT2Zc397BrUY5s5nHKDPb9jYmJkz549aj9vWL9+vQqQsLWXpd6NHj16qCD9zz//dPlRRCKyDWxhN2nZkRwBPeA65Fzg9u71yomXI+biExERkdP688BlOXAhRgJ9veTVHrXsXR1yQE4xVFi3bl012j58+HDZuXOnbN26VUaNGiWPPfaYfuX7S5cuSZ06ddTtWkB/zz33SGJionz77bfqMubf44QRdiKivNoZfl2uxCabvR2BPW5HOSIiIqKCkpSaLh/8dUydH9mlhpQpxkFKctKF8mDRokUqkL/77rtVuu9DDz0kn3/+udG8hePHj+sXFNi7d69aGR9q1Khh9Fjh4eFStSq3gCCivImMNx/QG7oae9PmdSEiIiL3gdXur8YlS8XiReSp9tXsXR1yUE4xUg9Y6f6HH35Qc+SxNd3cuXPVAnUaBOlYhKxz587qMv7isqkTA3rr4L18+eWXbfocWMCwMPZdtzVkhSBDxNaQpfLpp5/a/HkoS157xaf8dUx+2nleUtO5vUxhTo3Yfvqa/LH/kvqLy0RERK7gSuxN+XpT1hZ2Y++rI/4+rr91MLl4UE8FB4sTYmGz5557Lsdtzz//vLrNcAHD3377TW0HSLmvAv/222/LxIkTja7/5Zdf1NQQrOnQuHFj+fvvvy0+Dt5vrJSJ9SaCgoLUmhKrV682KjN+/Hh57733VAcX2V6raiWlfLC/mjtvDlajjYxPkTG/HZLOH2+Q77afleQ0TvWx9eKF7T9cL4/P3iEv/bRf/cVlXE9EROTsndVv/HJQbqZlSIvQEtKroeMsxkuOh0G9m8I2fT/99JPaX9wwKEU2RJUqVXJkSWD3ALJsyZIlKghv166d/rpt27bJ448/Lk899ZTs27dP+vbtKwMHDpT//vvP7ONgv3cE9StXrlSLQ3bp0kX69Omj7q9p0KCBVK9eXRYuXGjz10WiFr+b2KeeOp89sPe4dZrWv4mM71VXShfzk8uxyTLhj8PS4aMNMmfzGTUfjgoWAvcRC/fmWOvgamyyup6BPREROXtn9eZT0er6bvXKusxWyGQbDOptIDU11ewJW/HltSzWCchL2fxo1qyZCuwxKqzBeQT0TZs2NZt+f+zYMQkICFDBv+bnn39W+5QfOYL9u0XtVPD000/rR5q7du0qBw4cMHrMDz74QMqWLas6CxDwokPBkhs3bsiAAQPUY+K5atasKfPmzVO3nT17Vh3oUH8EwKgfRsS3b99u9Bi//vqr1K9fX+2GgCkYhunrX375pQqUNb///rt6zFmzZumv69atmxohNwedJAi+DU2fPl0t8vj666+rBR/fffddVbcZM2aYfZxp06bJG2+8IS1btlSv8/3331d/ly1bZlQOz4XnpMKB7eqwbV25YONUfFzG9X2bVJSnO4TJ5je6yLt960uFYH+Jik+R/604Ku0/3CAzNpyS+GTj7zTZZjcCwO1MxSciImfvrIYP/zrGzmpyjYXynMmUKVPM3obg7IknntBf/uSTT3IE75rQ0FCjNHgEiNpCgIayp3vn1ZNPPqkCYwTLgHUKhg0bJhs3bjR7H6SRo84jR46U9u3bq0ULkcb/4YcfSr16WSOZjzzyiAq8//rrLwkODpavv/5aLXB44sQJNeqPTgDMoUdgi8f4/vvv1aKHYWFhZp8Xae3oNMBjhoSEyKlTp4yyDOCtt95SdcN7jPMYIUc5b29vNeL96KOPquft37+/GkHHayhVqpR6jzt16iQvvviiREVFqY6DTZs2qefBe4HXh88InQRjxowxW8ctW7bIoEGDjK7DfUaPHm10HTo5Vq1aJXmFrRuxlgTeO0OtWrVSKfgpKSmqo4IKJ7DHtnVY5R6L52GuPVLzDbexw3y3wW2qymMtq8jSfRdlxobTcv56kny8+rh8vem0DGtXTYa1qyrFA3zt+lrcZTeCNtXvfC9fdA5Y+syJiPKKxxOytrNaw61zyRIG9W4MaeBjx46Vc+fOqcvYKhAjv5aCekAwjNRw3N/X11eNKL/wwgv6wBbbCkZGRuoDTQTaGPlGevozzzyjRqIxOo8T/O9//5O1a9daHK0/f/68yiBo0aKFumxqscPXXntNevXqpc5PmjRJjcojqEdHxNSpU1XHAjoHoFatWqqT4OOPP1ZBPUbpETQjmH/44YfVe/Dqq6+qjhTAa0Jg37ZtW5P1Q3YC5rdrWyxqsIUiMhIModMA1+cV3r+EhATVKWEIz4VMDTwWOoCocOCfaV4CRV9vT+nfsoo81KySLDt4Wb5cf0pORyXK9HUn5dst4TKoTahaxTakKDtkbLUbgbly1vyoxsgIfkgZdiJgfQVMx0AnDxFRXvF4QuasPxZRqJ3V5HoY1NsAAmVzMLKdPRA1J/vcmZdeekkKEoJLBMHz589XuwLgPEan8wKj+giM8XoOHz6sryvS7BGAYgTcEEbVT58+rc4fPXo0xyJ9WAxuw4YNZp9vxIgRahtDbFWI1eX79euXI8Bu1KiR/nz58ln/HNG5gKAez4n57IYw9x0dDBkZGeLl5SUdO3ZUwTzS7BHwo/Pio48+UlMOEOyj8wKp/aZoWQNYDK8gYZoDOij++OMPKVOmjNFtyIYAU9kb5Di8vTzlgaaV5P7GFWXVf1fli/Un5djVePlq42mZtzVcBrQOlWc6hknZIO47W9C7EZgqZ82Pai0VMvvIiTZvH9Mu+EOciPKCxxMydCMxVf4Nvy47zlxTJ/wuKMhObXI/DOptAKPX9i5rTQr+qFGj1HlL87yzQ/CemJiogvorV67og2gE9DhvarS/ePHi+a7nvffeqzIKkCGwZs0aNeqOlfoxiq3x8fHRn9c6GZC6nldYO+Cbb76RzZs3q6wArAegBfoI6pGibw46MfCcmPtvqFy5chIREWF0HVL8cX1ukDWBtQmwej46GrK7fv26vnOGHB9Ggns1Ki/3Nign645FquD+4MVYNWr//Y5z0r9FZXm2U5hUKmG644hy7kZgblTD49ZaByiX3x/Vuc3bx3MwFZKI8oLHE8pvEJ/fTm1yPwzq3RwWcUMKNwLSHj165Ok+CCaRso556wjoMScfI+gYOcYCfEgHxzx2UynygAXj/v33Xxk8eLD+uh07duT6vAhehwwZok4dOnRQi88ZBvWW4DkxvcAQLiPbAKP0gKAdCwIiiEaAD/iLqQEoi3R8Sx0uWFMAI/yG+9QjA2HdunX6hQYBGQl33XWXxfr++OOPqsMFgb02pSA7rKBfqVKlPGdXkGPw9PSQ7vXKSre6ZeSfk9HyxbqTsvvcDRXY/7jzvDzYrKKM7FxDqoYE2ruqDgs/eu9rWF51iJiDkXfDH8d5WVzv9SUH1Q+t9AydnLuWyFRIIrpj6RmZsmD7WR5P3Mz1xFTZGY4A/rrZIL5mmaJyV1gpdWoeWkIemLlVdTLrrOisJtIwqHdzCGiRmq6dzwukzmPlfKwEj0XaMKqNaQQY6ceIMgJZpMcjdR1B8+XLl2XFihXywAMPqDnxmEaATgGcRwr8okWLVAq/pYXyJkyYIM2bN1fz5PGcy5cvV4F6XiEgR/r85MmT1UJ5WMAOK97PnDnTKH2/RIkSKuUdj68F9Xht6PQw3KrOFHSKYE0BwwAerxWdBVhpH8E5gvX9+/fLnDlzjKZrXLp0Sb777jt1Gc+PjgvM52/durV+/j06TbDwoAYZBYYdCORc0KY61SotHWuGqH/6GLnfdvqa/Lz7oizZc1Gtpv98l+pSo4zjbCfpKAs8XbyRJD/vvqDOF/XzkoSUDKPbkfGQPY01t8X1ID45XaatPWlVXb7+57RcS0yRZlVKqOwBbjnkGhylrbsqV39/Y5JSZdOJKFl3NFI2Ho+UuOS8bWvK1Gr3CeLR5rEFbvbOaGSN4ZtgGNh7mOmsJjLEoJ5UmnleIfBECjz2TMdoPE7YKx2r2Pfu3VulyeN2jOJjJX0t1Rxp7NqCcQiqMb8e27ZhcTzMlcec+dWrV1scCUfwi+3rENxipN6a7dyQQYBV99E5gMAeUwSwvZzh7gL4MY7HRQcEXo8W6OP9qV27tgQGWh45xcJ/6KjAgnla8I15/wjS0QEybtw4tTI/3i/D7fOQ7YCFADWYAoCtDzG9ACcNAn2sfwB437D4oDWr6JNjQrvDyAxOe85dVwvqbTgeJUv3XZLf91+S+xqUl+e71JB6FfL+PXXlBZ4w6vXyT/tVAN60SnH5cfhdsu98jPoxvPZopCw7cFn+uxSX7x/L+Bxqly0mUfHJsuJQ7gtabjwepU5QNshPmlYuIc1Ci0vTKiWkYcVgtRsCORdHaeuuyhXfX6xLdDoqQR2D1h+NlN3nrovhjpqmOh9NmbP5jJQu6qeOQ+wgdP4gvlbZotK6mvkg3tzWudm/H+Wc/PtBhcNDhyMRmRUXF6cCNARmWBXdcN42Aqvw8HCpVq1agS+QRs4J2/mhA8HcYomY4482hY6C7IsmWuOrr76SpUuXyt9//y32xu9BwTt0MVaN3P995PZ6DN3qlpUXutaQxpUtr02BXRrQsXbfffcZHa/uhLm56NpPTlsv8GQ4qrf5RLQs2XtRivp5y8oXO0iVUrfXILhwPUk6fbxB/ZheO7qjUZbD9tPR8vjsf3N9LnQS4Ac1nrP9h+vNpkJC8SI+ap2EAxdj5OiVeHUfQ96eHqozBqP46IBAwF+5ZBGn/rFui/blSOzd1l2dNe+vo7e1lPQMdVzCaPz6Y5Fq61JD6By8u24ZdWpYsbg6Nlk6nhhqVClYnu1YXXo24Bx7W7KmjV1LSFGfd9ac+OtyPMJ0EG84Ep/f3W1cPZPFXaQVwDFMi0MxYJjbICxH6okKELbIW7Zsmc2fBweHL774wubPQ/bRsFKwfDO4hRy7Gqf2uV9+8LKsPRqhTh1rlVbBfcuqhTOvLi9z0Sf+eVg61CwtgX7ehTKqB4+0qGgU0EPlkgFyd92ysuZIhHy86rjc16i8+kHUsmoJ9d5Zkn2+In5A5ZYK+cFDDfUByM3UDDl4MUb2XYiRveduyN7zMRKdkKIWQ8Rp/ras+4QU9ZUm2mh+5RLSuHKwBPjyX7Ej4GJmtuUK7y++0xuORapAfvPJKElMvT367uvlqToEEcR3qV1GHY8M5XY8ebdvAzkZGS+Ld11Qx4znf9groaUCZHiHMHm4eSVm/RSywgzi87t1LpEh/pIgKkBYHPCFF16w+fNgVXxyfXXKBckXjzeVl7vVlJkbTqt0/H9ORKlT62ol5cW7a0pbG6dp5mUuekRcitSfuFr8vD0luIiPFA/wUX9xCrr11/Bk6nY/b688j+rB/K3nVFpj9lHTeuWDVFC/+tYJAny9JMngx3de5ytakwpZxNdLWoeVUidAEtylmJsquN93PivIP3I5VqITUvUdNIDnw4geRvK1Ef1qIYFOPZrvrHJr69piZkPm7pQaZYrq2zH+Fi/iK8Ha5VttGttZkvXvryMtFofv8ZErcSqlHruWICvHML8V6dRda2eNxrerEWKxYzOvx5OX7q4pC7afk++2n5Vz15Jk/O//ybS1J2RYu2oysHWoamdknyAex+q7wkrqg/hSBRTEExUEBvVERA6ueumi8umjjdWPva82nZYley6orXEGzPlXmlUpLi90rSmda5e2SSBozcJNKemZEhmfok7W8vfJ6hDICop8pZi/t1o40FKqavZRPXQCfL4u50J3WkD/RKsq0rFWiFXzFXEdnsPaVEh8FtieEKf7G1dQ1yWnZcjhy7FqDQCc9p6/oeqBoAGnRf9mra2BILFpZS3IzxrNL+bPH/K2EpecJltPRsuCbWfzVH7LqWh1yk0xP2/VaaWC/luBv/6yUeeWr74MLhfx8XKqTp28pAonpqTLr3sv5unxIuPsu1gcvqfbTkfr0+qzd0RgrYyudbIC+QYVgtWOJnmVl+MJAsXR3WvJsx3D1IKgczaHqw7Cj1cfl5kbTsnjrarIk+2rSYXiRQr0dbtjEL/1ZKQsOeMpX36xVU5GJuYowyCenAmDeiIiJ4F08ykPNlTp99/8c0ZtgYcR4GHzd0mDikEyqktN6VKzYNPy87on7pwhLdQPoNibaTlOMUlZf+NM3IaACiNfyWmZkpyWokb98yL7qJ6l1F7NhuORMrlfA6uD9IJKhUT6bPPQkuqkuRJ781aQnzWaf+hSrHq/sFgiToD4rlYZ49F8dPRYE0yQ8ejriYgE1R6QSr3n3A1Jz7YegiUDWleW4gG+6nOKQTu+1b5jbqaq67CII8SnpKsTAjJrII0bo7HaiH9WsJ8V+BtmuqAOWmYArkPHT2Gnree26B1e+3fbzqpjVV5XgP9yw0kp6u8tHaqXkMKCue7rjkWoEfmtp6PV8ciww7F9jdIqiEcwXzboztaOyevxBKP+anT+rlBZcfCKzNp0Wi3GNmdLuMzfdlbtkIKdPmqVzblDCudkm546cXsk/po6BmRBRk1WQM8gnpwZg3oiIieDEZp37q8vI7tUV6M4C3ecUyu+P7dwj9QqU1TaBHtIj0ydFMTYrvphE+gr1xJTLc5FxxxS/GisbOXjZ2bqVBCUPdjHFIPFt7aty0smQV6mCRh2AjhKem/54CJSvmERua9hVpZAanqmHL0Sp0bxtdH8izduqlRQnH7alfWeIJOhiX40P2t+PtNyzcNI8dZT0aqjZNPxSLmcra2ElQ5UW0z+sf+y3EhMtbhP9Lt9G1oMkLBDA9o0An5sbRZzq0NLdQKojoDUrHZ+q1MAZWJv4juQKmkZOknNyJSo+BR1sgY6f5AdgGDfcIqLfnoALuuzBG51CNwql5/52uamxyBAfm7hXrVuxIELsfoFJKuWClArhuO9sdSFghHTpxbsltpYOTzIQ+7JyJSCXicPx52Dl2Jl/dEIlVZ/+LLxjhkVgv3V+hxd65aRNmGl7Dqf3cfLU/o1rSh9m1RQ2+QhuEd6ODIfcLq7Thl5rnN1/Torrri7QH6D+H9vrUyP08lILYi/DW2srEecPNq5qbStWUZKBvrapa5EBYFBPRGRk8IIzLj76spznarL3C3hKn34RGSCnIj0kn8+3yrPd62pfgjiR2F+pWdmmr1/Qeydi9FmNRc5W0CKH1d5Ceq1TIK8ThNw9H2gfb091Q4HOA1rd7vOWso+RvSxiBYCo80no9VJU710oErX1wJ9jOC56+ictsUYthvEiDw6cxAwa7D+A9aj6FKnjHSuVUa/6CLWqrjTfaIxl75EoK86iQRaVWdMFdFG/1Xwr2UD6DNesjoEtA4CrRMsISVdZbxgNByn89fzN/1FWxvAMEtA6wAw7CBApxIWyLS0gObeczHqL97np9pXUx1/fx+5muvik2eiE2Xh9nNyPCJBjkd4yYbpW2VE5xryUPOKRmtvWDsajfdoy0mk1Ueozh0Effrn9xA15QWBPEbkMWLraFMgUJ/Otcuo0/4LMfL1ptOy6vBV1SmBU/PQEmpR0K83nTHZ0YL33ZV3b8hLEF+nXLFbC9uVlFbVSkkxXw+1OnmP+mUdcocFImswqCcicnIIgF/rUVuGdwyTuZtPy+x/Tkn4tSR57ZcDMn3dCRnRKecP4ryatfGMXI1LVj/iseCcYXq8LffOxQ90jC6Z2wIq+2r1eZ0mkNdyjgR17lG/nDppI8FIw0WArwL9CzESHp0op6OyTkv2ZM1dDvT1Up0DWto+RvZdOZ00KTVdtp++pg/kkeFgCCuJI7DE+hN3mRl9tec+0QjakHKNU0Ur50unZWQaTXVB8J898NeyBgyzBHAewbG101/y6sOHGkr/llWsfn9Hdqoh87aekdmbTsqFGzdl3NJDarE4rAT/ROsqauX5vIxGY5tLBPEIehHwIQtCg20xkZ2BlHq0CWf6buC7/NXA5nImKkFmbw6XX/dcVNNIcBIn3l3A1kF89pF4bDlG5CoY1BMRuQiMpI3qUl0qJByX6OJ1Zd62c3LhetYPYux7j4WXHmtVJc+ppKci42XGhlPq/PsPNFQp4oU1TzMvW8oZjppa2wngzDAS3KBisDoNapN1HdKa91+4oUZH9124IfvPx6jttrDYIE4apEBjNF8L9GuXK3ZHmRz2hs4MzIvfeCJK/bDH9AXDrAeMvCOQx4g8dhXIi/wujmhP+AyxnZa1W2ohOwAj2Frwr/97q1PgdgeAweWbaSqgMsx8MMdcx0lu7y+yBZ7vHCYV4o9JbEh9mbv1nArg31t5VD5be8JoN4vso9Gv3lNLElIyVDCfPdBDx87ddbJG45GujjbizMJKF1XrrLzSvaa8t/yo/HHgslPtLmANTEf5NzwrgMf0g1MWg/isOfFMpyd3wqCeiMjF+HuJPNOhmjzZvrpaoOrrf06rH3PvLDsiX244Lc90rCYDWofm2H7JMJ21dFE/+WzNCTWy1aV2aendqLwaSSzMH4PWjJpa2wngavDjtWudsuqkfZbY81rNyz93Q43m40fw2WtJ6rR03yV92nWjSrdH8/HXkbMZsDI5ftRjNH7j8Uj1WgxhhLtLndIqkEdbDfDN388cd9knGt9pLLCHkzXrYSAj4vHZO3ItZ64t5fX99fMSGdomVIa0DZPf912SrzZmZSGZon3nP/n7hNHztAgtId1uzY8Pc9HtIvE+4/VZCuqdZQpSfoN4dOBlTXkhck8M6onyqXPnztKkSROZNm2axXIdO3aU5557Tp544gmb1SU6Olrq1asne/fulUqVKtnseci5YO90bH004K4q8svui/LVxtNqNer3Vx5T5zHPdXDbqhLk72NycSXASBZWjLfXD2FrRk3tmTrtaPD+1CkXpE7YAgsw4rr/4u0gH+n7mJuP9xYnTaUSRbJG87EQX2gJqVc+yK4jmkih1laq337mmtHK5D5eHqo9YF48gnnsCuCKQZujKezMGLS/R1tWlooliqitPHPTrnop6d+qinSqWdptFpB09ilIOL5r6fTYstVUEF+3fNDt1emrMognMsSg3g0NHTpUFixYoL9csmRJadmypXz00UfSqFEju9bN1fz5558SEREhjz32mP665ORkGTdunCxevFhSUlKkR48eMnPmTClbNmuEzVyK5MSJE2X27NkSExMj7dq1k6+++kpq1qypbg8JCZHBgwerMt9++22hvDZyHphLj22R+resrEZosdcxRjgxovX1P2ekQ40Q+eu/qyZ/nCOd+b9LsWq/dXuxZtTUGVOnCwuCG8whxklbAfxMdILaSk+bn48V9jEXHadlt0b9EFBhb24E+Y0qFpMYK6ddW7ugWUp6huwKv5EVyB+PlDNRxvtHI5jMWjCstLSrEaLmRlPhsldmjOHidpagA+D+xhXEnTjbFCTDIB4nrAeSHYN4orzjf0I31bNnT5k3b546f/XqVRk/frz07t1bzp8/n+/HTE1NFV/fwjngIsjNyMgQb2/HbsKff/65DBs2TDw9b49yIaBfu3at/PLLLxIcHCyjRo2SBx98ULZu3Wr2cdDhgsdCZ0y1atXk7bffVp0BR44cEX//rF53PE/z5s3l448/Vh01RKbm3T7aorI82LSirDh0Rb5cf0rNOV3531Wz93HGxZXcJXX6TmHngRpliqkT2gXEJ6ep1fUR5GvB/o2ktGyLcHnLzFOb1Grb2EoPW5fVrxBscv50XrfXQgYJ0uk3HIuSbaejjeZMaynUaqX62qUdcmVyd2SPzBhnH422JUefgpRbEI+vdN1yQdKaQTxRvjh2ROSkEhNz9jZqvLy89EFYbmURCBYpUiTXsoGBed8yR+Pn5yflymWtpIy/Y8aMkQ4dOkhUVJSULp01ivPmm2/K0qVL5eLFi6rMgAEDZMKECfptP9555x35/fffVVD63nvvyblz5yQzM1ONJL/22mvyxx9/qJHoFi1ayGeffSaNGzc2W59t27bJyJEj5dixY9KgQQPVyfDAAw/Ivn37VIr7xo0bpUuXLmrrEdx26NAh+fvvv6VNmzby+uuvy08//SRxcXH650LmAcyfP19efvllVScN6ozHRseA4et49dVXVbB848YNuffee9WoeLFixfTv/YgRI+S3335T1+H15Qbv5fr162X69On662JjY2XhwoXq1LVrV3UdOlfq1q0rO3bskLvuuivH46CeSPHH6+7bt6+67rvvvlMj+6i3lgVQv359qVChgvrMnnrqqTy0AnJXWGitb5OK0qdRBbWA3mdrT7rs4kpkHcytxug3TtrxB1kd2kj+nnPX5diVOLVK+spDV9VJS4OvVyFrNF+bn48Mj5GLTO9jjsDjpW415WZqhhqNPxFhnGpbppifCuAxN75dzRA1RYQcT2FnxjjbaLSjdLQgS+eDBxsW6hSkyLhk2RF+O4g/YyaIv706fUm1fSIRuXhQf/36dXnhhRdk2bJlKth96KGHVLBUtGhRs/d59tln1Yjo5cuXVbm2bdvKhx9+KHXq1LFpXS3V6b777pMVK1boL5cpU0aSkkwv+tKpUycVzGqqVq2q5k5npwWn+ZWQkKCCzBo1akipUrd/tCN4RVCMQBFB9PDhw9V1b7zxhr7MqVOn5Ndff1XBLjos4JFHHlGdEX/99Zcaif7666/l7rvvlhMnTpgcQUYw3qdPH/Xe/PDDD6pzAIG4Keh8+OSTTyQsLExKlCih6oLnxwh2aGioGtHGCDbqZc1o9enTp1WAvHz5chXUP/roo/LBBx+ozgpAx8GmTZtURwU+M4y2Y/46OhzM2bJliwQEBKiAXbNnzx61hUq3bt3016E9VqlSRbZv324yqA8PD1fZFIb3wfvaunVrdR/D1P5WrVrJ5s2bGdRTnkdqq+ZxRXBnWVyJChZGxLFqPE4PNqukjl9Ll62UCg3ukkOXkbqPYP+GRCekyoELMeo0f1vWfRHXWdrHfJpBZxLKoiNAG43HPH6OxjuHwsyMcfTRaEfraJm/NVxWH4mQuuWK2TygZxBPZF9OE9RjlPjKlSuyZs0a9aMCqcbPPPOMCgLNQSoy7oeACZ0CGJG95557VJCkBaDuCsGr1vmAUejy5cur6wzTxDEybNihgNFpjIgbBvVIuceosTa6j0B2586dEhkZqbIBAEE4AuYlS5aozyw7fIb48YaRcWQxYMG3S5cuqU6E7N59913p3r27vt6YV46OB4ysAx4DbQTzyhGI5xUyDPA42sj8oEGDZN26dSqoR6cHHg8dH+icAHQi5LYgHTonMJpu+J4iOMcUheLFixuVRTncZop2ffY596bugw4YZDcQ5RXTWclaWJEcK023r1VW37GMOfhZAX5Wyv6hS7GSmYf+5g41Q1TqP/7yBz7lBRfEzHtHS+WSReTvoxGy/cx1tf1jXrd1zIsIBPG3VqbHKvWmgnh0zum3mKta0m0WLSSyB6cI6o8ePSqrVq2SXbt2qfRq+OKLL9TILgJGBDKmGAaQCEr/97//qRTws2fPSvXq1W1WXwSB5mTvTEDwa45hMAiod0FBKjsCYsDINBZqQ2CMgBwj3oCF3DCPG6PYeE3p6ekSFBRk9DgoqwX0cODAAVXWcMQfbt68qR7HlOPHj6sF+gynJWDE2RTt8wc8Hjp4sGicBlMDcF+0GWugfWgBPaCTQ/ts8DzovMDIuAZZALVr17b4mHjNhq+pMCBDwlzmB5EpTGelO4VO2colA9QJ0zpgyZ4L8tovB3O978PNK0kfN1vQjO4cF8TMGyxw2rlWadlwPEp+2nVext57O3PwjoJ4jMRHM4gnciROEdQjxRgjm4YBHVKREfT++++/an50bjCqi7nLWGSscmXzu7FiDjhOhqnhGgSQhnAZIxQY5cVJYzgP3pSCLmtYJi9QZ6SFI4Vd880336h0dvydPHmyes+R5aBlNyDdG0H+1KlT9c+Hx8F8fsPnj4+PVwEx5pJnh8/QVF216QOGt2nntfdWu4z3wPA2wzKGj6d9Ltplw9u1z9fwdnQGZK9b9ufO/jymHtsQAn90mBjejtF1dBDgerzfGqyQj9tMPRbS/QGZKoaj9bgPOqkM73Pt2jW1Er61beJO4LnwPuD74O4ZMI5AO05lP15Z8ta9teWFnw6YTWfF7ZkZ6ZJ5e+0yclN5bV/liuVt1L1UgLdVbZXcR17aWosqGGjIGmzgMcq0R5tXVEH9z7vOy11Vi6uFL7FmBRaftNQJgiB+59kb8m/4DdV5En4tyUQ6fTGVtdO6aglpUbWEBBcxDuId/budn/+XRIXZvqy5r1ME9Ugx1gIbDVY9R9BkLmVZgxFopIsjqMfIKlKzLa3QPmXKFJk0aZLJ23Df7HXAAnIYmUag5izQQDDqbthhgcAMnSRYyA3Xb9iwQXV+YBE8DeapI3jT7ofgGCvQGz4O3mN8Jti2DdMesjMsq0E5pLZjYTktZR/zwgGfG+6jjT6j00DLYECGAD5LrJuAefzaa0NGB/aFx/3Q6YD7ICDWFhRENoJhXUy9DtQf7wmuw/Mg6Mf6Bv369VO3Y+E9rBGAOfCmXhPUqlVLvRfYUUBLt8d1eCysq3D//fer606ePKnKNGzY0ORjIesBwTwWCdQ6YlAOHVrYxs7wPsiUaN++vdk62QLaPrIS/vnnH9WuyDFkP17lZlgtD/ntrKfEpN7+kRfsq5MHq2ZKxrk9svKcDSpJLtu+kHpf3NdLYtS/RlOBg06K+4pEHdkhK61LrCI3Y+2xjIxl6ESKeHnJjaR0GbZgr/764reO741LZXXlxqaKnIz1kFNxWaeoZOPvrYfopGKgSI0gndQM0klYkE4CvG+I6G5ISrjI1nBxWmxj5Kjty5rsW7sG9Vj0DAvXWWJtGnV2GG3GHGwEdUjVxwJo2DrMXFr02LFjZfTo0frLCI60kX08jrbyuxb4XbhwQc1NL+w06zuB14AgVmsoGDWeMWOG6pzA1mpIsccK9Fj1HoEkVpLHXwSiSLPUUvARgGNk1jAlH4EqVqRHsImF5hDEYqFC3B8BsWG2hebJJ59Uc9cxBx4r7iPARWcM4L3F4yOzAJAirz0f/iJ4RzZBxYoVVecAtnNDgImV9HE7phngvmhnWGgRgTDWBdDub+514PNE5wGuwwl1xPNgHj06mLDeAG5Hp0L2KQkaBNcYNT948KDaLlCr/8CBA9Uq+6gz7vvSSy+p90ybrw9YVwDviZaFgoUD8RrwuSDbBLsQYNrJ448/rm97+DwR1ON9N1cnW8D3ABkUHTt2dKrvgatCxxb+gWQ/XuXmPhF5I1Mnu8/dkMj4lDyN5JD7saZ9+VSNUBkgYjIDxEP+92Bj6VHfeK0Qojs9lpGx1Ycj5OaOrO+hIXTgzj3hJe2ql5RLMclql4uc6fTFpHVVLGpXQlqGlpCgbCPxzo5tjBy9fVkzSGfXoB5biA0dOtRiGYxMYjQ8+9xzjAhi8TttWzZzkDaOU82aNdWoKlKeseUXgiFTEOBpo8XZ4QMx/FAQGCPIRXCXff67I0OdV69erYJKLdDECuzYN13bZg0B+CuvvCIvvviiGsnu1auXCkQR2GqvVVuZOPtrRwD/1ltvqRXYMfqOzwgBH9LyTb1PGMXGrgbYMq5Zs2ZqxBpB6xNPPKECcsP3N/t7jUAX2QNDhgxRI/LoNMBr0+b0I6hGFgA6DObMmaMCZ7wGrLdg6XVkvw4dQsgawJZyeL/QdvFF0z5/U3A9FnT88ccf9aPyGP1///331YJ/yC7Ae4vV+tGJYfg4WGfAMCsBnR0I2tGJgSwBdBhgnQmtswPwHqJjA7smFCbUEe9D9u8H2Vd+Pg+U1hY/I7rT9tW7SSXx9vbigmZ0R/i/Jf8yMnXy3l/HLZbZevq6+ov+2/oVglU6PebEt6xWMkc6vatiGyNHbV/W3M9Dd6f7oRUCjNZj5HL37t1qRXvAHuU9e/ZUo8nmFsrLDgEUgnoEULl1JmgQuKFTACu0P/zwwzlG6rGSPkZOOUJZsBYtWqQCYkwHyG3dAUeG9HvsH4/t77CooJbSj5H0gu4IQqcVOmHQGVKY+D1wvJ5hdKxhIVH+SCFHaF8ILLigGVmLx7I7t/30NXl89o5cy73eo5YMvKuq2wTxGrYxcvT2pcWhiIdyy8J1ijn12OcbATy2OJs1a5Z6kzDXG/tzawE9tkDDKCy2V8Pq52fOnFELu2GRN8yJRvCPtGQEiHhzybHgc0NWBrIHkEKOkWlMlXDmgB6QpYDt8DClQNtVwBaio6PV1AlzGShERO6wjzkR3YaOtLyuku9uAT2Rq3GKoF4buUUgj8AdI5wPPfSQ2m5Ng0AfKcvaPHGMGGKxtWnTpqk541hoDCng27Zty7HoHjnGiDZS7vEXafpITcecclegLa5nS5hmgAUhiYiIiACZMQVZjogcl9ME9VjpHinwlvYZN5xJgBF8pDyQc0BAyqCUiIiIqGBgqkv5YH+5GptstFilxuPWGhcoR0TOzXlWdyMiIiIiojxPfcGilJB9FQvtMm7nGhdEzo9BfQFwgrUGiWyG7Z+IiMgxYZeJrwY2UyPyhnAZ13MXCiLX4DTp945IW8kQ8/idfUE3ovzS1rHgyrFERESOB4F793rluAsFkQtjUH8HvLy81B7rkZGR6jL2DNf2NycyBVvapaamqm3gCnpLO3uM0COgR/vH9wDfByIiInI83IWCyLUxqC+ALctAC+yJcguEb968qTI7XKUDCAG99j0gIiIiIqLCxaD+DiEwwxZs2CYP2+oRWYI28s8//6jtFV0hXR2vgSP0RERERET2w6C+gCCwYXBDuUEbSU9PF39/f5cI6omIiIiIyL6ce1IvERERERERkRtjUE9ERERERETkpBjUExERERERETkpzqnPw2rlgK274uLiOA+a7nihPLYlsiW2MbIlti8qLGxrZGtsY+To7Qv3NYxHLfHQ5aWUG7t48aJUrlzZ3tUgIiIiIiIiN3PhwgWpVKmSxTIM6nORmZkpx48fl3r16qk3NCgoyN5VIieGHjd0ErEtka2wjZEtsX1RYWFbI1tjGyNHb18I0+Pj46VChQri6Wl51jzT73OBN7BixYrqPD4QfumpILAtka2xjZEtsX1RYWFbI1tjGyNHbl/BwcF5KseF8oiIiIiIiIicFIN6IiIiIiIiIifFoD4P/Pz8ZOLEieov0Z1gWyJbYxsjW2L7osLCtka2xjZGrtS+uFAeERERERERkZPiSD0RERERERGRk2JQT0REREREROSkGNQTEREREREROSkG9URERHRHPDw85Pfff7d3NYiIiNySWwT1U6ZMkZYtW0qxYsWkTJky0q9fPzl+/LhRmeTkZHn++eelVKlSUrRoUXnooYckIiJCf/uBAwfk8ccfl8qVK0uRIkWkbt26Mn36dKPH2LJli7Rr1049BsrUqVNHPvvss0J7nWQfQ4cOVT9on3vuuRy3oU3hNpQhyg+2L7IHtCn8ryQqaFFRUTJixAipUqWKWhW6XLly0qNHD9m6dau9q0Yu9D/zgw8+MLoenY64nsgRY0hDOBZ6e3tLkyZNrKqrWwT1mzZtUm/2jh07ZM2aNZKWlib33HOPJCYm6su88sorsmzZMvnll19U+cuXL8uDDz6ov33Pnj3qw1y4cKEcPnxY3nrrLRk7dqx8+eWX+jKBgYEyatQo+eeff+To0aMyfvx4dfrmm28K/TVT4cIX9aeffpKbN28afcl/+OEH9cPlTqC9knuzZfsiIipM+MG7b98+WbBggZw4cUL+/PNP6dy5s1y7ds3eVSMX4e/vLx9++KHcuHHD3lUhJ7epkGJITUxMjAwePFjuvvtu6yurc0ORkZHYxk+3adMmdTkmJkbn4+Oj++WXX/Rljh49qsps377d7OOMHDlS16VLF4vP9cADD+gGDhxYgLUnRzNkyBBd3759dQ0aNNAtXLhQf/2iRYt0jRo1UrehDPz111+6du3a6YKDg3UlS5bU9erVS3fq1Cn9fcLDw1W7++mnn3QdO3bU+fn56ebNm2eX10Wu175wvHr++edzHA9x/Fu7dm0hvipylnYHoaGhus8++8zo9saNG+smTpyov4zj1tKlSwu9nuRcbty4odrKxo0bLZZ56qmndCEhIbpixYqp49b+/fv1t6Pdof3NmjVLV6lSJV2RIkV0jzzyiPotR4RjV+/evXV16tTRvf766/rrcXwyDHuWLFmiq1evns7X11cd4z755BP9bWPHjtW1atUqx2Pjf+6kSZMK4VWQu8aQ/fv3140fP15/nLOGW4zUZxcbG6v+lixZUt+Dgp6Xbt266csgdR4jYNu3b7f4ONpjmIKe6G3btkmnTp0KtP7kmJ588kmZN2+e/vLcuXNl2LBhRmXQszd69GjZvXu3rFu3Tjw9PeWBBx6QzMxMo3JjxoyRl156SWV8IC2RqCDa19NPP61G91NSUvT3Qc9xxYoVpWvXroX4aojIHSE1FSekQhsehww98sgjEhkZKX/99Zf6fdasWTM1anX9+nV9mVOnTsnPP/+sRsdWrVqlfm+NHDmyEF8JOTIvLy95//335YsvvpCLFy/muB3t6tFHH5XHHntMDh06JO+88468/fbbMn/+fHX7gAEDZOfOnXL69Gn9fTDCevDgQXniiScK9bWQ+8SQ8+bNkzNnzsjEiRPzVTe3C+rx4/bll19Wc98bNGigrrt69ar4+vpK8eLFjcqWLVtW3WYKgvXFixfLM888k+O2SpUqqXliLVq0UCkb+CFNrm/gwIFqXYVz586pE+bE4LrsaYdIyalRo4aaK4PADP9Qjhw5YlQObRTlqlWrJuXLly/kV0Ku2r60dLA//vhDfx/8iNHmIBIR2RLmieKYg9R7/ObCb7Fx48apYAlwjEMwhTRW/IaqWbOmfPLJJ6rskiVLjKYffffdd+o417FjRxW8YYqSud9s5H7QoY32YSpAmjp1quooQiBfq1Yt9T8Q02c//vhjdXv9+vWlcePGqhNcs2jRImndurX6/0ruKdOGMeTJkyfVgB4GWnCczA+3C+oRZP/333/q4J9fuH/fvn3VgQLzKrLbvHmzGimbNWuWTJs2TX788cc7rDU5g9KlS0uvXr3UDxb0tuF8SEiIURl8abFYRlhYmAQFBUnVqlXV9efPnzcqhx8zRAXdvjDPcNCgQSrYh71796rjGRfaI6LCgs5HzDnFXPqePXvKxo0b1Wg8jm1YUCohIUG/4JR2Cg8PNxo1xSgYMow0bdq0UT+4sy9gRe4N8+rRgYSsR0O4jMDMEC7jf2hGRoZ+tF4L6jHDCL/lcR25r+dtFEOizSEDZNKkSaqTKb/y1xXgpNALt3z5crWQHUbTNVh5NTU1VS1OYNjTgpULcZshjHihdw+9K1gEzxSMrkLDhg3VYyCtBz+0yT1SpNHOYMaMGTlu79Onj4SGhsrs2bOlQoUK6kcIevvQ/gxh0UUiW7QvZA5h9AIpiegcQNo97kNkDqZxZE2bv40LeNKdQAdj9+7d1QmjpTgu4UcuUuiRnYZAP7vsI2FEuUEWB6YwYlEyazuv8bv9zTffVJ3fWKT2woUL0r9/f5vVldw3hoyPj1eDwZhGpP3Gw+83/N/FqP3ff/+dpymSbhHU40154YUXZOnSpeofhRZ0a5o3by4+Pj5qDip6kAG9vRjdQu+v4XwavKlDhgyR9957L0/PjQ/F3Lwxcj0YdcCXG6nM2efCY2VftCsEXB06dNCnGhIVZvtCZyMyQVAOoxCmVl8lyp4lcuXKFf3luLg4NXJKVFDq1aun5tljxB4pq/ghq2UamYLfZxjtR+clYGVqdD7Vrl27EGtNzgBb26Ej27BtYEux7Fso4jJGSTEfHxC4YU0spN0jqEcHFFYwJ/eiK4QYEpmVmCppaObMmbJ+/Xo17Sj7c7p1UI90Cfx4xTxS7DOozXEIDg5W+wXi71NPPaUWmMKiBXhz8QHiw7jrrrv06RL4MPBDGuW0x8CXHz94tJEzpIRhgQRAbw7mgr344ot2e+1UuNAetDQv7R+DpkSJEiqlEFscYiQCX3jMnyEq7PaFUTH0BiMjBPMOiSzB/z6kRiMTBCMREyZMyNH+iPICnY9YCA9ZR40aNVK/yTBC9dFHH6mUVCw2hd9e2Asa1yHIQvC+YsUKdazSpqZhpB8/jvEbC51M+J2Fhc+yj4wRoSMbafOff/65/rpXX31V7T0+efJkNfqOBc3QwY1AyhDuhwwSdKZ/9tlndqg9uUMM6enpqZ+jr0EHEo5z2a+3SOcG8DJNnQy3Crt586baXqBEiRK6gIAAtRXdlStX9LdjawFTj4FtMDSff/65rn79+ur+QUFBuqZNm+pmzpypy8jIKPTXTPbZ+skUwy3H1qxZo6tbt67aqg5bo2BbH8OtoLQt7fbt21do9Sf3aV+a+Ph4dZzCMY/IlEGDBukeeughdT42NlZts4P/a5UrV9bNnz+fW9pRviQnJ+vGjBmja9asmdp6E8eh2rVrqy2ckpKSVJm4uDjdCy+8oKtQoYLaKgptbsCAAbrz58+r27WtnvD7CmX8/f11Dz/8sO769et2fnXkqP8z8dsKW9eZ2tIObaxKlSq6jz/+2OT2ivh/inaK/5vkfqSQYsjs8rOlncetChMRkZs4e/asVK9eXXbt2qXSXYlMTfXAKs+cnkGOBusUIVV///799q4KEZHDcLvV74mI3BUWN0PaFxZoQVoYA3rK7saNG2oxIMwdNNx3l4iIiByXW8ypJyKirIWAunTpouapGu75TKTBXGdkcGDOKeY4ExERkeNj+j0RERERERGRk2L6PREREREREZGTYlBPRERERERE5KQY1BMRERERERE5KQb1RERERERERE6KQT0RERERERGRk2JQT0REREREROSkGNQTERGRRUOHDhUPDw918vHxkbJly0r37t1l7ty5kpmZmefHmT9/vhQvXtymdSUiInI3DOqJiIgoVz179pQrV67I2bNn5a+//pIuXbrISy+9JL1795b09HR7V4+IiMhtMagnIiKiXPn5+Um5cuWkYsWK0qxZMxk3bpz88ccfKsDHCDxMnTpVGjZsKIGBgVK5cmUZOXKkJCQkqNs2btwow4YNk9jYWP2o/zvvvKNuS0lJkddee009Nu7bunVrVZ6IiIhyx6CeiIiI8qVr167SuHFj+e2339RlT09P+fzzz+Xw4cOyYMECWb9+vbzxxhvqtrZt28q0adMkKChIjfjjhEAeRo0aJdu3b5effvpJDh48KI888ojKDDh58qRdXx8REZEz8NDpdDp7V4KIiIgce059TEyM/P777zlue+yxx1QgfuTIkRy3LVmyRJ577jmJjo5WlzGi//LLL6vH0pw/f17CwsLU3woVKuiv79atm7Rq1Uref/99m70uIiIiV+Bt7woQERGR88LYAFLpYe3atTJlyhQ5duyYxMXFqbn2ycnJkpSUJAEBASbvf+jQIcnIyJBatWoZXY+U/FKlShXKayAiInJmDOqJiIgo344ePSrVqlVTC+hh0bwRI0bIe++9JyVLlpQtW7bIU089JampqWaDesy59/Lykj179qi/hooWLVpIr4KIiMh5MagnIiKifMGceYy0v/LKKyoox/Z2n376qZpbDz///LNReV9fXzUqb6hp06bqusjISOnQoUOh1p+IiMgVMKgnIiKiXCEd/urVqyoAj4iIkFWrVqlUe4zODx48WP777z9JS0uTL774Qvr06SNbt26VWbNmGT1G1apV1cj8unXr1AJ7GL1H2v2AAQPUY6BDAEF+VFSUKtOoUSPp1auX3V4zERGRM+Dq90RERJQrBPHly5dXgTlWpt+wYYNa6R7b2iFtHkE6trT78MMPpUGDBrJo0SIV9BvCCvhYOK9///5SunRp+eijj9T18+bNU0H9q6++KrVr15Z+/frJrl27pEqVKnZ6tURERM6Dq98TEREREREROSmO1BMRERERERE5KQb1RERERERERE6KQT0RERERERGRk2JQT0REREREROSkGNQTEREREREROSkG9UREREREREROikE9ERERERERkZNiUE9ERERERETkpBjUExERERERETkpBvVERERERERETopBPREREREREZE4p/8DBtocZ8nA++4AAAAASUVORK5CYII=", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "ndsi_timeseries = ndsi_masked.mean(dim=[\"y\", \"x\"])\n", + "ndsi_timeseries.plot(marker=\"o\", figsize=(12, 4))\n", + "\n", + "# Reference lines\n", + "plt.axhline(0.4, color=\"blue\", linestyle=\"--\", label=\"Snow threshold (0.4)\")\n", + "plt.axhline(0.2, color=\"gray\", linestyle=\"--\", label=\"Mixed snow (0.2)\")\n", + "plt.axhline(0.0, color=\"black\", linestyle=\"--\", label=\"Bare ground (0.0)\")\n", + "\n", + "plt.title(\"Mean NDSI Over Time\")\n", + "plt.xlabel(\"Date\")\n", + "plt.ylabel(\"NDSI\")\n", + "plt.legend()\n", + "plt.grid(True)\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "70f98e20", + "metadata": {}, + "source": [ + "> 🔍 **Note on Gaps in the NDSI Time Series**\n", + ">\n", + "> Breaks in the line plot occur when scenes are fully masked due to clouds or shadows — resulting in no valid NDSI values for that date.\n", + "> These time steps appear as `NaN` in the data and are rendered as gaps in the plot.\n", + ">\n", + "> This is expected behavior and helps highlight when cloud-free observations were not available for the selected area and time range." + ] + }, + { + "cell_type": "markdown", + "id": "a293cde1", + "metadata": {}, + "source": [ + "### ❄️ Interpreting NDSI Values\n", + "\n", + "The **Normalized Difference Snow Index (NDSI)** is a powerful tool for identifying snow cover in satellite imagery by comparing reflectance in the green and shortwave infrared (SWIR) bands. \n", + "Snow is typically **bright in the green band** and **absorptive in the SWIR band**, leading to high NDSI values.\n", + "\n", + "#### Thresholds and Uncertainty\n", + "\n", + "Interpreting NDSI values isn't always straightforward — thresholds can vary depending on:\n", + "- **land cover** (e.g., forested vs. open terrain),\n", + "- **illumination and sensor angle**,\n", + "- **scene conditions** (e.g., snow under cloud shadow or mixed with vegetation).\n", + "\n", + "That said, commonly used guidelines include:\n", + "\n", + "| NDSI Value Range | Interpretation |\n", + "|------------------|--------------------------------------------|\n", + "| > **0.4** | Likely **snow-covered** surface |\n", + "| 0.2 – 0.4 | Possibly **mixed snow** or patchy areas |\n", + "| 0.0 – 0.2 | Bare ground or dry soil |\n", + "| < **0.0** | Water, vegetation, or cloud shadow |\n", + "\n", + "These are not hard rules — think of them as **decision aids** rather than strict classifiers.\n", + "\n", + "In this notebook, we include horizontal reference lines on the NDSI time series plot to help guide interpretation. \n", + "However, if your AOI is forested, coastal, or mountainous, it may be worth adjusting the snow threshold or inspecting example scenes visually for validation." + ] + }, + { + "cell_type": "markdown", + "id": "a5464d8d", + "metadata": {}, + "source": [ + "### 📤 Exporting Results for Use in Power BI\n", + "\n", + "To make the NDSI results available in external tools like Power BI, we can export a simple table summarizing:\n", + "\n", + "- the **date** of each satellite scene,\n", + "- the **mean NDSI** value for that date,\n", + "- a **classified result** indicating snow presence (e.g., `\"snow\"`, `\"mixed\"`, or `\"no snow\"`).\n", + "\n", + "This table can be saved as a CSV file and uploaded to Power BI or other data analysis platforms for further visualization and reporting." + ] + }, + { + "cell_type": "code", + "execution_count": 139, + "id": "1131555f", + "metadata": {}, + "outputs": [], + "source": [ + "# Create a DataFrame from the NDSI time series\n", + "ndsi_df = ndsi_timeseries.to_dataframe(name=\"ndsi\").reset_index()\n", + "\n", + "# Add a simple classification based on threshold\n", + "def classify_ndsi(value):\n", + " if pd.isna(value):\n", + " return \"no data\"\n", + " elif value > 0.4:\n", + " return \"snow\"\n", + " elif value > 0.2:\n", + " return \"mixed\"\n", + " else:\n", + " return \"no snow\"\n", + "\n", + "ndsi_df[\"classification\"] = ndsi_df[\"ndsi\"].apply(classify_ndsi)\n", + "\n", + "# Export to CSV\n", + "ndsi_df.to_csv(\"ndsi_summary.csv\", index=False)" + ] + }, + { + "cell_type": "markdown", + "id": "5077e211", + "metadata": {}, + "source": [ + "### 📊 Visualizing Results in Power BI\n", + "\n", + "To explore the results interactively within Microsoft’s ecosystem, you can upload the `ndsi_summary.csv` file to [Power BI](https://app.powerbi.com/). \n", + "This allows you to create custom dashboards and visualizations—such as time series plots or spatial summaries—directly from your NDSI data, all within the familiar Microsoft suite." + ] + }, + { + "cell_type": "markdown", + "id": "a59fd87a", + "metadata": {}, + "source": [ + "## ✅ Next Steps\n", + "\n", + "This notebook introduced the core concepts of time-series monitoring using Sentinel-2 imagery and NDSI to detect snow cover over a user-defined area and date range.\n", + "\n", + "To continue exploring site monitoring workflows, check out the next notebook:\n", + "\n", + "👉 **[site-monitoring-hls.ipynb](./site-monitoring-hls.ipynb)** \n", + "This notebook uses HLS (Harmonized Landsat and Sentinel) data to demonstrate cross-sensor monitoring and highlight additional techniques like combining indices and refining temporal analysis.\n", + "\n", + "### 💡 Other ideas to try:\n", + "- Use NDVI to monitor vegetation changes before/after snow cover\n", + "- Compare results across different years to study seasonal patterns\n", + "- Apply thresholds to detect persistent snow vs. transient snowfall\n", + "- Export your xarray dataset for visualization in other GIS tools\n", + "\n", + "---\n", + "\n", + "## 📎 Supporting Materials\n", + "\n", + "- [Sentinel-2 STAC collection on Planetary Computer](https://planetarycomputer.microsoft.com/dataset/sentinel-2-l2a)\n", + "- [NDSI overview – USGS](https://www.usgs.gov/landsat-missions/normalized-difference-snow-index)\n", + "- [STAC specification](https://stacspec.org/)\n", + "- [ODC-STAC documentation](https://odc-stac.readthedocs.io/)\n", + "- [Xarray documentation](https://docs.xarray.dev/)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": ".venv", + "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.10.8" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}