{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "44ef49d8",
   "metadata": {},
   "source": [
    "# Calibration and error propagation\n",
    "\n",
    "Previously, we compared `phoptic` to [`Source-Extractor`](https://sextractor.readthedocs.io/en/latest/index.html) and showed that the two programs produce [very similar results](sextractor_comparison.ipynb). However, we neglected image calibrations in that example to keep things simple.\n",
    "\n",
    "When performing image calibrations, `phoptic` will automatically propagate the errors; in this notebook, I will prove that `phoptic` does this correctly. Additionally, I will also verify that calibrations applied by `phoptic` are correct. \n",
    "\n",
    "## Generating data\n",
    "\n",
    "First, we need to generate some data that contain only noise. We will include typical values for the bias and dark current, and simulate long exposures to allow for a significant dark noise contribution. We will also generate flat-field images, which introduce some additional noise. These images will be generated as if they were produced by OPTICAM-MX, but we'll only simulate data from a single camera for simplicity.\n",
    "\n",
    "Let's define the properties of these generated data: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "fdc354ea",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:13.031913Z",
     "iopub.status.busy": "2026-06-16T12:21:13.031832Z",
     "iopub.status.idle": "2026-06-16T12:21:13.034260Z",
     "shell.execute_reply": "2026-06-16T12:21:13.033874Z"
    }
   },
   "outputs": [],
   "source": [
    "X, Y = 512, 512  # image dimensions\n",
    "\n",
    "gain = 1.  # gain of OPTICAM-MX\n",
    "read_noise = 1.1  # nominal read noise of OPTICAM-MX\n",
    "dark_curr = 0.1  # nominal dark current of OPTICAM-MX\n",
    "texp = 30.  # exposure time of science images (in seconds)\n",
    "flat_texp = 3.  # exposure time of flat-field images (in seconds)\n",
    "fltr = 'g'  # image filter\n",
    "\n",
    "binning = f'{2048 // X}x{2048 // Y}'  # binning mode\n",
    "\n",
    "bias_value = 400.  # typical bias according to Handbook of CCD Astronomy, Howell\n",
    "\n",
    "flat_level = 40000.  # mean flux of flats\n",
    "\n",
    "# dark noise contributions\n",
    "dark_signal = dark_curr * texp\n",
    "flat_dark_signal = dark_curr * flat_texp\n",
    "\n",
    "n_calibration_images = 30  # number of images to generate for each calibration step\n",
    "\n",
    "science_flux = 1000.  # mean flux in science images (before noise)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e1eec252",
   "metadata": {},
   "source": [
    "We can compute the analytical variances for each calibration:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ecb6071d",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:13.035158Z",
     "iopub.status.busy": "2026-06-16T12:21:13.035082Z",
     "iopub.status.idle": "2026-06-16T12:21:13.080750Z",
     "shell.execute_reply": "2026-06-16T12:21:13.080201Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "true_bias_var = read_noise**2 / n_calibration_images\n",
    "\n",
    "true_dark_var = np.pi / (2 * n_calibration_images) * (dark_signal + read_noise**2 + true_bias_var)\n",
    "\n",
    "true_flat_var = np.pi / (2 * n_calibration_images) * (flat_level + read_noise**2 + true_bias_var + flat_dark_signal)\n",
    "effective_flat_var = true_flat_var * science_flux**2 / flat_level**2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "70542fa0",
   "metadata": {},
   "source": [
    "I have included $\\frac{\\pi}{2 N}$ factors in the dark and flat variances since `phoptic` computes master darks and master flats using the median (to remove outliers due to, e.g., cosmic rays), while the master bias is computed using the mean. Note that the dark and flat-field variances include a bias contribution, and the flat-field variance further includes a (small) dark-noise contribution. The flat variance also has an additional term because it acts multiplicatively, whereas the bias and dark corrections are additive.\n",
    "\n",
    "Now let's create the directories for storing our simulated data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "af934725",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:13.086758Z",
     "iopub.status.busy": "2026-06-16T12:21:13.086627Z",
     "iopub.status.idle": "2026-06-16T12:21:13.088936Z",
     "shell.execute_reply": "2026-06-16T12:21:13.088719Z"
    }
   },
   "outputs": [],
   "source": [
    "from pathlib import Path\n",
    "\n",
    "out_dir = Path('out')\n",
    "\n",
    "if not out_dir.joinpath('biases').is_dir():\n",
    "    out_dir.joinpath('biases').mkdir(parents=True)\n",
    "\n",
    "if not out_dir.joinpath('darks').is_dir():\n",
    "    out_dir.joinpath('darks').mkdir(parents=True)\n",
    "\n",
    "if not out_dir.joinpath('flats').is_dir():\n",
    "    out_dir.joinpath('flats').mkdir(parents=True)\n",
    "\n",
    "if not out_dir.joinpath('flat_darks').is_dir():\n",
    "    out_dir.joinpath('flat_darks').mkdir(parents=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "97a40484",
   "metadata": {},
   "source": [
    "We'll use `numpy`'s `default_rng()` Generator to generate random noise. For reproducibility, I'll pass a fixed seed:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "33f7ef94",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:13.089710Z",
     "iopub.status.busy": "2026-06-16T12:21:13.089627Z",
     "iopub.status.idle": "2026-06-16T12:21:13.099811Z",
     "shell.execute_reply": "2026-06-16T12:21:13.099455Z"
    }
   },
   "outputs": [],
   "source": [
    "rng = np.random.default_rng(42)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8aeb4644",
   "metadata": {},
   "source": [
    "Let's now generate our calibration images.\n",
    "\n",
    "### Biases\n",
    "\n",
    "Bias images should be taken with the telescope shutter closed and have an exposure time of 0.0 seconds:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b5186356",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:13.100596Z",
     "iopub.status.busy": "2026-06-16T12:21:13.100494Z",
     "iopub.status.idle": "2026-06-16T12:21:16.972234Z",
     "shell.execute_reply": "2026-06-16T12:21:16.971821Z"
    }
   },
   "outputs": [],
   "source": [
    "from astropy.io import fits\n",
    "\n",
    "for i in range(n_calibration_images):\n",
    "    bias = rng.normal(bias_value, read_noise, size=(X, Y))\n",
    "    \n",
    "    hdu = fits.PrimaryHDU(bias)\n",
    "    \n",
    "    hdu.header['EXPOSURE'] = 0.0\n",
    "    hdu.header[\"BINNING\"] = binning\n",
    "    hdu.header[\"GAIN\"] = gain\n",
    "    hdu.header['FILTER'] = fltr\n",
    "    hdu.header['INSTRUME'] = 'OPTICAM'\n",
    "    \n",
    "    path = out_dir / 'biases' / f'image_{i}.fits.gz'\n",
    "    if not path.is_file():\n",
    "        hdu.writeto(path, overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6a7bbb48",
   "metadata": {},
   "source": [
    "Let's pass these bias images to a `BiasCorrector`: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "08f054a1",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:16.973803Z",
     "iopub.status.busy": "2026-06-16T12:21:16.973681Z",
     "iopub.status.idle": "2026-06-16T12:21:19.359852Z",
     "shell.execute_reply": "2026-06-16T12:21:19.359525Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:   0%|                                                                                                                              |[00:00<?]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  30%|████████████████████████████████████▌                                                                                     |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  60%|█████████████████████████████████████████████████████████████████████████▏                                                |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  90%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████▊            |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] 30 1 bias images.\n"
     ]
    }
   ],
   "source": [
    "import phoptic\n",
    "\n",
    "bias_corr = phoptic.BiasCorrector(\n",
    "    out_directory=out_dir / 'reduced',\n",
    "    data_directory=out_dir / 'biases',\n",
    "    instrument=phoptic.OPTICAM_MX(),\n",
    "    )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1a70369c",
   "metadata": {},
   "source": [
    "### Darks\n",
    "\n",
    "Dark images should again be taken with the telescope shutter closed, but now have an exposure time equal to that of the science images. Like all images, dark images also include a bias offset. We'll generate one set of darks for the science images, and another set of darks, using a different exposure time, for the flat-field images:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "4c644f00",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:19.360945Z",
     "iopub.status.busy": "2026-06-16T12:21:19.360803Z",
     "iopub.status.idle": "2026-06-16T12:21:27.390966Z",
     "shell.execute_reply": "2026-06-16T12:21:27.390585Z"
    }
   },
   "outputs": [],
   "source": [
    "\n",
    "for i in range(n_calibration_images):\n",
    "    dark = rng.poisson(dark_signal, size=(X, Y)).astype(np.float64)  # dark noise\n",
    "    dark += rng.normal(bias_value, read_noise, size=(X, Y))  # add bias\n",
    "    \n",
    "    hdu = fits.PrimaryHDU(dark)\n",
    "    \n",
    "    hdu.header['EXPOSURE'] = texp\n",
    "    hdu.header[\"BINNING\"] = binning\n",
    "    hdu.header[\"GAIN\"] = gain\n",
    "    hdu.header['FILTER'] = fltr\n",
    "    hdu.header['INSTRUME'] = 'OPTICAM'\n",
    "    \n",
    "    path = out_dir / 'darks' / f'image_{i}.fits.gz'\n",
    "    if not path.is_file():\n",
    "        hdu.writeto(path, overwrite=True)\n",
    "    \n",
    "    dark = rng.poisson(flat_dark_signal, size=(X, Y)).astype(np.float64)  # flat-field dark noise\n",
    "    dark += rng.normal(bias_value, read_noise, size=(X, Y))  # add bias\n",
    "    \n",
    "    hdu = fits.PrimaryHDU(dark)\n",
    "    \n",
    "    hdu.header['EXPOSURE'] = flat_texp\n",
    "    hdu.header[\"BINNING\"] = binning\n",
    "    hdu.header[\"GAIN\"] = gain\n",
    "    hdu.header['FILTER'] = fltr\n",
    "    hdu.header['INSTRUME'] = 'OPTICAM'\n",
    "    \n",
    "    path = out_dir / 'flat_darks' / f'image_{i}.fits.gz'\n",
    "    if not path.is_file():\n",
    "        hdu.writeto(path, overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d9b6728e",
   "metadata": {},
   "source": [
    "Let's pass these dark images to `DarkNoiseCorrector`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "6cecdbf6",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:27.392551Z",
     "iopub.status.busy": "2026-06-16T12:21:27.392470Z",
     "iopub.status.idle": "2026-06-16T12:21:28.414819Z",
     "shell.execute_reply": "2026-06-16T12:21:28.414397Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:   0%|                                                                                                                              |[00:00<?]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  33%|████████████████████████████████████████▋                                                                                 |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  67%|█████████████████████████████████████████████████████████████████████████████████▎                                        |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] 30 1:g dark images.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:   0%|                                                                                                                              |[00:00<?]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  30%|████████████████████████████████████▌                                                                                     |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  60%|█████████████████████████████████████████████████████████████████████████▏                                                |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  90%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████▊            |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] 30 1:g dark images.\n"
     ]
    }
   ],
   "source": [
    "# science image corrector\n",
    "dark_corr = phoptic.DarkNoiseCorrector(\n",
    "    out_directory=out_dir / 'reduced',\n",
    "    data_directory=out_dir / 'darks',\n",
    "    instrument=phoptic.OPTICAM_MX(),\n",
    "    bias_corrector=bias_corr,\n",
    "    )\n",
    "\n",
    "# flat-field image corrector\n",
    "flat_dark_corr = phoptic.DarkNoiseCorrector(\n",
    "    out_directory=out_dir / 'reduced_flat_dark',\n",
    "    data_directory=out_dir / 'flat_darks',\n",
    "    instrument=phoptic.OPTICAM_MX(),\n",
    "    bias_corrector=bias_corr,\n",
    "    )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7bed2f55",
   "metadata": {},
   "source": [
    "Note that I have also passed our previously defined `BiasCorrector` instance to both dark noise correctors to automatically bias-correct the dark images.\n",
    "\n",
    "### Flats\n",
    "\n",
    "Flat-field images should be taken using a bright, uniform field with an exposure time that prevents the detector from saturating. Flat-field images contain a bias offset **and** a dark noise contribution: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "1889d203",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:28.415669Z",
     "iopub.status.busy": "2026-06-16T12:21:28.415595Z",
     "iopub.status.idle": "2026-06-16T12:21:33.188929Z",
     "shell.execute_reply": "2026-06-16T12:21:33.188286Z"
    }
   },
   "outputs": [],
   "source": [
    "for i in range(n_calibration_images):\n",
    "    flat = rng.poisson(flat_level, size=(X, Y)).astype(np.float64)\n",
    "    flat += rng.normal(bias_value, read_noise, size=(X, Y))  # add bias\n",
    "    flat += rng.poisson(flat_dark_signal, size=(X, Y)).astype(np.float64)  # add dark noise\n",
    "    \n",
    "    hdu = fits.PrimaryHDU(flat)\n",
    "    \n",
    "    hdu.header['EXPOSURE'] = flat_texp\n",
    "    hdu.header[\"BINNING\"] = binning\n",
    "    hdu.header[\"GAIN\"] = gain\n",
    "    hdu.header['FILTER'] = fltr\n",
    "    hdu.header['INSTRUME'] = 'OPTICAM'\n",
    "    \n",
    "    path = out_dir / 'flats' / f'image_{i}.fits.gz'\n",
    "    if not path.is_file():\n",
    "        hdu.writeto(path, overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bbe97b5f",
   "metadata": {},
   "source": [
    "Let's pass these flats to a `FlatFieldCorrector`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "e2d5b953",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:33.189880Z",
     "iopub.status.busy": "2026-06-16T12:21:33.189798Z",
     "iopub.status.idle": "2026-06-16T12:21:33.689035Z",
     "shell.execute_reply": "2026-06-16T12:21:33.688445Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:   0%|                                                                                                                              |[00:00<?]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  33%|████████████████████████████████████████▋                                                                                 |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory:  67%|█████████████████████████████████████████████████████████████████████████████████▎                                        |[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\r",
      "[OPTICAM] Scanning data directory: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████|[00:00<00:00]"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] 30 1:g flat-field images.\n"
     ]
    }
   ],
   "source": [
    "flat_corr = phoptic.FlatFieldCorrector(\n",
    "    out_directory=out_dir / 'reduced',\n",
    "    data_directory=out_dir / 'flats',\n",
    "    instrument=phoptic.OPTICAM_MX(),\n",
    "    bias_corrector=bias_corr,\n",
    "    dark_corrector=flat_dark_corr,\n",
    "    )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0d65e95a",
   "metadata": {},
   "source": [
    "Note that I have passed our previously define `BiasCorrector` instance **and** the `DarkNoiseCorrector` instance that I created to subtract the dark noise from the flat-field images.\n",
    "\n",
    "## Propagation\n",
    "\n",
    "Now that we have defined our calibrators, we can use them to calibrate some synthetic images and evaluate how accurate `phoptic`'s error propagation is. `phoptic` uses a helper function, `phoptic.utils.helpers.propagate_errors()` to propagate errors, so this is what we will use here:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "d99176d4",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:33.690084Z",
     "iopub.status.busy": "2026-06-16T12:21:33.690007Z",
     "iopub.status.idle": "2026-06-16T12:21:37.802374Z",
     "shell.execute_reply": "2026-06-16T12:21:37.801887Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] 1 master bias image not found. Attempting to create.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master bias image(s) created.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master bias image(s) saved to out/reduced/master_bias.fits.gz.\n",
      "[OPTICAM] 1:g master dark image not found. Attempting to create.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master dark image(s) created.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master dark image(s) saved to out/reduced/master_darks.fits.gz.\n",
      "[OPTICAM] 1:g master flat-field image not found. Attempting to create.\n",
      "[OPTICAM] 1:g master dark image not found. Attempting to create.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master dark image(s) created.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master dark image(s) saved to out/reduced_flat_dark/master_darks.fits.gz.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master flat-field image(s) created.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[OPTICAM] Master flat-field image(s) saved to out/reduced/master_flats.fits.gz.\n"
     ]
    }
   ],
   "source": [
    "from phoptic.utils.helpers import combine_variances\n",
    "\n",
    "calibrated_images = []\n",
    "bias_vars = []\n",
    "dark_vars = []\n",
    "flat_vars = []\n",
    "errs = []\n",
    "\n",
    "noiseless_image = np.ones((X, Y)) * science_flux\n",
    "\n",
    "for i in range(100):\n",
    "    image = rng.poisson(noiseless_image).astype(np.float64)\n",
    "    image += rng.normal(bias_value, read_noise, size=(X, Y))  # add bias\n",
    "    image += rng.poisson(dark_signal)  # add dark noise\n",
    "    \n",
    "    # apply bias correction\n",
    "    bias_corrected_image, bias_var = bias_corr.correct(image, camera='1')\n",
    "    bias_vars.append(bias_var)\n",
    "    \n",
    "    # apply dark noise correction\n",
    "    dark_corrected_image, dark_var = dark_corr.correct(bias_corrected_image, key='1:g')\n",
    "    dark_vars.append(dark_var)\n",
    "    \n",
    "    # apply flat correction\n",
    "    flat_corrected_image, flat_var = flat_corr.correct(dark_corrected_image, key='1:g')\n",
    "    flat_vars.append(flat_var)\n",
    "    calibrated_images.append(flat_corrected_image)\n",
    "    \n",
    "    # error propagation function used by phoptic\n",
    "    var = combine_variances(\n",
    "        flat_corrected_image,\n",
    "        bias_var=bias_var,\n",
    "        dark_var=dark_var,\n",
    "        flat_var=flat_var,\n",
    "        background_rms=0.,\n",
    "        read_noise=read_noise,\n",
    "        )\n",
    "    errs.append(np.sqrt(var))\n",
    "\n",
    "calibrated_images = np.asarray(calibrated_images)\n",
    "bias_vars = np.asarray(bias_vars)\n",
    "dark_vars = np.asarray(dark_vars)\n",
    "flat_vars = np.asarray(flat_vars)\n",
    "errs = np.asarray(errs)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d0aaf722",
   "metadata": {},
   "source": [
    "If `phoptic`'s error propagation is correct, then the average propagated error should be close to the empirical error. We can compute the empirical error by taking the standard deviation of the calibrated images. Let's check that the average ratio between the propagated and empirical errors is close to 1:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "70eb3028",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.803698Z",
     "iopub.status.busy": "2026-06-16T12:21:37.803624Z",
     "iopub.status.idle": "2026-06-16T12:21:37.837841Z",
     "shell.execute_reply": "2026-06-16T12:21:37.837491Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9930406244538639\n"
     ]
    }
   ],
   "source": [
    "empirical_err = np.std(calibrated_images, axis=0)\n",
    "\n",
    "measured_err = np.mean(errs, axis=0)\n",
    "\n",
    "print(np.mean(empirical_err / measured_err))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "25a20e24",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.838709Z",
     "iopub.status.busy": "2026-06-16T12:21:37.838639Z",
     "iopub.status.idle": "2026-06-16T12:21:37.840183Z",
     "shell.execute_reply": "2026-06-16T12:21:37.839915Z"
    },
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "# check errors are within 1%\n",
    "assert np.isclose(np.mean(empirical_err / measured_err), 1., rtol=0.01, atol=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e1fae7c",
   "metadata": {},
   "source": [
    "That's reassuring! We can also compare each calibrator's variance to its analytical value:\n",
    "\n",
    "#### Biases:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "f3a44c8f",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.841003Z",
     "iopub.status.busy": "2026-06-16T12:21:37.840938Z",
     "iopub.status.idle": "2026-06-16T12:21:37.858390Z",
     "shell.execute_reply": "2026-06-16T12:21:37.857996Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.999480261881968\n"
     ]
    }
   ],
   "source": [
    "print(np.mean(bias_vars / true_bias_var))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "b65db048",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.859573Z",
     "iopub.status.busy": "2026-06-16T12:21:37.859440Z",
     "iopub.status.idle": "2026-06-16T12:21:37.877156Z",
     "shell.execute_reply": "2026-06-16T12:21:37.876782Z"
    },
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "assert np.isclose(np.mean(bias_vars / true_bias_var), 1., rtol=0.01, atol=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "56a144f2",
   "metadata": {},
   "source": [
    "#### Darks:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "cc5db80e",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.878377Z",
     "iopub.status.busy": "2026-06-16T12:21:37.878309Z",
     "iopub.status.idle": "2026-06-16T12:21:37.895414Z",
     "shell.execute_reply": "2026-06-16T12:21:37.895065Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9997573182415519\n"
     ]
    }
   ],
   "source": [
    "print(np.mean(dark_vars / true_dark_var))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "c9f8070d",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.896199Z",
     "iopub.status.busy": "2026-06-16T12:21:37.896132Z",
     "iopub.status.idle": "2026-06-16T12:21:37.913616Z",
     "shell.execute_reply": "2026-06-16T12:21:37.913272Z"
    },
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "assert np.isclose(np.mean(dark_vars / true_dark_var), 1., rtol=0.01, atol=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6d9d61e0",
   "metadata": {},
   "source": [
    "#### Flats:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "367b6ffb",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.914776Z",
     "iopub.status.busy": "2026-06-16T12:21:37.914707Z",
     "iopub.status.idle": "2026-06-16T12:21:37.932150Z",
     "shell.execute_reply": "2026-06-16T12:21:37.931820Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.0004358917515763\n"
     ]
    }
   ],
   "source": [
    "print(np.mean(flat_vars / effective_flat_var))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "283832e7",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.933124Z",
     "iopub.status.busy": "2026-06-16T12:21:37.933055Z",
     "iopub.status.idle": "2026-06-16T12:21:37.950541Z",
     "shell.execute_reply": "2026-06-16T12:21:37.950064Z"
    },
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "assert np.isclose(np.mean(flat_vars / effective_flat_var), 1., rtol=0.01, atol=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8280df48",
   "metadata": {},
   "source": [
    "All variances are well within 1% of their analytical values!\n",
    "\n",
    "Let's compare the propagated errors to their analytical/empirical values graphically:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "1c92eae4",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:37.951763Z",
     "iopub.status.busy": "2026-06-16T12:21:37.951658Z",
     "iopub.status.idle": "2026-06-16T12:21:38.253541Z",
     "shell.execute_reply": "2026-06-16T12:21:38.253243Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABCoAAAQuCAYAAAAaxNzcAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlcelbwAAAAlwSFlzAAAPYQAAD2EBqD+naQAA17hJREFUeJzs3Xd809X+x/F3S2lpoAMaNg1bCIhAQYaAXEDZuBGVgBf1igoquPG6RXFclQuOq9dNnbhRUQRxIDJkKCOAslI2YXSQUjq+vz/49XtJBxRIm2/S1/Px6OM2J59vcpJbk8P7e77nRBiGYQgAAAAAAMACIoPdAQAAAAAAgEIEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJZBUAEAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAMggoAFWrXrl2aOnWqNm/eHOyuAAAAlGrv3r2aOnWqNmzYEOyuAJVOVLA7ACC0HThwQG+99ZZfW40aNeRwOPS3v/1N0dHRfvdt2bJFEydOVIsWLdS0adOK7CoAAAhBOTk5eumll8zbVapUUXx8vJKTk9WxY0fVrFmzXJ53+/btmjhxourVq6czzjijXJ4DQMkIKgCclt27d2vixInq0qWLunfvLklKT0/X5MmTdeTIEb355pvq37+/WV+/fn3deuutatasWbC6DAAAQkh2drYmTpyotm3b6rzzzpNhGMrIyNCaNWu0cuVKDR06VP/6178YWwBhhKACQECcf/75mjx5snk7OztbHTp00OjRo7Vjxw5FRh690qxx48aaOnVqkHoJAABCVbdu3YqNIdasWaMRI0aoU6dOWrhwoZxOZ3A6ByCgCCoAlIvY2Fide+65evXVV3XgwAElJSVJOrpGxfvvv68LL7zQvPQjIyNDr7/+unlsTEyMmjVrpj59+hS7dOTw4cOaN2+etm/frtq1a6tr165q0KBBxb0wAABgGW3bttWXX36pVq1aady4cfr+++8lSfn5+Zo+fbpZFx0dLYfDob59+8pms5nte/fu1TvvvKPBgwerefPmmjNnjjZu3KiLLrqo1Odcv369Zs+erXbt2qlfv37l9tqAyoygAkC5KCgo0IoVK1SvXj3VqlXLbC9pjYq8vDxt2bLFrDlw4IAee+wxRUZG6vvvv1eLFi0kSStWrNCAAQNUv359de3aVQcPHtSdd96pcePGaeLEiRX6+gAAgDU0adJEgwYN0ueff66dO3eqfv36MgzDb2yRnp6uZ599Vunp6frmm2/UqVMnSf9bhyImJkbXX3+9ateurZ07d+rMM89UYmJiseeaN2+eLrvsMp1zzjm69tprK+gVApUPQQWAgFiyZIk5HTMjI0Nz5sxRRkaG3nnnHUVERBz32Fq1ahWbypmTk6NevXrppptu0pw5cyRJ//znP9W4cWMtXrzYvJTkyJEjWrBgQcBfDwAACB0dO3bU559/rj/++EP169dXVFRUsbFFfn6+hg4dqjFjxuiPP/7wu+/ZZ5/V559/rjZt2qigoEBer1c7duzwq3n11Vd144036oYbbtDUqVNVpUqV8n5ZQKVFUAEgINLT080zFz6fT4cOHVK1atV05MiRMh1vGIaWLVumVatWKSMjQ4ZhKDExUQsXLpRhGIqIiNCBAweUn5+v3NxcxcTESDo6lbNv377l9bIAAEAISEhIkHR0PHKsP/74QytWrFB6eroKCgpUvXp1c6wRHx9v1vXt21dt2rSRJEVGRqpOnTpmUFFQUKA777xTzz33nKZOnarx48dX0KsCKi+CCgABUXQxTUm6/vrrNXToUK1cuVJnnnlmqcfu2bNHQ4cO1YYNG9SnTx81aNBAVatWVVZWlg4dOiSfz6fq1avrlltu0ejRo9W0aVMNGzZMvXv31vnnn6/atWuX98sDAAAWlpWVJUmKi4uTJGVmZuqSSy7RwoUL1bdvXyUnJys6OloHDx6UdHRtimODirPOOqvUx7777ru1fft2ffHFFxo6dGj5vQgApshgdwBA+Pr73/+u/Px8vfvuu8ete+CBB7Ru3Tr9/vvv+vTTT/XCCy9o6tSpOvvssyUdnW0hSVdeeaXWrl2r8ePHa+vWrRo7dqwaNWqkJ598stxfCwAAsK7Vq1dLkjkr4plnntH8+fO1aNEizZo1Sy+++KKmTp1qLn5ZOLYoVLjod0l69OihqlWr6vXXX1dOTk45vQIAxyKoAFBuCnfsKDx7UZoVK1aoffv2aty4cbH2olq2bKl7771X33zzjfbs2aNhw4bpnnvu0a5duwLWbwAAEDq8Xq++/PJLnX322eZYYsWKFWratKnatWvnV1vS2OJELrroIs2aNUtz5szR0KFDdejQoYD0G0DpCCoAlJtPP/1UktSrV6/j1jVp0kR//fWXsrOzzbbZs2dr1apVfnVLly71ux0bG2uu2n348OFAdBkAAIQQr9eryy67TPn5+Zo2bZrZ3qRJE+3YsUP79u0z25YuXWpuX3qy+vfvr++++07Lli3TeeedpwMHDpx23wGUjjUqAATEsbt++Hw+LV26VF988YWuvfZajRgx4rjH3nvvvfr666917rnn6sILL9TmzZu1fv163XjjjZoyZYpZN3nyZKWlpalnz55q2LChNm/erLffflvjx49XkyZNyvHVAQCAYFuzZo051sjMzNTq1av11VdfqVWrVpo/f766detm1t52221699131atXL1155ZXasWOHlixZojvuuEOTJk06pefv3r27fvzxR/Xv31+9e/fWnDlzVK9evUC8NABFRBhFL9ACgJOwZ88ePf74435t1apVU3Jysvr06WNeK1po69ateu6553T99df73bd9+3Z9+OGH2rNnj9q2bavhw4frp59+0ldffaWnnnrKvIxk1apV+uGHH7R9+3bVrVtX/fv3V9u2bcv/hQIAgKDIzs72CxeqVKmiGjVqyOFwqGvXrqUu2L1v3z6999572r59u1q0aKERI0Zo9erVev/993XffffJbrdr+/btevrpp3X11VerY8eOfseXdt/GjRv1/PPPq3nz5uwAApQTggoAAAAAAGAZrFEBAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwDIIKAAAAAABgGQQVAAAAAADAMggqAAAAAACAZUQFuwPBVlBQoB07diguLk4RERHB7g4AABXOMAxlZmaqQYMGiozkHEZ5YcwBAKjsyjrmqPRBxY4dO5ScnBzsbgAAEHRpaWlq1KhRsLsRthhzAABw1InGHJU+qIiLi5N09I2Kj48Pcm8Qcg4dkho0OPr7jh1S9erB7Q8AnIKMjAwlJyeb34koH4w5cFoYcwAIA2Udc1T6oKJw6mV8fDyDBpy8KlX+93t8PIMGACGNyxHKF2MOnBbGHADCyInGHFyICgAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJZBUAEAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAyCCoAAAAAAIBlEFQAAAAAAADLIKgAAAAAAACWQVABAAAAAAAsIyrYHQCsxOPxyOv1mrftdrscDkcQewQAAMIBYwwAKDuCCuD/eTweOZ1O+Xw+s81ms8ntdjOQAAAAp4wxBgCcHIIK4P95vV75fD6lpqbK6XTK7XbL5XLJ6/UyiAAAAKestDHGzz//LKfTWeIxzLgAUJkRVABFOJ1OpaSkBLsbAAAgzBSOMex2u2w2m1wuV6m1zLgAUJkRVAAAAAAVyOFwyO12+61ZcSxmdQKo7AgqgBNwu93m70zDBAAAgeBwOBhTAEApCCpQqR27AvexgYSkEqdlMg0TAAAAAMoXQQUqrdJW4Lbb7ZKKT8tkGiYAAAAAlD+CClQqRWdQHLsCt1T80g6mZQIAAABAxSKoQKVR2gyKXr16EUYAAAAAgEWcUlBx8OBBrVu3Ti1atDCnyRe1d+9ebdu2TU2bNlViYmLQa4Cie5hLLI4JAAAAAFYTeTLFf/75p6699lo5nU51795d33zzTbEawzA0btw4NWrUSFdccYXq1aunBx98MGg1QFGFe5inpKScUkjhdru1fPlyLV++XGlpaeXQQwAAAMYcACqvk5pRsWrVKnXv3l3Tp09X9erVS6x5+eWXlZqaquXLl6tt27ZasGCB+vTpo44dO+qiiy6q8BogUEraBcQeG6u9QewTAAAIP4w5AFR2JzWj4pJLLtF1110nm81Was2rr76qyy67TG3btpUk9ezZU/369dNrr70WlBogUAp3AVm2bJmWLVum1NRU+bKzg90tAABgQR6Px5wNUXQL9BNhzAGgsgvoYpr5+fn6448/dN111/m1d+3a1QwPKrKmJDk5OcrJyTFvZ2RknOSrRCgpusvH6WIXEACoWAsWLNAvv/yiPn36qEuXLsXuz8/P1+zZs/Xnn3+qcePGGjp0qKKjo4NWA0gn3gK9LBhzAKjMAhpUZGZmKjc3V7Vq1fJrT0pK0r59+yq8piRTpkzRww8/fGovECElEIMEAEBw/PbbbxozZozi4uL0+++/q0qVKsWCipycHA0aNEibN29W//799fLLL2vy5Mn64YcfFB8fX+E1QCEW8AaA0xPQoKJq1aqS5DdjQZKys7PN+yqypiSTJk3SbbfdZt7OyMhQcnJyGV4dQk1FDxJWrlypWsnJDEIAIABiY2P17rvvql27dqpXr16JNf/5z3+0cuVKrV27VvXq1VN6errOPPNMPfHEE3r88ccrvAYoqnABbwDAyTmpNSpOpHr16kpKStKOHTv82nfs2KHGjRtXeE1JYmJiFB8f7/eD8Ha6u3yUxm63yxYba97u0bOnnE6nPB5PwJ4DACqrtm3bql27dsetmTlzpoYOHWoGGQkJCRoxYoRmzpwZlBoAABAYAQ0qJKlfv36aNWuWeTs/P19fffWV+vXrF5QaoLw4HA4tX77cvP3aq6/K5/OZa2IAAMqX2+1Wq1at/NpatWqljRs3Kjc3t8JrisrJyVFGRobfDxAIK1eu5MQIgLB2UkFFRkaGFi1apEWLFkmSNm7cqEWLFmnTpk1mzf33368VK1bopptu0uzZszVq1CgdPHhQd9xxR1BqgPJ07GVDrVu3DmJPAKDyycrKUkJCgl9bYmKiDMPQoUOHKrymqClTpighIcH84VJTnCpmcQKobE4qqFi/fr0mTJigCRMmqGvXrpo9e7YmTJigDz74wKw588wz9csvvygzM1NPPPGEbDabfv31VzVq1CgoNQAAIDxVr1692CyF9PR0STK3Uq/ImqImTZqk9PR08yctLe3kXiDw/5jFCaCyOanFNM8++2xzNsXxdOjQQTNmzLBMDQAACD9nnHGG/vrrL7+2v/76S02aNDG3Da3ImqJiYmIUExNz6i8QOAazOAFUJgFfowIAAKAiXHzxxZo1a5b2798vSfL5fPrwww91ySWXBKUGAAAERkC3JwUAAAgEr9erV199VZJ06NAhzZ8/X3l5eWrRooUuu+wySdLNN9+sjz76SL169dLQoUM1b948Va1aVffee6/5OBVZAwAAAoMZFQAAwHLy8/N18OBBHTx4UOPGjVO7du108OBBv4UrbTabFixYoEmTJqlatWq68cYbtWzZMiUlJQWlBgAABAYzKgAAgOXUrVtXTzzxxAnrYmJi5HK5LFMDAABOHzMqAAAAAACAZRBUAAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIPFNAEAAIDT5PF45PV6JUlutzvIvQGA0EZQAQAAAJwGj8cjp9Mpn89nttlsNtnt9iD2CgBCF0EFAAAAcBq8Xq98Pp9SU1PldDolSXa7XQ6Ho1yf99iZGxXxfABQUQgqAAAAgABwOp1KSUkp9+dJSkqSzWaTy+Uy22w2m9xuN2EFgLBAUAEAAACEkOTkZLndbr81MVwul7xeL0EFgLBAUAEEGNMwAQAIf8FePNPhcDDGABC2CCqAAGEaJgAAlQOLZwJA+SKoQNg49syGVPFnN5iGCQBA5RCsxTMBoLIgqEBYKOnMhlTxZzeYhgkAQOVRUYtnAkBlQ1CBsFDSmQ2JsxsAAAAAEGoIKhBWrHhmg8U1AQAAAKDsCCqAcmK321lcEwAAAABOEkEFUE4cDgeLawIAAADASSKoAMoRi2sCAAAAwMmJDHYHAAAAAAAAChFUAAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsAy2JwUAAADCgNvt9rttt9vZJh1ASCKoAAAAAEKY3W6XzWaTy+Xya7fZbHK73YQVAEIOQQUAAAAQwhwOh9xut7xer9nmdrvlcrnk9XoJKgCEHIIKAAAAIMQ5HA4CCQBhg8U0AQAAAACAZRBUAAAAAAAAyyCoAAAAAAAAlsEaFQhZHo/HXDSq6HZcAAAAAIDQRFCBkOTxeOR0OuXz+cw2m80mu90exF4BAIBwdOzJEYkTJABQ3ggqEJK8Xq98Pp9SU1PldDolHd1DnNWuAQBAIJV0ckTiBAkAlCeCCoQ0p9OplJSUYHcDAACEqZJOjkicIAGA8kRQAQAAAJwAJ0cAoOKw6wcAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAyWKMCqGBFtzRjMS4AAAAA+B+CCqCC2O122Ww2uVwuv3abzSa3201YAQAAAAAiqAAqjMPhkNvtltfrNdvcbrdcLpe8Xi9BBQAACLhjZ3IyixNAqCCoACqQw+FggAAAAMpdSTM5mcUJIFQQVAAAAABhpuhMTmZxAgglBBUAAABAGGImJ4BQxfakAAAAAADAMggqAAAAAACAZRBUAAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsIyoYHcAKCuPxyOv1ytJcrvdQe4NAAAAAKA8EFQgJHg8HjmdTvl8PrPNZrPJbrcHsVcAAAAAgEAjqEBI8Hq98vl8Sk1NldPplCTZ7XY5HI4g9wwAAAAAEEgEFQgpTqdTKSkpwe4GAAAAAKCcEFQAFnDsmhvMFAEAAABQmRFUAEFkt9tls9nkcrnMNpvNJrfbTVgBAECQsIA3AAQXQQUQRA6HQ263228w5HK55PV6CSoAAAgCFvAGgOAjqACCzOFwEEoAAGARLOANAMFHUAEAAAAUwQLeABA8kcHuAAAAAAAAQCGCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGWwmCYAAABQSbjdbvN3djMBYFUEFQAAAECYs9vtstlscrlcZpvNZpPb7SasAGA55RJUfP755/r444+1Z88eNWjQQCNHjlS/fv38atauXaupU6dq69atatmype688041bty4XGoAAACAyszhcMjtdsvr9Uo6OrPC5XLJ6/USVACwnICvUTF9+nRdccUVat++vSZMmCCHw6H+/ftr5syZZs26devUrVs35efn64YbbtDOnTvVpUsX7dy5M+A1AAAAAI6GFSkpKUpJSZHT6Qx2dwCgVAEPKmbOnKkRI0bo9ttv18CBA/XQQw/pvPPO00cffWTWPProo2rbtq1ee+01XXzxxfrggw9ks9n0zDPPBLwGAAAAAACEjoAHFR06dNDq1auVk5MjScrIyNCGDRvUsWNHs2bOnDm64IILzNtRUVEaOnSo5syZE/AaIBS53W4tX75cy5cvl8fjCXZ3AMCy5s2bp2uuuUYDBw7U6NGj9fnnnxer2bp1q8aPH68BAwbo+uuv19q1a8utBgAAnL6ABxVPP/20unXrJofDoe7du6tp06YaPXq07rrrLkmSz+eT1+tVw4YN/Y5r2LChtm7dGtCakuTk5CgjI8PvB7CKYxe66tSpkzp16iSn00lYAQAlmDlzpgYMGKAmTZpowoQJOuusszRixAhNnz7drCm8LHTXrl264YYblJubq27dumndunUBrwEAAIER8MU0P/vsM73zzju655571K5dOy1atEj//ve/1bt3b/Xt21dHjhyRJMXGxvodZ7PZzPsCVVOSKVOm6OGHHz6NVwiUHxa6AoCymzlzps477zw98MADkqSBAwdq1apV+uijj3TzzTdLkv71r38pLi5OH3zwgapUqaKLL75Ya9as0eTJk5WamhrQGgAAEBgBn1Fx66236uabb9bdd9+twYMH65FHHtHw4cN1++23S5Li4uIUFRWlffv2+R23b98+1axZM6A1JZk0aZLS09PNn7S0tNN+zUAgsdAVAJRNhw4dtGHDBmVmZkqSDh8+rNWrVxe73HTIkCGqUqWK2XbhhRcWu5Q0EDUAACAwAjqjIj8/X/v371dycrJfu8Ph0Ndffy1JqlKlitq1a6fly5f71SxdutQcWASqpiQxMTGKiYk5tRcIAAAs4+6779bhw4fVtGlTtWzZUhs3btSIESP01FNPmTVbt24t8TLRvXv3Kjs7W7GxsQGrKSonJ8dcs0sSl5sCAFBGAZ1RUaVKFXXt2lUzZsyQz+eTJB08eFAffPCBevToYdaNGTNGM2fO1IYNGyRJS5Ys0dy5czVmzJiA1wAAgPD0/fffa9q0abrpppt0//336/bbb1dqaqrfgppHjhwp8TLRwvsCWVPUlClTlJCQYP4UPZEDAABKFvA1Kl577TVdccUVcjgcat68udavX6+zzjpL//73v82am266SStXrlT79u11xhlnaP369ZowYYIuu+yygNcAAIDwdMcdd+jyyy/XI488IkkaPHiwsrKydMstt2j48OGSpJo1a5Z4mWjVqlVVo0aNgNYUNWnSJN12223m7YyMDMIKAADKIOBBxRlnnKFly5Zp69at2r17txo2bKhGjRr51VSpUkWvvfaaHn74YaWlpalZs2aqW7duudQAAIDw5PV6S7zcdN++fcrPz1eVKlXUsWPHEi8TbdeunbneRKBqiuJyUwAATk3AF9OUpIiICDVp0kRdu3YtFlIcq1GjRurevftxw4VA1QAAgPByzjnn6IMPPtDBgwclHd26fMaMGerWrZsZHowZM0Zz5szR0qVLJUnr16/Xxx9/XOxS0kDUAACAwAj4jAoAAICKMG3aNI0YMUJNmjRRq1attHHjRjkcDr3//vtmzfDhw7V06VL16tVLrVq10oYNG3TVVVfppptuCngNEIrcbrf5u91uZzt0AJZAUAEAAEJS/fr19dNPP2nbtm3avn276tatq8aNGysiIsKv7qmnntLtt9+uTZs2KTk5ucTZnoGqAUKF3W6XzWaTy+Uy22w2m9xuN2EFgKAjqAAAACGtUaNGJwwN6tate8JLRANVA4QCh8Mht9str9cr6ejMCpfLJa/XS1ABIOgIKgAAAIBKyOFwEEoAsKRyWUwTAAAAAADgVDCjAgAAAJWax+PxuwQCABBcBBUAAACotDwej5xOp3w+n9lms9lkt9uD2CsAqNwIKmBZnN0AAADlzev1yufzKTU1VU6nUxLbdAJAsBFUwJI4uwEAACqS0+lUSkpKsLsBABBBBSyKsxsAAAAAUDkRVMDSOLsBAAAAAJUL25MCAAAAAADLIKgAAAAAAACWQVABAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwDHb9AEKA2+02f2ebVgAAAADhjKACsDC73S6bzSaXy2W22Ww2ud1uwgoAAAAAYYmgArAwh8Mht9str9cr6ejMCpfLJa/XS1ABAAAAICwRVAAW53A4CCUAAAAAVBospgkAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAy2PUDAAAAgKSjW6EXstvt7DwGICgIKgAAAIBKzm63y2azyeVymW02m01ut5uwAkCFI6gAAAAAKjmHwyG32y2v1yvp6MwKl8slr9dLUAGgwhFUAAAAoFLxeDx+/yDHUQ6Hg1ACgCUQVAAAAKDS8Hg8cjqd8vl8ZpvNZpPdbg9irwAAxyKoAAAAQKXh9Xrl8/mUmpoqp9MpiUUjAcBqCCoAAABQ6TidTqWkpAS7GwCAEkQGuwMAAAAAAACFCCoAAAAAAIBlcOkHEIKOXaGc62oBAAAAhBOCClgGW4WdmN1ul81mk8vlMttsNpvcbjdhBQAAAICwQFABS2CrsLJxOBxyu91+gY7L5ZLX6yWoAAAAABAWCCpgCWwVVnYOh4P3BQAAAEDYIqiApbBVGAAAAABUbuz6AQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJZBUAEAAAAAACwjKtgdAAAAAGBNbrfb/N1ut8vhcASxNwAqC4IKAAAAAH7sdrtsNptcLpfZZrPZ5Ha7CSsAlDuCCgAAAAB+HA6H3G63vF6vpKMzK1wul7xeL0EFgHJHUAEAAACgGIfDQSgBIChYTBMAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAyCCoAAAAAAIBlEFQAAAAAAADLIKgAAAAAAACWERXsDgA4fW632/zdbrez5zkAAACAkEVQAYQwu90um80ml8tlttlsNrndbsIKAAAAACGJoAIIYQ6HQ263W16vV9LRmRUul0ter5egAgAAAEBIIqgAQpzD4SCUAADgODwej1+oDwCwNoIKAAAAhC2PxyOn0ymfz2e22Ww22e32IPYKAHA8BBUAAAAIW16vVz6fT6mpqXI6nZJYeBoArI6gAgAAAGHP6XQqJSUl2N0AAJRBZLA7AAAAAAAAUIigAgAAAAAAWAaXfiBoWIEbAAAAAFAUQQWCghW4AQAAQs+xJ5dYlBRAeSm3oCI/P1+rVq1SRESE2rVrp8jI4leZbN26VR6PRy1atFD9+vVLfJxA1cBaWIEbAAAgdNjtdtlsNrlcLrPNZrPJ7XYzfgMQcOWyRsWcOXPUpEkTXXbZZbrhhht09tlna+PGjeb9+fn5uvrqq9WmTRtNmDBBzZo10+233+73GIGqgbUVrsCdkpLClxwAAIBFORwOud1uLVu2TMuWLVNqaqp8Pp95GS8ABFLAZ1SsXr1aF1xwgR555BHdddddkqS1a9dq7969at68uSTp+eef1xdffKHff/9dLVq00G+//aYePXqoa9euuvzyywNaAwAAAOD0ORwOTiwBqBABn1HxxBNPqFWrVmZIIUlt2rRRt27dzNtvvvmmLr/8crVo0UKS1LlzZ/Xv319vvvlmwGsAAED48nq9uu6662S32xUfH6/Ro0dr3759fjWfffaZWrduraioKDVr1kxvvPFGsccJVA0AADh9AQ8qvv/+ew0ePFhZWVn65Zdf9Oeff6qgoMC8Py8vT6tXr1ZKSorfcSkpKVq5cmVAa0qSk5OjjIwMvx8AABB6srOz1adPH23atEmLFy/Wvn37NGTIEH355ZdmzZIlSzR8+HCNHz9eBw4c0COPPKLrr79eX3/9dcBrAABAYAQ0qDAMQ7t27dLmzZvVpk0b3XHHHerVq5c6deqkv/76S5KUlZWlvLw81apVy+/YpKQkHThwIKA1JZkyZYoSEhLMn+Tk5NN+3QAAoOK98MIL2rJli2bOnKnmzZuratWqGjFihK6++mqz5t///re6dOmi8ePHKy4uTi6XS4MHD9azzz4b8BoAABAYAQ0qIiIiVLVqVc2dO1e//PKLfv31V23evFnx8fG6/vrrJUnR0dGSjp4FOZbP5zPvC1RNSSZNmqT09HTzJy0t7VRfLgAACKLPP/9c/fv3V1JSUqk1CxcuVO/evf3a+vbtq0WLFskwjIDWAACAwAj4pR9NmjTReeedZ85UiI2N1ZVXXqmFCxfKMAzZbDbVrl1b27Zt8ztu27ZtatKkiSQFrKYkMTExio+P9/sBAAChZ/PmzapXr54uvfRSxcTEqG7durrmmmv81qjYvXu3ateu7XdcnTp1dOjQIWVlZQW0piguNwUA4NQEPKgYNGiQtm7d6te2ZcsW1alTRxEREZKk/v376/PPPzfvz83N1Zdffqn+/fubbYGqAQAA4ckwDP3nP//RoEGDlJ6erh9//FG//fabRo0addzjCtfOKhyXlFcNl5sCAHBqAh5U3H333dq6dauuv/56ffHFF5oyZYqmTZum++67z6y5//77tW7dOl199dWaOXOmhg8frpycHN1+++0BrwEAAOGpfv366ty5s6677jpVq1ZNrVu31qRJk/Ttt9/q0KFDZs2ePXv8jtu7d69q1Kih6tWrB7SmKC43BQDg1AQ8qKhfv76WLl2quLg4vfLKK9q0aZO++uorc40KSWrVqpUWL14sm82mN998U02aNNHixYtVr169gNcAlZHb7dby5cu1fPlyeTyeYHcHAMpFr169/HYWk6T8/HxFREQoMvLoEKdHjx6aP3++X828efPUvXt3cyZEoGqK4nJTAABOTVR5PGjDhg31zDPPHLemdevWeumllyqkBqgs7Ha7bDabXC6X2Waz2eR2u+VwOILYMwAIvFtuuUX//e9/9fzzz2vUqFHavHmzHn/8cV1yySWKjY2VJE2YMEHdunXTv/71L1199dX64osvNHv2bL9tRQNVAwAAAiPgMyoABI/D4ZDb7dayZcu0bNkypaamyufzyev1BrtrABBwTZs21Xfffaf33ntPDRs21MUXX6wBAwbotddeM2tSUlL02WefKTU1VU2aNNG//vUvvf32237rWQWqBgAABEa5zKgAEDwOh4PZEwAqje7du+uXX345bs3gwYM1ePDgCqkBAACnjxkVAAAAAADAMggqAAAAAACAZXDpBwAAAMKKx+Mx12dyu91B7g0A4GQRVAAAACBseDweOZ1O+Xw+s81ms8lutwexVwCAk0FQAQAAgLDh9Xrl8/mUmpoqp9Mp6ej23Sw0DQChg6ACAAAAYcfpdColJSXY3QAAnAIW0wQAAAAAAJZBUAEAAAAAACyDSz8AAAAAnJJjd1VhLRAAgUJQAQAAAOCk2O122Ww2uVwus81ms8ntdhNWADhtBBUAAAAATorD4ZDb7ZbX65V0dGaFy+WS1+slqABw2ggqAAAAAJw0h8NBKAGgXLCYJgAAAAAAsAxmVKDCeDwev+mBAAAAAAAURVCBCuHxeOR0OuXz+cw2m80mu90exF4BAAAAAKyGoAIVwuv1yufzKTU1VU6nUxJbWAEAAAAAiiOoQIVyOp1KSUkJdjcAAAAAABbFYpoAAAAAAMAyCCoAAAAAAIBlEFQAAAAAAADLYI0KoBI4djtYFjEFAAAAYGUEFUAYs9vtstlscrlcZpvNZpPb7SasAAAAAGBJBBVAGHM4HHK73fJ6vZKOzqxwuVzyer0EFQAAAAAsiaACCHMOh4NQAgAAVAguNwUQCAQVAAAAAE4Ll5sCCCSCCgAAAACnhctNAQQSQQUAAACA08blpgACJTLYHQAAAAAAAChEUAEAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAyCCoAAAAAAIBlEFQAAAAAAADLIKgAAAAAAACWQVABAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwjKhgdwBAxXO73ebvdrtdDocjiL0BAAAAgP8hqAAqEbvdLpvNJpfLZbbZbDa53W7CCgAAAACWQFABVCIOh0Nut1ter1fS0ZkVLpdLXq+XoAIAAACAJRBUAJWMw+EglAAAAABgWSymCQAAAAAALIMZFQAAAADKBQt4AzgVBBUoNx6Px28tBAAAAFQOLOAN4HQQVKBceDweOZ1O+Xw+s81ms8lutwexVwAAIBxxcsR6WMAbwOkgqEC58Hq98vl8Sk1NldPplMR0PwAAEHicHLEuFvAGcKoIKlCunE6nUlJSgt0NAAAQpjg5AgDhh6ACAAAAIY+TIwAQPtieFAAAAAAAWAZBBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDJYTBOA357zrJQOAAAAIJgIKoBKzG63y2azyeVymW02m01ut5uwAgAAAEBQEFQAlZjD4ZDb7ZbX65V0dGaFy+WS1+slqAAAAAAQFAQVQCXncDgIJQAAAABYBotpAgAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIM1KgAAAABUiGO3RJfYFh1AyQgqAAAAAJSrkrZEl9gWHUDJCCoAAAAAlKuiW6JLbIsOoHQEFQAAAADKHVuiAygrFtMEAAAAAACWQVABAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwjHINKr777jt169ZN9957b7H7Vq5cqdGjR6t379667rrr9Ndff5VbDSqGx+PR8uXLtXz58mJ7ZAMAAAAAUBblFlTs3r1b1113nfbv368NGzb43bd69Wr17NlTNWrU0D333KNDhw6pW7du2rZtW8BrUDE8Ho+cTqc6deqkTp06yeVyyWazyW63B7trAIBKIDc3V99++60WL15c4v1paWmaP3++Nm3aVOpjBKoGAACcnnIJKgzD0KhRozRhwgSdccYZxe5/9NFH1aFDB7344osaNGiQUlNTlZiYqGeeeSbgNagYXq9XPp9PqampWrZsmZYtWya3280WVACACnHvvfdqyJAhuv3224vdd+utt6pVq1a65557dNZZZ2nkyJHKy8srlxoAAHD6yiWoePLJJ2UYhiZMmFDi/fPmzdOwYcPM21WqVNGQIUM0d+7cgNegYjmdTqWkpCglJYWQAgBQIb799lt9+eWXGj58eLH73nvvPb3yyitauHChFi9erN9//11ff/21XnjhhYDXAACAwAh4ULFo0SJNnTpVb731liIiIordf+jQIe3bt08NGjTwa2/QoIG2bt0a0JqS5OTkKCMjw+8HAACEpl27dunaa6/VjBkzVL169WL3v/XWWxowYIA6dOggSWrevLmGDx+uN998M+A1AAAgMAIaVKSnp+vKK6/Uiy++WCxAKJSbmytJiomJ8WuPjY017wtUTUmmTJmihIQE8yc5ObmsLw8AAFhIQUGBRo0apXHjxqlz584l1qxcuVIdO3b0a+vYsaNWr16t/Pz8gNYUxckRAABOTUCDih9++EE7duzQU089pW7duqlbt2765ZdfNH/+fHXr1k07d+5UXFycqlatqv379/sdu2/fPiUlJUlSwGpKMmnSJKWnp5s/aWlpgXjpAACggj3xxBPKycnR3XffXWrNwYMHVatWLb+2pKQk5eXlKSsrK6A1RXFyBACAUxMVyAc799xz9eOPP/q13XHHHYqOjtbjjz+uWrVqqUqVKjrrrLO0dOlS3XDDDWbd4sWLzTMVgaopSUxMTLFZGAAAILRs3LhRDz30kKZPn645c+ZIkrZt26aDBw/qm2++Uffu3ZWQkKDo6GhlZ2f7Hevz+SRJ0dHR5v8GoqaoSZMm6bbbbjNvZ2RkEFYEiMfjkdfrlSS2RAeAMBTQoKJmzZrq1q2bX1tiYqKqVavm137ttdfqnnvu0e233642bdpowYIFmjdvnmbOnBnwGgAAEH5yc3PVt29fffrpp2bbmjVrlJWVpalTp6pJkyZKSEhQkyZNis2e3LZtm+rUqaPY2FhJClhNUZwcKR+FW6IXBkWS2BIdAMJMQIOKsrrhhhu0evVqdezYUU2aNNHWrVs1adIkXXTRRQGvAQAA4ad169b65ptv/Nquu+46rVu3zq994MCB+uijj/Tss8+qatWqMgxDn376qQYOHBjwGlSMY7dEdzqdkiS73c5uYyHs2Fkx/H8JQJIiDMMwyvMJ1q9fr8jISLVs2bLYfXv37tW2bdvUpEkT1axZs8TjA1VTmoyMDCUkJCg9PV3x8fEndSz+Z/ny5erUqZOWLVumlJSUYHen4hw6JNWocfT3rCyphBXnQ0ml/f8RqOTC5buwMKhYsGCB2bZ792517NhRnTt31pVXXqlZs2bp66+/1tKlS82xSaBqTiRc3udgq7TfVWE25pBKnx3jdrsJK4AwVdbvwnKfUdGqVatS76tdu7Zq16593OMDVQOg7Ipe78vZDQCh4Mwzz1RiYqJfW926dbVkyRJNnTpV77//vhwOR7FwIVA1AE6Ow+GQ2+32W2/E5XLJ6/Uy7gAquaBc+gHAmux2u2w2m1wul187ZzcAhIIJEyaU2N6oUSP961//Ou6xgaoBcHIcDgfjCwDFEFQAMBU9syFxdgMAAABAxSKoAOCHMxsAAAAAgiky2B0AAAAAAAAoRFABAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwDIIKAAAAAABgGQQVAAAAAADAMggqAAAAAACAZRBUAAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALCMq2B0AEBrcbrf5u91ul8PhCGJvAAAAAIQrggoAx2W322Wz2eRyucw2m80mt9tNWAEAAAKOkyMACCoAHJfD4ZDb7ZbX65V0dPDgcrnk9XoZOAAAgIDh5AiAQgQVOGUej8fvH68IXw6HgwECAAAoV5wcAVCIoAKnxOPxyOl0yufzmW02m012uz2IvQIAAEAo4+QIAImgAqfI6/XK5/MpNTVVTqdTEtcQAgAAAABOH0EFTovT6VRKSkqwuwEAAAAACBORwe4AAAAAAABAIYIKAAAAAABgGQQVAAAAAADAMggqAAAAAACAZbCYJgAAACzN4/HI6/VKktxud5B7AwAobwQVAAAAsCyPxyOn0ymfz2e22Ww22e32IPYKAFCeCCoAAABgWV6vVz6fT6mpqXI6nZIku90uh8MR5J4BAMoLQQUAAAAsz+l0KiUlJdjdAABUABbTBAAAAAAAlkFQAQAAAAAALINLPwCckmNXXedaYQAAAACBQlAB4KTY7XbZbDa5XC6zzWazye12E1YAAICA4+QIUPkQVAA4KQ6HQ263228/e5fLJa/Xy8ABAAAEDCdHgMqLoALASXM4HAwQAABAueLkCFB5EVQAAAAAsCROjgCVE7t+AAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGVEBbsDAMKD2+02f7fb7XI4HEHsDQAAAIBQRVAB4LTY7XbZbDa5XC6zzWazye12E1YAAAAAOGkEFQBOi8PhkNvtltfrlXR0ZoXL5ZLX6yWoAAAAAccsTiD8EVQAOG0Oh4NBAgAAKFfM4gQqD4IKAAAAAJbHLE6g8iCoAAAAABASmMUJVA5sTwoAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAM1qhAmXk8Hr/FiwAAAMoDYw4AqNwIKlAmHo9HTqdTPp/PbLPZbLLb7UHsFQAACDeMOQAABBUoE6/XK5/Pp9TUVDmdTklH97Jm1WUAABBIjDkAAAQVOClOp1MpKSnB7gYAAAhzjDkAoPIiqABQLo69ppgzYQAAAADKiqACQEDZ7XbZbDa5XC6zzWazye12E1YAAICA4+QIEH4IKgAElMPhkNvt9lut3eVyyev1MnAAAAABw8kRIHwRVAAIOIfDwQABAACUK06OAOGLoAIAAABASOLkCBCeIoPdAQAAAAAAgEIEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJbBrh8AAAAAwobb7TZ/t9vt7AoChCCCCgAAAAAhz263y2azyeVymW02m01ut5uwAggxBBUAAAAAQp7D4ZDb7ZbX65V0dGaFy+WS1+slqABCDEEFAAAAgLDgcDgIJYAwUC5BRXZ2tn7//Xfl5uaqbdu2qlWrVol1a9eu1datW9WyZUu1aNGiXGtwcjwej5lGS/7X+gGngutFAQAAAJRFwIOKxx57TC+++KKSk5NVpUoV/f7773r00Uc1ceJEsyY3N1dXXnml5s6dq/bt22vZsmUaOXKk/vOf/ygiIiKgNTh5Ho9HTqdTPp/Pr91ms8lutwepVwhVXC8KoDytXLlSc+fOVVZWljp06KALLrhAkZH+m5plZ2frnXfe0Z9//qnGjRtr5MiRSkhIKJcaAABw+gIeVCQkJGjt2rXmF/fHH3+syy67TL169VLnzp0lSVOnTtUPP/ygP/74Qw6HQ6tWrVKXLl107rnnauTIkQGtwcnzer3y+XxKTU2V0+k02zkLjlPB9aIAysu1116rlStXqk+fPqpRo4YmTpyoZ555RnPnzlVMTIwkKSsrS7169ZJhGLrwwgs1Y8YMPfPMM1q0aJFq164d0BqcGmZxAgCKMcpZQUGBER0dbbz88stmW7t27YybbrrJr+7CCy80+vfvH/CaE0lPTzckGenp6WU+JtwtW7bMkGQsW7Ys2F2xvqwsw5CO/mRlBbs3IYG/L8B6QvW7cPny5X63d+zYYURHRxuvv/662fbYY48ZderUMQ4ePGgYhmH4fD6jRYsWxi233BLwmhMJ1fe5PG3dutWw2WyGJL8fm81mbN26NdjdsxbGHCeNMQdgPWX9LowsJb8ImEWLFunIkSPmmfm8vDytXbtW7du396vr0KGD/vjjj4DWlCQnJ0cZGRl+PwAAIPR07NjR73b9+vVlt9u1fft2s+3zzz/XsGHDzJmesbGxGj58uD777LOA1+DkHTuLc9myZeYPlwcCQOVWrkHFwYMHNWbMGA0aNEi9evWSdHTqZH5+vmrWrOlXm5SUpIMHDwa0piRTpkxRQkKC+ZOcnHx6LxIAAFjCd999px07dujcc8812zZs2KDmzZv71TVv3lwej0c5OTkBrSmKkyNl53Q6lZKSYv4QUgBA5VZuQUVWVpYGDx6s6tWr69133zXbo6OjJanYQo1ZWVnm9aSBqinJpEmTlJ6ebv6kpaWdyssDAAAWsmXLFo0ePVp///vf/YIKn8+nuLg4v9r4+HhJRxfHDGRNUZwcAQDg1JRLUJGVlaVBgwbp8OHD+u6775SYmGjeZ7PZVLdu3WIBQVpampo2bRrQmpLExMQoPj7e7wcAAISubdu2qV+/furSpYteeeUVv/vi4uKKzbQ8cOCAIiIiVL169YDWFMXJEQAATk3Ag4pDhw5p8ODBOnTokObOnatatWoVqxk0aJA+/fRTFRQUSDo6NXLWrFkaNGhQwGsAAED42r59u/r06aM2bdpo5syZqlq1qt/9bdq00bp16/za1q1bp5YtW5q1gaopipMjAACcmoAHFRdccIFWrFih8ePH66efftJnn32mzz77zO/L/f7779eWLVt0+eWX64033tCwYcMUERGh2267LeA1AAAgPO3YsUN/+9vf5HQ69fHHH5uXhR7r8ssv15dffmkusLl//3598MEHuvzyywNeAwAAAiPgQUVSUpL69eunL774Qm+++ab5s3LlSrOmWbNmWrZsmZo2bapvvvlGXbp00dKlS2W32wNeAwAAwtMll1wij8ejRo0a6a677tKECRM0YcIEzZw506wZO3asunXrpq5du+rqq6/W2WefrUaNGunuu+8OeA0AAAiMqEA/4IcfflimuiZNmujpp5+ukBoA1uN2u83f7XY7K7wDOGnXX399iTtpJCUlmb9XrVpVs2fP1rx58/Tnn3/q8ssv14ABAxQVFRXwGgDWxJgDCD18uwKoUHa7XTabTS6Xy2yz2Wxyu90MHACclGuuuaZMdRERETrvvPN03nnnlXsNAOtgzAGELoIKABXK4XDI7XbL6/VKOnqWw+Vyyev1MmgAAAABw5gDCF0EFQAqnMPhYIAAAADKHWMOIDQFfDFNAAAAAACAU0VQAQAAAAAALIOgAgAAAAAAWAZrVECS5PF4/BYaAgAAAAAgGAgqII/HI6fTKZ/PZ7bZbDbZ7fYg9gqVDXucAwCAisCYA7A+ggrI6/XK5/MpNTVVTqdTEh/aqDjscQ4AACoCYw4gdBBUwOR0OpWSkhLsbqCSYY9zAABQERhzAKGDoAJA0LHHOQAAqAiMOYDQQFABAACACsMC3gCAEyGoAAAAQIVgAW8AQFkQVAAAAKBCsIA3AKAsCCoAAABQoVjAGwBwPJHB7gAAAAAAAEAhggoAAAAAAGAZXPoBwJKOXQme65cBAEB5YcwBWA9BBQBLsdvtstlscrlcZpvNZpPb7WbgAAAAAoYxB2BdBBUALMXhcMjtdsvr9Uo6epbD5XLJ6/UyaAAAAAHDmAOwLoIKAJbjcDgYIAAAgHLHmAOwJhbTBAAAAAAAlsGMCgAAAAD4f8curimxwCYQDAQVAAAAACq9khbXlFhgEwgGgopKyuPx+C0cBAAAAFRmRRfXlFhgEwgWgopKyOPxyOl0yufzmW02m012uz2IvQKOjz3OAQBAeWNxTcAaCCoqIa/XK5/Pp9TUVDmdTkn8ww/WxR7nAAAAQOVCUFGJOZ1OpaSkBLsbwHGxxzkAAABQuRBUALA8pmECAAAAlQdBBQAAAAAcB2tlARWLoAIAAADlhp3GEMpYKwsIDoIKAAAAlAt2GkOoY60sIDgIKgCEpKJn5ZiGCQDWw05jCAeslQVUPIIKACGlpCmYEtMwAcDK2GkMAHAyCCoAhJSiUzAlpmECAAAA4YSgAkDIYQomAAAAEL4IKgCEDbYOAwAAAEIfQQWAkMfWYQAAAED4IKgAEPLYOgwAAFQkZnEC5YugopLweDx+/4gDwg3rVgAAgPLGLE6gYhBUVAIej0dOp1M+n89ss9lsstvtQewVAAAAEFpKm8X5888/y+l0SmKGBRAIBBWVgNfrlc/nU2pqKh+gAAAAwGk4dhYnMyyA8kFQUYk4nU6lpKQEuxsAAABAWGCdLKB8EFQAAAAgYFgXC5UN62QBgUdQAQAAgIBgXSwAQCAQVAAIW2wdBgAVi3WxgKMYgwCnh6ACQNgpbWGrTz75RLVr1zZrGDQAQPlgXSxUViyuCQQGQQWAsFN0Yau9e/fqkksu0cCBA80aBg0AACDQWFwTCAyCCgBhqejCVgwaAABARWBxTeD0EVQAqBQYNAAAAAChgaACAAAAAMoRi2sCJ4egIkyxhzkAAAAQXCyuCZwagoowxB7mQNlwdgMAAJQnFtcETg1BRRhiD3Pg+Di7AQAAKgrrZAEnj6AijLGHOVAyzm4AAIBgYlYncHwEFQAqpZLObjBoAAAA5YlZnUDZEFQAqPQYNADAqWMBb6DsmNUJlA1BBYBKr7RBw88//8w6LwBwHCzgDZw8ZnUCJ0ZQAQDyHzQwwwIAyoYFvIHTw5gDKBlBBQAUwbRMADg5LOANnBrGHEDJCCrCBNeHAoHFtEwAAFARGHMAxRFUhAGuDwXKF9MyAQBARShtzPHJJ5+odu3aZg3jD4Q7goowwPWhQPliWiYAAKgIRccce/fu1SWXXKKBAweaNZwsQWVAUBFGuD4UKD9MywSAo7jcFChfRcccJ9qZTGIcgvBDUAEAJ4lpmQAqKy43BSreiXYmk5hlgfBDUBGiOJsBBA/TMgFUVlxuCgRX0TGIxCWpCE8EFSGIsxlA8JVlWiYDBgDhistNgeAp6XJUiUtSEV4IKkIQZzMA6znRGhYS/50CCE3M4gSsjUtSEY4IKkIYZzMAazre9aMMGgCEEmZxAtZX1ktSGYMglBBUhAjOZgCho6TrRxk0AAhFzOIEQsPxLkllLS2EIoKKEMDZDCD0lHYpCGc7AFhdSSdHmMUJhJZT2eL0WIxBEGwhH1T8+uuvevrpp7V161a1bNlS9913n84888xgd+u0FR0kcDYDCH2ncrbj2OCiKD4HgIq1evVqTZ48WX/++acaN26sO++8U927dw92twKKkyNAeCrLFqfHOtEYpCjGJAi0kA4qVqxYoT59+uiWW27RxIkT9dZbb6lnz55auXKlmjRpEuzunbLSBgm9evXiAwAIIycbXBRVdNrmsQGnxKABCKTNmzerZ8+eGj58uKZOnaovvvhCffv21cKFC9WxY8dgd++0cHIEqFxKukT1WGUZgxTFrFAEWoRhGEawO3GqLr30Uu3fv1/z58+XJBmGodatW+v888/X888/X6bHyMjIUEJCgtLT0xUfH1+e3T2uooMEl8vFICEUHDok1ahx9PesLKl69eD2B2GlaPBwrKKfE4WDiqIBJ4MGlIVVvgutbNy4cZo/f77WrFmjiIgISdK5556rOnXq6KOPPirTY1jxfS7t5AjXrlsQYw5UoOONQYoqyxjkVDBuCU9l/S4M6RkV8+fP16RJk8zbERERGjx4sObMmRPEXp1Y0f/wS/uPmxkUQOVW2j7pUulbkX3zzTeqXbv2KV1KUh4YZCBczJ8/X4MGDTJDCkkaMmSInnzyySD2qmxOFHoygwJAUccbg5TkZGeFlsXpjlv4LAttIRtUHDp0SAcOHFD9+vX92uvXr6+0tLRSj8vJyVFOTo55Oz09XdLRZCdQdu3apV27dpV4n9frlcvlUnZ2tl97bGysPv74Y/Ma0KSkJCUmJga0XygHhw797/eMDCk/P3h9QaWSmJioJUuWaN++fWZbUlKSkpOTJUktWrTwu7/ws+d0Bw0nKzY2VqmpqVzfHmD16tVTvXr1AvZ4hd81ITzJstylpaWVOOY4cOCAfD6fbDZbsWOCPeaQSh93HCs2NlYdOnQwPz8C3UcECGMOWFhiYqISExMlFR+DnIpAjFsYgwRGsMYcIRtU5OXlSZKio6P92mNiYpSbm1vqcVOmTNHDDz9crP3YL+dgyM7O1qWXXhrUPuA0NWgQ7B4AlsNnW2jJzMxUQkJCsLthSXl5eSWOOSSVOu6w6pijqOzs7LBYiLxSYcwBnBBjEGs70ZgjZIOKuLg4RUdHF0vq9u3bd9zUbNKkSbrtttvM2wUFBdq/f7+SkpL8pnNWFhkZGUpOTlZaWpplrpcNNbyHp4/3MDB4H09fZX0PDcNQZmamGvCPn1IlJSWVOOaIiYlRXFxciccw5ihZZf3vLBB4704P79/p4f07Pbx/R5V1zBGyQUVkZKQ6dOigxYsX68YbbzTbFy5cqE6dOpV6XExMjHkGpFDhNKXKLD4+vlL/BxMIvIenj/cwMHgfT19lfA+ZSXF8nTp10uLFi/3aFi5cqA4dOigyMrLEYxhzHF9l/O8sUHjvTg/v3+nh/Ts9vH9lG3OU/M0aIsaOHauPPvpIy5YtkyTNmTNHP/zwg66//vog9wwAAISTsWPHat68eZo7d64k6bffftOnn36qsWPHBrlnAACEn5CdUSFJ11xzjdavX68ePXooKSlJBw4c0JQpUzR48OBgdw0AAISRwYMH67HHHtMFF1ygmjVrat++fbr11ls1ZsyYYHcNAICwE9JBhSQ9+eSTuu+++7Rr1y41bNiwxFW3UbqYmBg9+OCDxaamoux4D08f72Fg8D6ePt5DHM9dd92l8ePHa/v27apXr16pa1Pg+Pjv7NTx3p0e3r/Tw/t3enj/Tk6EwV5kAAAAAADAIkJ6jQoAAAAAABBeCCoAAAAAAIBlEFQAAAAAAADLIKgIYzk5OXr00UfVvXt39ezZU08//bTy8vKOe0ynTp105pln+v28+eabFdNhi5o/f74uv/xynXnmmfrll1/KdMy8efM0bNgwde7cWaNHj9bGjRvLuZfWtmPHDj3yyCNKSUnRuHHjTliflpZW7O/wZN7/cLR+/XrdfPPN6tGjhwYOHKjnnntOOTk5JzxuxowZOu+883T22Wfr1ltv1f79+yugt9aUnZ2t6dOna8iQITrnnHM0duxY/fnnn8c95rvvvivxb3Hfvn0V1GsgtOTl5enpp59Wz5491b17dz366KMn/Kz629/+Vuy/senTp1dQj61l165dmjx5slJSUnTDDTeU6ZgDBw5o4sSJ6tKli/r166e33nqrnHtpXW63W7feeqvatWtXpr8hPuP/Z/fu3br//vvVp08f9enTR//85z/LNGaYP3++LrjgAnXu3FmjRo064fdquFqzZo3GjRunc845R4MGDdL06dN15MiR4x4zevToYn97d955ZwX12PpCftcPlG7MmDFatGiRpk2bpiNHjmjcuHHasmWLXnjhhVKPWbNmjZ599lmde+65Zlv9+vUroruW9MADD+inn37SFVdcoZkzZyozM/OEx8yfP18DBw7UAw88oN69e+uFF15Qjx49tGrVKtWuXbsCem0t27dvV/fu3TV69GjVrFlTW7duPeExOTk5WrNmjb799ls1aNDAbG/SpEk59tS6Nm3apEsuuUQ33XSTrrzySu3cuVN33HGHFixYoI8//rjU41588UXdcccdmjZtmpo2baoHHnhA559/vpYsWaIqVapU4CuwhpEjRyo5OVk33XSTatasqRdffFFnn322fvvtN7Vo0aLEY9LT0+XxeLRw4UK/9oSEhIroMhBybrnlFn322Wd64YUXVLVqVd18881at26d3nnnnVKPWbdunSZMmKChQ4eabXXq1KmI7lrK7t271aVLF40aNUp2u11btmw54TH5+fkaMGCAoqKiNGXKFG3dulU33nijMjIydPPNN5d/py3k22+/1cSJE3X99dcrMzNTu3fvPuExfMb/T+/evTVq1Cg9+OCDys3N1QMPPKAvvvhCixcvLnVXxZ9++kn9+/fXfffdpzvuuEMvvfSSevbsqT/++EN169at4FcQPG63W1deeaVuvPFGuVwupaWl6fbbb9eiRYuO+9m3adMmnXvuubrpppvMtsTExArocYgwEJbWrl1rSDLmzp1rtn3wwQdGZGSksX379lKPi4mJMWbNmlURXQwJPp/PMAzD2Lt3ryHJmD179gmPOffcc43LLrvMvH3kyBGjTp06xoMPPlhe3bS03NxcIzc31zAMwxgxYoQxZMiQEx7z559/GpKMzZs3l3PvQsORI0eM/Px8v7bPP//ckGTs3LmzxGPy8vKMOnXqGA899JDZlpaWZkRERBgzZ84s1/5a1eHDh/1u5+fnG40aNTruf5szZ840EhISyrdjQJjYtm2bERkZaXz00Udm27fffmtIMtxud6nH1a1b15gxY0ZFdNHScnNzjSNHjhiGYRgjR440BgwYcMJjPvnkEyMiIsJIS0sz2x5++GHDbreb372VRXZ2tvl7q1atjH/+858nPIbP+P8p+h25devWE459+/bta1x00UXm7dzcXKN+/fpleu/DSU5OjlFQUODX9sEHHxiSjP3795d6XI8ePSrtvw/Kgks/wtT8+fMVGxurv/3tb2bbkCFDVFBQoB9//PG4x95///3q3LmzLrvsMs2aNauce2ptsbGxJ1Wfk5OjhQsXavDgwWZb1apV1b9/f82fPz/Q3QsJUVFRioo6tclbV111lbp27apRo0bpt99+C3DPQkfVqlUVGen/cV34nhYUFJR4zJo1a7Rnzx6/v8VGjRqpffv2lfZvsei+5ZGRkYqMjCz1PSzk8/l07rnnqmfPnrrxxhvLNCsIqIx+/PFHFRQUaNCgQWZb3759Va1atRN+7jz11FPq3LmzLrnkEn344Yfl3VVLioqKUtWqVU/qmPnz56tdu3Zq1KiR2TZkyBB5vV6tWrUq0F20tGrVqp3ScXzGH1X0O7Lwb7G078i8vDwtWLDAb5wRFRWlAQMGVLpxRnR0tCIiIvzaTjROK/Tuu++qU6dOGjp0qF566SXl5+eXWz9DDUFFmEpLS1OdOnX8pndXr15dcXFxSktLK/W4zp0766677tL06dOVkpKiyy+/XM8++2xFdDks7Nq1S3l5ecUul6lfv/5x33cUN2TIEN1zzz3617/+pVq1aqlbt26aPXt2sLtlCXl5eXrsscfUvXt3v0tjjlX498bfYulef/11bdu2TRdffHGpNZGRkRozZowefPBBPfjgg9qxY4fOPPPMSr/uDFCStLQ0JSQk+E0Tj4qKUu3atY/7uVN4Xfbzzz+vHj166Nprr9UDDzxQEV0OeWlpaSV+zhfeh+PjM750Dz74oOrUqaNevXqVeP/u3bt15MgRxhklOHLkiKZMmaI+ffooKSmp1Lrk5GRNnDhRL7zwgi666CI9/PDDGjVqVAX21NpYoyJM5eXlKTo6ulh7TEyMcnNzSz3uhx9+MBPA7t27yzAM3Xfffbr55ptPOuWvjAoXKy363p/ofYe/Jk2a6MsvvzRv9+rVS3v27NGkSZP8ztRVVjfccIM2bNigRYsWlVpzvL/Fw4cPl2v/QsHPP/+scePG6fHHH1fHjh1Lrbvgggt0ySWXmLf79u2rDh066NFHH630Cw0DRZ3q2OObb74xxx7dunVTtWrVdOutt+q2227jeu0TKOk9LzwzzrjjxPiML9n06dP15ptv6uuvv1ZcXFyJNYx5S2YYhq699lpt27btuOM06eiC58d+9iUnJ2vgwIG666671KFDhwrorbUxoyJMJSUlFVuxuKCgQAcOHJDdbi/1uKJT9Hv37q3s7GyS5TIqTE2Lvvf79u077vsOfyVdKtK7d2+tXbv2hFPowt348eP12Wefad68eWrevHmpdfwtlm7hwoUaMmSI7rjjDt19993HrS36t1ilShX16NFDq1evLs8uAiEpKSlJ+/fvl2EYfu0n+twpaeyRn5+vdevWlUs/w0lJ473C25X9s74s+Iwv7pVXXtEdd9yhDz/8UOedd16pdbVq1ZLEOONYhmHo+uuv15w5czRv3jw1btz4uPUlffZJqtR/f8ciqAhTnTp10v79+/XXX3+ZbUuWLFF+fr46depU5sfZtm2bJCk+Pj7gfQxHiYmJatasmRYvXuzXvnDhwpN631Hctm3bVKNGjWJrNVQmt9xyi9577z3NnTtXZ5111nFrzzrrLFWtWtXvb/Hw4cNasWJFpf5b/PXXXzVw4EDdfPPNevTRR0/pMbZt28ZnIlCCTp06KT8/X8uWLTPb1q9frwMHDjD2KCedOnXS77//7jdTbuHChYqKijrh9wRKVpk/41999VXdfPPNeu+993TRRRcdtzYuLk5nnHEGY97/ZxiGxo4dqy+++ELff/+92rRpc9KPwWdfEUFdyhPlJjc312jZsqVxxRVXGHl5eUZOTo4xYMAAIyUlxa+uffv2xhtvvGEYhmF8/fXXxpdffmmuWvvnn38aLVu2NPr27VvR3bec4+368dhjjxmDBg0ybz/55JNGUlKSsW7dOsMwjq76GxERYSxatKjC+mtVpe36sWnTJqNt27bGggULDMMwjNdff91YsmSJef8vv/xiJCQkGLfcckuF9dVqJkyYYCQlJRkrVqwotebvf/+733s0atQoo02bNsa+ffsMwzCMBx54wIiLizN27dpV3t21pMWLFxvx8fHHXY189uzZRtu2bQ2v12sYhmE89dRTxqZNm8z7Z8yYYURERBivv/56ufcXCDUFBQVGSkqKMXjwYCMnJ8fIy8szhg8fbpxxxhl+O1D07NnTmDZtmmEYhvHDDz8YH330kbmz0datW4327dsbnTp1KraKfmVS2q4fHo/HaNu2rfHjjz8ahmEYe/bsMeLj443777/fMAzD2L9/v9G2bVvjqquuqtD+Wk1pu37wGV+61157zYiOjjY+/vjjUmueeOIJo3///ubtZ555xqhZs6axdu1awzAM46OPPjIiIiLM8VxlcsMNNxh16tQxVq9eXWrNVVddZdxxxx2GYRjGunXrjBdffNHc6Wffvn3GgAEDjAYNGhhZWVkV0merI6gIY6tXrzbatGljJCQkGDVq1DBSUlKMjRs3+tVUqVLFePrppw3DOPrlN3z4cCMhIcFo3LixER0dbYwePdrYu3dvMLpvCYVfaK1btzYkGY0bNzbatm1rPPfcc2bNuHHjjObNm5u38/LyjLFjxxrR0dFG/fr1jbi4OOOVV14JQu+to2fPnkbbtm2NhIQEIy4uzmjbtq1x9tlnm/e73W6/IGjFihVG7969jaSkJKNBgwZG9erVjbvvvrvY1lmVxdKlSw1Jht1uN9q2bev3c2yg07t3b+PCCy80bx88eNAYNGiQERMTY9SuXduoX7++MWfOnCC8Amvo1KmTUaVKlWLvYeGgwTAM47333vPb9vWrr74y2rRpY9SvX99ISkoy6tSpY7zwwgvBegmA5W3cuNHo2LGjUaNGDSMhIcFo06aNsWbNGr+apKQk8x+Ru3btMkaNGmXEx8cbTZo0MaKjo43LL7/8uFuph7PevXub35c1atQw2rZt63eSqXD77mO3kp87d67RoEEDo3bt2kZMTIwxYMAA48CBA0HofXDt2rXL/FyPjo42ateubbRt29YYNWqUWcNnfMlycnKMyMhI82/u2J9jtw6+9dZbjcaNG5u38/PzjZtuusmIjo42GjRoYNSoUcN46aWXgvAKguvnn382JBl16tQp9v79/vvvZl3Xrl2NESNGGIZhGFlZWcaECROMhIQEo2nTpkZMTIzRr18/M/SBYUQYRpELCRF2tm7dqsjISCUnJxe7b82aNapXr57firTZ2dnatWuXkpOTT3lbyXCRnp5e4srFtWvXVt26dSVJO3fuVFZWllq2bOlXc/DgQe3du1fJycmnvGVWuFi3bp256FKhKlWqyOl0Sjq6OvKGDRvUpEkT1ahRw6xJT09XRkaGGjZsWKkv+TjeOjFNmzZV9erVJUmbN29WlSpV5HA4/Gr27NmjzMxMNWnSxG8noMpm48aNys7OLtaemJhobu1X+N9869at/T7/du/eLcMwVK9evQrrLxDK0tLSVFBQUOI12uvWrVOtWrVUp04dsy0nJ0c7duxQw4YNS1yQs7JYv359sYUIIyMjzWnkhd+XjRs39lvkMD8/X1u2bFGNGjXM8Ullk5eXV+K6JtWrV1fTpk0l8RlfGsMwtGbNmhLvq1+/vvnvhF27dikjI0NnnHGGX016err27NlTace8hw4d0ubNm0u8r3nz5oqNjZUkbdq0SdHR0X7bCefl5SktLU316tUz63AUQQUAAAAAALCMynuKEgAAAAAAWA5BBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJZBUAGgQiQmJmrBggXB7kbIWbBggRITE4PdDQAAQkarVq30ySefBLsbIefPP/9UYmKi9uzZE+yuAAQVQHm7+uqrlZiYqHvuuafYfXfffbcSExN19dVXB6FnFSs9PV15eXnHrSkoKFCzZs3022+/VVCvysdXX32lRo0aBeSx8vLylJ6eHpDHAgCEt7vuukuJiYkaM2ZMsfueeeYZJSYmatCgQUHoWcVKT0/XkSNHTljXpUsXffXVVxXQo/Lz22+/KTExUYcOHTrtx8rPz1d6eroKCgoC0DPg9BBUAOXs0KFDql69ul577TW/L82cnBy99tprql69ekC+XMLBokWLlJmZqZSUlGB35bTk5ubq4MGDwe4GAKCS8fl8stlsev/993XgwAG/+15++WVVr15dmZmZQeqdtWzevFm///67evfuHeyunJbCExqGYQS7K0BAEVQAFaB79+6qVauWvvjiC7Pt008/Vd26dXX22WcXq3/ppZfUsWNH1a9fXz179tSsWbP87h87dqwSExNVs2ZNtWzZUn//+9+1a9cuv5rvv/9evXv3VoMGDdStWze9/fbb5n0pKSl64403/OqvuOIKTZo0ybx9/vnn65///KeuueYaJScn64orrpAkHTx4UDfffLOaNWumxo0ba/jw4dq4caPfY/35558aOHCg6tevrx49euizzz4r0/v05ZdfatCgQYqMLPmj6UTPPX78eHXr1k05OTmSJMMwdNFFF+nCCy+UYRj65JNP1Lx5c7399tvq2rWr6tevr4svvljbt2/3ex63261LL71UDRo0kNPp1MSJE4sN7Nxut4YPH65GjRqpXbt2mjZtmgzD0KJFizRy5EgdOnRIiYmJSkxM1P3331/mx/3mm2/UqVMnJScna/DgwVq7dm2Z3jsAACSpWbNmOvvss/XOO++YbT/88IMyMjJKnE3x7rvvqmvXrqpXr566du2q1NRUv/vvv/9+8/usadOmuvzyy7V582a/miVLlqh///5q2LChOnXqpOnTp5v/cB40aJCeeOIJv/pbbrlF1157rXl75MiRGjdunG655RY1adJEffr0kSRlZ2frnnvu0RlnnKHk5GQNHTpUf/zxh99j7dixQ5deeqnq16+vLl266M033yzT+zRr1iz17t1bNWrUKPH+Ez33I488ojZt2vjNevzHP/6hXr16KTc317x089NPP1WPHj1Uv359DRgwQBs2bPB7nq1bt8rlcik5OVktW7bUP/7xD3m93mI1o0ePVuPGjdWqVSs9+uijysvL0+bNm3X++edLkho1aqTExESNHTu2zI+7cOFCnXPOOWrYsKH69eunJUuWlOm9AyoCQQVQQa655hq9/vrr5u3XXnvN70u60MMPP6xp06bpmWee0bJly3Trrbfqqquu0vfff2/WPPfcc9qyZYs2b96szz//XOnp6Ro+fLh5/4EDBzR06FBddNFFWrZsmV544QX98MMP5pdjRkaG+Y/5QllZWcrOzjZvZ2Zm6oknnlD79u21aNEivfbaa8rLy1P//v3l9Xo1a9Ys/fTTT2rSpIl69+5t/oM7NzdXgwYNUlxcnH766Sc9/fTTuuuuu8r0Hn355ZcaOnRoifeV5bkff/xx7d692wxcpk2bph9//FHTp09XRESEjhw5ok2bNmnq1Kl66aWX9P333yszM9MMMiRpy5Yt6tmzp84++2z9+uuv+uijj7RmzRqNHDnS7MvGjRvVvXt3xcfHa86cOfrggw+0ceNG/fbbbzr77LP1yiuvqHr16tqyZYu2bNmif/7zn2V63PXr1+uCCy7QpZdeqgULFmj06NFlfu8AACh07bXXFhtzjB49WlFRUX51L7/8su6++2498MADWrFihR588EHddtttevfdd82awu+wLVu2aM6cOapRo4YGDx5sXs6Zn5+vwYMHq0uXLlqyZIneeOMNbdq0yVyXKjMzU4cPH/Z7Xp/P5zeb9NChQ3rxxRdVs2ZN/fTTT+YJmksvvVS///67PvjgA/3666/q2bOnzj33XO3cudM89uKLL1Z6errmzZunl19+WdOmTdPu3btP+B4db8xRlue+5557FBMToxtvvFGS9MEHH2jGjBl64YUXVLVqVXOmw1133aWnnnpKP/30kxISEjRgwABzhu2+ffvUs2dPNWjQQD/88IO+/vprpaena9iwYeblF3v37tU555yjzMxMffHFF5o1a5Zyc3M1Z84cNW7cWJ9++qkkac2aNdqyZYuee+65Mj3unj17NGDAAJ1zzjlasGCBJk6cqNtuu+2E7xtQYQwA5erSSy81Lr30UmPnzp1GTEyMkZaWZmzevNmIiYkx9uzZY1x44YXGpZdeahiGYWRkZBixsbHGjz/+6PcYEydONGtK4vV6DUnGli1bDMMwjJUrVxqSjIMHD5ZY37x5c+Oll17yaxsyZIhx6623mre7du1qDBs2zK9m5syZRp06dYzc3Fy/9mbNmhkzZswwDMMw3nvvPSMuLs7IyMgw7//qq68MScb8+fNLfQ1bt241qlatWmqfy/LchmEYCxYsMKKiooxnnnnGiImJMd577z3zvvfee8+QZCxdutRs27ZtmxEVFWV89913hmEYxrhx44yLL764WN8kGWlpaYZhGMYNN9xgtG/f3igoKCixr59++qlRvXp1v7ayPO7YsWON3r17+9XceeedBh/VAICyGDdunNGjRw/j0KFDRnx8vLF8+XLjwIEDRmxsrOF2u42xY8caPXr0MAzDMPLz8406deoYH3zwgd9jPPnkk2ZNSXJzc41q1aoZixYtMgzDMHbt2mVIMtatW1difY8ePYwHH3zQr+3aa681RowYYd6+8MILjbPPPtuv5ueffzaio6P9xhOGYRjdu3c3nn76acMwDGP+/PlGlSpVjG3btpn3L1u2zJDk9/1fVGZmphETE2Ns3LixxPvL8tyGYRhr1641YmNjjccee8xITEw0nnvuOfO++fPnG5KMWbNmmW0ZGRlGfHy8OW559NFHje7du/s9R1ZWlhEVFWWOVR566CEjOTnZOHz4cIl9/fXXXw1JRmZmptlWlsd98MEHjdatW/uNZZ555hlDkrFz584SnwuoSFElxxcAAq1evXoaMGCA3nrrLR0+fFhDhw5V7dq1/Wr++OMPZWdn6+KLL1ZERIQMw5BhGDp8+LBatWpl1q1fv14PPfSQfvvtN3m9XnM2wNatW9W4cWO1adNGnTp1Uvfu3eVyudSnTx917dq11EsqStO+fXu/24sXL9b+/ftVr1498zkNw1BmZqZ5CcaqVavUtm1bxcXFmcd17979hM/15ZdfqmfPnkpISCjx/rI8tyT16NFD99xzj26//XaNGjXKvGSlUHR0tDp16mTebtiwoRwOh1avXq3zzjtPixcv1po1a2S32833v/D5Nm7cqEaNGmnJkiXq06ePIiIiTvi6ju3/iR531apVOvfcc/2OK8t7BwDAsWw2m6644gq9/vrrat26tVJSUtS6dWu/ms2bN2vPnj0aO3asbrrpJvN7KScnR/Hx8WZdWlqaHnzwQS1cuFB79uxRQUGBcnJytHXrVnXt2lV169ZV//79NXDgQF199dXq06ePevToUWz2xomUNObIy8tT06ZN/b73s7KyzNpVq1YpOTlZDRs2NI/r2LGjqlWrdtznmjNnjpo1a6ZmzZqVeH9ZnluSnE6nnn76aY0fP179+/fXrbfeWuyxjv0ej4uL05lnnqnVq1ebz7N8+fJiY4P8/Hxt3LhRnTt31pIlS9SjRw/FxMQc9zUV7f+JHnfVqlXq1q2b31iGMQeshKACqEDXXHONbrvtNh05ckQvv/xysftzc3MlST/99JPfl64kValSRdLRSyDOP/989e/fX5988onq1aungoIC1atXz5xKWLVqVS1cuFCzZs3S3LlzNWLECMXFxenbb78tdTcKo4RFmIp+Kebm5qpt27b64YcfitUWDgry8/PNvhYqy2DlRFMwy/LchdauXauIiIgSF7SMjIwsFjBERUWZU1hzc3N1zTXXaPLkycWOLbyOtaCgQFWrVj3RSyrW/xM97qm+dwAAFHXNNddo0KBBatSokSZMmFDs/sIxx8cff1xsEetjT2xceOGFatasmd577z01bNhQ0dHRaty4sd8C4bNnz9bs2bM1Z84cXXfddcrNzdXXX3+tNm3alNi3so456tSpI7fbXWptSd+bERERxdqKKsuY40TPXWjNmjWKiIhQZmamCgoKij13Sd/rx445hg0bpv/+97/Fnqd69eqSTn3McaLHZcwBq2ONCqACDRkyRD6fT4ZhaMCAAcXub9u2rSIjI82tpo79KZyhsHHjRqWlpWny5Mlq166dateurU2bNhV7rOjoaF166aV66aWXtHHjRuXl5em1116TJCUkJBTb8rLowlglOeuss7Ru3Trl5OQU619hWHDGGWfI7Xb7rYGxYsWK4z6uz+fT/PnzjztoKMtzS0evw50zZ45mz56t+fPnFwuEDh8+rHXr1pm39+/fr61bt5ozVs466ywtXLhQCQkJxZ6n8Au8Xbt2+vXXX0vta1RUVLGtvcryuGeccYZWrlzpd9yJ3jsAAErStWtXNWjQQJs3b9bll19e7P5mzZqpevXqWrp0abHvpcIZFenp6VqxYoXuu+8+dezYUXXq1NHBgweVkZHh91iRkZEaMmSI/v3vf2v9+vVq1KiRpk+fLun0xhy7du3S7t27i/UvNjZW0tHvzbS0NO3fv988bsOGDcfdTc0wDM2ePfuEY44TPbd0dDvyV199VV9++aU2bdpU4smIY7/Xc3JytHbtWr8xx5IlS2Sz2Yo9T2E40a5dOy1evLjULUMLxxDH3l+Wx2XMAasjqAAqUFRUlP766y+tX7++xLS/du3auv7663XPPfdo3rx5ys/P14EDB/TWW2/phRdekCTVr19f1apV08yZM2UYhjZt2qTx48f7Pc7PP/+sBx54QNu2bZN0dIHIAwcOmLMpunXrpnfffVder1e5ubl6+umnSzxrUNQVV1yhhg0b6qqrrjIHGVu2bNFdd91lfrmNGDFCUVFRuuuuu5Sdna0dO3bojjvuOO7jzp07V8nJyTrjjDNO67n/+usvTZgwQf/+9781YMAATZs2TbfddpvWr1/v91i33HKL9u3bp6ysLN1yyy1q1KiRuRL67bffrrVr1+r2229Xenq68vLytHz5cl111VXm8RMmTNDixYs1efJkHTp0SJmZmXruuefML/zk5GRlZ2frr7/+Mo8py+OOHz9es2fP1nvvvaf8/HwtWrRI06ZNO+57BwBAaZYuXaqdO3eWuLNFdHS07rzzTj3++OP6/PPPlZubq8zMTH300Ud6/PHHJUnx8fFKSkrSRx99pIKCAu3cuVPXXHON3+OsXbtWEydONE+a7Ny5U7t37/Ybc3zyySdKS0tTfn6+3nzzTc2fP/+EfR8wYIA6deokl8tl7oC1fft2TZ48WfPmzTNrmjRpoptvvlmZmZnav39/sTFRSe9JTk6OevTocVrPvXv3bl1zzTV66KGHNHjwYL3xxhuaPHmyFi1a5PdYd955p7Zv367s7GzdfffdMgzDvCx13LhxSk9P1z/+8Q/t3btXBQUFWrt2rf7+97+bYdCNN96oHTt26LbbblN6erqys7P1+uuvm4usF77Pq1atMp+zLI97ww036Pfff9f06dOVm5urNWvW6LHHHjvh/y9ARSGoACpY9erVzWl3JZk2bZr+8Y9/aOTIkYqNjVXLli31008/6cILL5R0dNDwxhtv6LHHHlO1atXUpUsXXXbZZX6P0blzZ0VHR+ucc85RbGysunfvrr///e8aM2aMJOmhhx5SvXr11KBBA9WvX1/r168v0z7iNptNP/74o2rWrKm2bduqWrVq6tevn+rUqaO2bduar++TTz7RN998o/j4eHXu3NnvH+MlOdEUzLI8d15enkaOHKn+/fubg6gxY8ZoyJAhGjlypDnFtVatWurRo4fatm2rxMRErVq1Sh9//LF5RqJ9+/aaO3euFi1apKSkJCUkJGjs2LF+Z6NSUlL09ddf6+OPP1Z8fLyaNm2qbdu2yel0mo8xZswYtWvXztyetCyP26VLF7344osaP368bDabrrnmGt10000n/P8FAICSxMbGlrr9piQ98MADeuCBB3TzzTfLZrPJ4XBo5syZ5k5iERERSk1N1YwZM1StWjW1bt1anTt39ltjq2XLlmrSpInOP/98xcbGqk2bNurXr595kmLixInq3LmzWrRoodq1a2v27NnmmOZ4qlSpom+//VZt2rRR165dFRMTo+7duys3N1fdunUzaz766COtXbtWiYmJatOmjbp3715sDbBjffnllxo4cOBxLw8py3Nfc801atWqle655x5JR7dhveGGG+RyuZSVlWU+1lVXXaUuXbooPj5eX3/9tT799FNzlmzjxo31448/avv27WrQoIFq1KihK664Qn379jVrmjZtqnnz5mnJkiWqVauWGjRooF9//dXc3r5evXq6++671b9/f3NsUZbHbdasmd599109/vjjql69ui644AJza1PACiKMki4SAxAwPp9P0tF/aJfk0KFDioiIKPH+7OxsvymGx7v/4MGDqlGjRrHrC3NyckpdgCk3N9ecAnjo0CFFRkaaj5eVlaWqVauWeqxhGDpy5MhxF3c6fPiweVlGaf2Tjp4NePvtt9W3b99SH+tEz11QUKCMjIxiz5GXl6esrCzFxcVp5syZGj9+vLmP+PHem8JjpeNfs5mTk6Po6OgSF9Y0DEMZGRmKiYnxuzzlRI9rGIZyc3MVHR1t9j8xMbHUPgAAIB0dF+Tn55caThzv/rKMOapVq6aIiAhlZGQoNja22NoJx/tezcvLU5UqVRQREaHs7GwZhmGOfYqOQYoqXOTzeItkHjvmKK1/0tGTDXfccccJT6Ic77kNw1B6erpsNpuio6PN9vz8fGVmZqpGjRpasGCB+vTpo9zcXEVFRZ1wzJGfn6/8/Hy/xyvqyJEjioqKKnVx9IyMDEVFRfmNKcvyuIV9KxxLJSQknNSC4UB5IKgAEFTLly9Xnz595PV6T3qxqJP1/vvv+wUVAACg8tixY4caN26s3bt3q1atWuX6XD/88INfUAHg5PBfDYCgateunbZv317uIQUAAKjc6tSpo3379vltvwrAmphRAaDSyM3NVXZ2NgMUAABQrrh0Ezg9BBUAAAAAAMAy2PUDAAAAAABYBkEFAAAAAACwDIIKAAAAAABgGQQVAAAAAADAMggqAAAAAACAZRBUAAAAAAAAyyCoAAAAAAAAlkFQAQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsAyCCgAAAAAAYBkEFQAAAAAAwDIIKgAAAAAAgGUQVAAAAAAAAMsgqAAAAAAAAJZBUAEAAAAAACyDoAIAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZBBUAAAAAAMAyCCoAAAAAAIBlEFQAAAAAAADLIKgAAAAAAACWQVABoEL9/vvvat26tX744YdgdwUAAOCkPPvsszrzzDOD3Q0g7EUFuwMAwsPmzZs1aNCgEu+LiYnR77//LknKzs7W+vXrlZWVddLPsWTJEo0ePVqvvvqqevbseVr9BQAA1jdr1izdeeedZaqdN2+eGjZsWOr9U6ZM0fvvv2+OSU7Fnj17tG7dulM+HkDZEFQACIicnBytX79eN9xwg2699Va/+yIjAzN5y+fznXLIAQAAQk/v3r312Wef+bW1b99eXbp00X//+1+/9jp16hz3sXbv3q3169cHuosAygFBBYCASkpKUuvWrYPdDQAAEAbi4+MVHx/v1xYREaHq1asz3gDCGGtUAAi67du3q3Xr1uZP+/btdfHFF+vTTz81az799FONHj1akvSPf/zDrH377beD1W0AAGARf/75p2688UadffbZSklJ0dVXX62VK1ea948fP15vv/22jhw54jfm2L59uyTp7bffNtucTqe6dOmisWPHau3atUF6RUDlxowKAEFXp04dv2mdBw4c0Ndff63LLrtMr776qsaMGaO+fftq8uTJuvrqq/Xggw+aa1TUrVs3SL0GAABWsHTpUvXr108dO3bU5MmTFRMTo6lTp6pr16769NNPNXjwYE2aNEmZmZl6//33/cYchZeLDBs2TF26dJEkFRQUKC0tTdOnT1f37t31xx9/qHHjxsF4aUClRVABIKD+85//6KOPPvJru/fee83ZECWpWrVqsemb3bt3165du/T4449rzJgxSkhIkMPhkCQ1atSI6Z4AAECSNG7cOMXFxenbb79VtWrVJB1d26JTp04aO3asNm/erIYNG6pmzZqKiIgocQxRs2ZN1axZ07zdpk0bnXfeeWrUqJGef/55Pf300xX2egAQVAAIsOHDhxdbTLMssx6WLFmiV155RatWrVJGRoYMw9C+ffu0f/9+HTlyRNHR0eXVZQAAEKL27NmjpUuXauLEiWZIIR1dx+Lqq6/WhAkTtHLlSnXu3Pm4j1NQUKB33nlHH3/8sbZs2aKcnBwZhqGDBw+yywcQBAQVAALqVBbT/OabbzR06FCNGDFCjzzyiBo0aKCqVatq+vTpevHFFwkqAABAiXbu3ClJ5qzLYxVerrFjx44TPs4tt9yi//73v3rkkUd01113KTExUZGRkbrooouUnZ0d2E4DOCGCCgBBN23aNDVu3FipqamKiIgw2w8fPhzEXgEAAKuLi4uTdHR9q6L27dvnV1OaQ4cO6eWXX9bNN9+su+++2+++vXv3qlGjRgHqLYCyYtcPAEGXmZmpevXq+YUUOTk5mj17tl9d4ZTOvLy8Cu0fAACwpqZNm6pevXqaN29esfvmzp2ratWqKSUlRdLRcUR+fn6xusOHDysvL0/169f3a//pp5+0f//+8uk4gOMiqAAQdL1799bSpUv1008/SZJ8Pp+uu+46NWvWzK+uadOmioyM9NtuDAAAVF4RERG6//779csvv+jJJ59UQUGBDMPQ22+/rQ8//FATJ05UQkKCJKl58+bKy8vT6tWr/R4jKSlJbdu21dtvv629e/dKkjZs2KAHHnhAzZs3r/DXBICgAoAF/POf/9Sll16qPn36qFGjRmratKm6d++ugQMH+tXVrVtX9957rx555BE1btxYrVu31ttvvx2kXgMAACu46aab9Pzzz2vq1KlKSEhQrVq1NGHCBD300EOaPHmyWedyudSzZ0916tRJLVu2VOvWrbV9+3ZJ0nvvvafIyEg1atRIjRo10sUXX6xnnnlGiYmJQXpVQOUWYRiGEexOAAh9R44c0aZNm2S322W320utO3z4sLZs2aLk5GRVr17d7z6fz6e9e/eai2nu27dPe/fuVatWrfwuC8nOztaOHTuUm5urunXr+m0nBgAAwtuGDRtks9mKrR1hGIa2b9+u/Px8NWrUSFWqVCnx+H379mnfvn0qKChQ8+bNVbVqVfO+vXv3+l0GsnXrVkVGRio5Odm8f//+/WrVqlU5vToAEkEFAAAAAACwEC79AAAAAAAAlkFQAQAAAAAALIOgAgAAAAAAWAZBBQAAAAAAsIyTCioOHz6sGTNm6JxzzlG1atX03nvvlVj34Ycf6qyzzlJCQoI6d+6sOXPmBLUGAAAAAACEhqiTKX755Ze1dOlSPfXUU+rVq5fy8/OL1cydO1cjR47UCy+8oGHDhumNN97QsGHDtGTJErVv377Ca06koKBAO3bsUFxcnN/2hwAAVBaGYSgzM1MNGjRQZCSTLcsLYw4AQGVX1jHHKW9PGhERoRkzZsjlcvm1DxgwQNWqVdPnn39utnXu3Flt27bVW2+9VeE1J7Jt2zZzX2QAACqztLQ0NWrUKNjdCFuMOQAAOOpEY46TmlFxIoZh6JdfftHkyZP92vv166ePP/64wmvKIi4uTtLRNyo+Pr7MxwEn7dAhqUGDo7/v2CFVrx7c/gDA/8vIyFBycrL5nYjywZgDOEmMnYCwU9YxR0CDiszMTB06dEh16tTxa69Tp4527dpV4TUlycnJUU5Ojl+fJSk+Pp5BA8pXlSr/+z0+ni9bAJbD5Qjlq/D9ZcwBlBFjJyBsnWjMUS4Xoha91iQyMlJFrzCpyJpjTZkyRQkJCeYPUzABAAAAALCOgAYVcXFxstls2rt3r1/7nj17VLdu3QqvKcmkSZOUnp5u/qSlpZ3aiwUAAAAAAAEX0KAiIiJCXbt21Y8//ujXPn/+fHXv3r3Ca0oSExNjTrlk6iUAAAAAANYS8Es/Jk6cqC+++EIffvihfD6fXnrpJf3222+65ZZbglIDAAAAAABCx0kFFfPnz1e1atVUrVo1SdI111yjatWq6cYbbzRrhg0bpueff1633367atSooWeeeUbvv/++unbtGpQaAAAAAAAQOiKM4608WURBQYGOHDlSrD0qKkpRUcU3ECkoKCi22GUwa0qSkZGhhIQEpaencxkIytehQ1KNGkd/z8pi5WoAlsF3YcXgfQZOEmMnIOyU9bvwpLYnjYyMNGdTlLXeSjUAAAAAAMDaTiqoAAAAqCi5ubn6/PPP9csvv+iCCy5Qnz59itVkZ2frnXfe0Z9//qnGjRtr5MiRSkhICFoNAAA4fUxDAAAAlvPzzz+refPmevfdd/Xaa69p2bJlxWqysrJ0zjnn6Pnnn1e1atU0Y8YMpaSk+G1dXpE1AAAgMAgqAACA5TRs2FCLFi3SJ598IpvNVmLNtGnTtGPHDv344496+OGH9f333ysyMlKTJ08OSg0AAAgMggoAAGA5zZo1U4MGDY5b8/nnn2vYsGHm5RexsbEaPny4Pvvss6DUAACAwCCoAAAAIWnDhg1q3ry5X1vz5s3l8XiUk5NT4TVF5eTkKCMjw+8HAACcGEEFAAAIST6fT3FxcX5thVudZWdnV3hNUVOmTFFCQoL5k5ycfHIvEACASoqgAgAAhKS4uDgdPHjQr+3AgQOKiIhQ9erVK7ymqEmTJik9Pd38SUtLO7UXCgBAJcP2pAAAICS1adNG69at82tbt26dWrZsqapVq1Z4TVExMTGKiYk59RcIAEAlxYwKAAAQki6//HJ9+eWX2r59uyRp//79+uCDD3T55ZcHpQYAAARGhGEYRrA7EUwZGRlKSEhQenq6ea0pUC4OHZJq1Dj6e1aWVMpUYQCoaFb8Lty1a5eeeOIJSdJ///tfderUSSkpKWrbtq3+8Y9/SJJyc3M1bNgwrV69Wv369dOCBQtkt9s1b9481fj/z9uKrDkRK77PgKUxdgLCTlm/CwkqGDSgovBlC8CirPhduG/fPs2YMaNYe7NmzXTBBReYtw3D0Lx58/Tnn3/K4XBowIABioryv7K1ImuOx4rvM2BpjJ2AsENQUUYMGnAyPB6PvF6vedtut8vhcJTtYL5sAVgU34UVg/cZlV3RcVRJ/MZWjJ2AsFPW70IW0wTKyOPxyOl0yufzmW02m01ut7vsYQUAAEAlVNI4qiSMrQBIBBVAmXm9Xvl8PqWmpsrpdMrtdsvlcsnr9Zpfpqc14wIAACBMFR1HlaRwbPXzzz/L6XQqMjtbHSq2mwAsgqACOElOp1MpKSnmbbfbLUnau3evLrnkEmZcAAAAlKLoOOpYdrtdNptNLpdLkmSTdKgC+wbAOggqgFNU9MtUOhpMfPPNN6pdu3aJMy4AAABQMofDIbfbbc5O3bBihXTddUHuFYBgIKgATlHRL1OJSz0AAABOh8PhMMdSkdnZQe4NgGAhqABOw7FfpgAAAACA0xcZ7A4AAAAAAAAUYkYFcBzH7uJRuGgmAAAAAKD8EFQAxzg2mChtFw+73X5Sj1kYcLDFFgAAAACcGEEF8P88Ho+cTmexYKJwFw/p5BbLZIstAAAAADh5BBXA//N6vfL5fEpNTZXT6ZR0ert4sMUWAABAYKxcuVIFsbHssAZUEgQVQBFOp1MpKSkBeSy22AIAADg1SUlJ5u89evaUT0dnu7rdbsIKIMyx6wcAAAAAy0lOTjZ//2XBAqWmpsrn85mzVQGEL2ZUAAAAALC0Dh06qCA2NtjdAFBBCCoAAAAAhIxjt4xnzQogPBFUAAAAALC8ojuqSaxZAYQrggoAAAAAlld0RzW32y2XyyWv10tQAYQZggoAAAAAp83j8fgtdFkel2Ucu6MagPBFUAEEQeFe4BLXVgIAgNDn8XjkdDrl8/nMNpvNpk8++US1a9f2W1cCAE6EoAKoICXtBS5xbSUAAAh9Xq9XPp9Pqampcjqd2rt3ry655BINHDjQrLHZbLLb7UHsJYBQQVABVJCie4EXxMZybSUAAAgrTqdTKSkpkuS3noTELFIAZUdQAQRBhw4dpOrVg90NAACAcsN6EgBOVWSwOwAAAAAAAFCIoAIAAAAAAFgGQQUAAAAAALAMggoAAAAAAGAZLKYJAAAAIGS53W7zd3YWAcIDQQUAAACAkGO322Wz2eRyucw2m80mt9tNWAGEOIIKAAAAACHH4XDI7XbL6/VKOjqzwuVyyev1ElQAIY6gArAApiwCAACcPIfDwbgJCEMEFajUPB6PXwpf0ZiyCAAAAAD+CCpQaXk8HjmdTvl8PrPNZrPJbrdXWB+YsggAAAAA/ggqUGl5vV75fD6lpqbK6XRKCs5lF0xZBAAAAID/IahAped0OpWSkhLsbgAAACAAWPsLCH0EFQAAAABOWrDX+iqKtb+A8EFQAQAAAOCkWGGtr6JY+wsIHwQVAAAAAE6KVdb6Koq1v4DwQFABAAAA4JSw1heA8hAZ7A4AAAAAAAAUIqgAAAAAAACWwaUfqFSstjo1AAAAAMAfQQUqDSuuTg0AAAAA8EdQgUrDqqtTAwAAAAD+h6AClQ6rUwMAAACAdRFUABZ07PoZzPoAAAAAUJkQVAAWYrfbZbPZ5HK5zDabzSa3201YAQAAAKBSIKgALMThcMjtdvvtTOJyueT1egkqAAAAAFQKBBWAxTgcDkIJAAAAAJVWZLA7AAAAAAAAUIgZFQAAAABOyOPx+F2eCgDlhaACAAAAwHF5PB45nU75fD6zzWazyW63B7FXAMIVQQUAAACA4/J6vfL5fEpNTZXT6ZTEFuoAyg9BBQAAAIAycTqdSklJCXY3AIQ5FtMEAAAAAACWQVABAAAAAAAsg6ACAAAAAABYBkEFAAAAAACwDIIKAAAAAABgGez6AQAAACBsud1u83e2VAVCA0EFAAAAgLBjt9tls9nkcrnMNpvNJrfbTVgBWBxBBQAAAICw43A45Ha75fV6JR2dWeFyueT1egkqAIsjqABCAFMWAQAATp7D4WDcBIQgggqENY/H45eihxqmLAIAAACobAgqELY8Ho+cTqd8Pp/ZZrPZZLfbg9irk8OURQAAAACVDUEFwpbX65XP51NqaqqcTqek0LxsgimLAAAAACqTyEA/oGEYeu211zR48GB17txZF1xwgT744INidb/++qsuueQSderUSVdccYVWr15dbjWo3JxOp1JSUpSSksI/+AEAAADA4gIeVDzxxBOaMGGCRowYof/85z/q37+/Ro0apTfeeMOsWbFihfr06aMWLVpo6tSpqlGjhnr+X3v3HR5Vlbhx/E2BJEOoGQxShl6GIhhYQBErKIoVVhAcEWEBFXZxVRTWgoou6y6urth2UVGMCig1FFlESqhC6DgUaUMVhg4DAZLz+4Mfd5kUCDDJzCTfz/PM82TuPXfm3Em5J+895aabtG3btoCXAQAAAAAA4SPgQcW0adPUqVMnPf7442rWrJn69eunNm3aaPr06VaZN998UzfccIP+/ve/q3Xr1hoxYoQSExM1bNiwgJcBAAAAAADhI+BBRatWrbRkyRIdOXJEkrRnzx6tWrVKrVu3tsrMnj1b99xzj/U8IiJC99xzj2bPnh3wMgAAAAAAIHwEfDLNv/71rzp16pQqVqyoihUrateuXXrppZf0xz/+UZJ04sQJHTp0SNdee63fcddee6127NgR0DI5SU9PV3p6uvX86NGjV36yAAAAAAAgoAIeVHz66acaNWqU3n//fTVq1EiLFy/Wyy+/rIYNG+qBBx7Q2bNnJUnFixf3Oy4mJkZnzpyRpICVycnQoUP1+uuvX8UZAsHndrutr8NxJRMAAAAAyE3Ag4qBAwfqueeeU8+ePSVJzZs318aNG/WXv/xFDzzwgEqWLKnixYvrwIEDfscdOHBAdrtdkgJWJieDBg3Ss88+az0/evSoqlSpcuUnDBQgu90um80ml8tlbbPZbHK73YQVAAAAAAqFgM5RkZGRIZ/Ply0oKF++vDXEIjIyUk2aNNGSJUv8yixcuFBNmzYNaJmcxMTEqFSpUn4PIFw4HA653W6lpaUpLS1NycnJ8vl88nq9wa4aAAAAAAREQIOKqKgo3XLLLRoxYoTV02HXrl0aNWqUbr/9dqtcnz599P333ystLU2S9N///ldz5sxR7969A14GKGwcDoeSkpKUlJQkp9MZ7OoAAIBCyuPxaPny5Vq+fLnfsFMAyG/5MkdFjx49VLlyZVWsWFE7d+7UPffco3fffdcq06NHD23YsEGtWrVSQkKCDh06pKFDh/qt4BGoMgAAAAAuj8fjkdPplM/ns7bZbLaLDrEGgEAJ+PKkVapU0cyZM7V//37997//1eHDhzVhwgSVK1fOr9zbb7+t/fv3a86cOfJ6vRowYEC21wpUGQAAUDh99NFHuvHGG1WtWjU1b95cb7/9tjXh9nkLFy5UmzZtVLVqVd18882aPn16ttcJVBmgsPB6vfL5fEpOTraGnDInFoCCEvCg4rz4+HjVrFlTcXFxuZYpWbKkateuLZvNlu9lAABA4TJ8+HA9//zzeuaZZzRnzhy9+uqrevvtt/Xaa69ZZdxut9q0aaNmzZrphx9+0D333KP7779fCxYsCHgZoDByOp3WkFNCCgAFJeBDPwAAAArC3Llzdfvtt6tTp06SpGrVqumhhx7S3LlzrTLvvPOO6tatq7/97W+Szv3T9dNPP+lvf/ubUlJSAloGQHhgmXcg9OVbjwoAAID81K5dOy1dulSbNm2SJO3cuVNz5871m6tq3rx5atu2rd9xd911l1JTUwNeBkBou3CZ96ZNm6pp06ZyOp3yeDzBrhqALOhRAQAAwtIf/vAH7d+/X06nUyVLltSRI0f05z//WYMGDbLK7Nq1S4mJiX7HJSYm6siRIzpx4oRKlCgRsDJZpaenKz093Xp+fql2AMFxfpn388u6u91uuVwueb1eelUAIYagAgAAhKUvvvhCQ4cO1ahRo9SyZUutWbNGvXr1UsWKFfXcc89Jkowxio72b+6cf56RkRHQMlkNHTpUr7/++tWcIoAAczgchBJAGGDoBwAACEtvvfWWevXqpa5du6pGjRp64IEH9Pzzz+uNN96wypQvX966e3re/v37FRsbq5IlSwa0TFaDBg3SkSNHrMeOHTuu+pwBACgK6FEBAADC0pkzZ7INuYiPj9eZM2eUmZmpyMhINW/eXPPnz/crM2/ePP3ud79TRESEJAWsTFYxMTGKiYm5qnMEAKAookcFCg2Px6Ply5dbjwtndAYAFD7t2rXT559/rl9++UWStG3bNg0fPlx33XWXIiPPNXH69u2r1NRUjR49WsYYzZgxQ5MnT1a/fv2s1wlUGQAAEBj0qECh4PF45HQ65fP5/LbbbDbZ7fYg1QoAkJ+GDRumjIwMNW/eXFFRUTpz5oweeOABvf/++1aZW2+9Vf/+97/1pz/9SU888YRiYmL05ptvWkuaBrIMAAAIDIIKFAper1c+n0/JyclyOp3WdtbGBoDCKz4+XiNGjNCIESN05MgRlSpVKsdhGD179tQTTzxhlYmKisq3MgAA4OoRVKBQcTqdSkpKCnY1AAAFrHTp0hfdHxkZqbJlyxZIGQAAcHWYowIAAAAAAIQMelQAAAAAkMfjsZbhZVJyAMFEUAEAAAAUcTlNTM6k5ACChaACAAAAKOJympicSckBBAtBBQAAAABJTEwOIDQwmSYAAAAAAAgZBBUAAAAAACBkEFQAAAAAAICQQVABAAAAAABCBpNpAoXAhWudM0M3AAAAgHBGUAGEMbvdLpvNJpfLZW2z2Wxyu92EFQAAAADCEkEFEMYcDofcbre8Xq+kcz0rXC6XvF4vQQUAAACAsERQAYQ5h8NBKAEAAHCFGEILhB6CCgAAAABFDkNogdBFUAEAAACgyGEILRC6CCoAAAAAFEkMoQVCU2SwKwAAAAAAAHAeQQUAAAAAAAgZBBUAAAAAACBkMEcFwpbH4/Gb/AgAAAAAEP4IKhCWPB6PnE6nfD6ftc1ms8lutwexVgAAAACAq0VQgbDk9Xrl8/mUnJwsp9Mp6dxa2MzaDAAAkDf0TgUQqggqENacTqeSkpKCXQ0AAICwQu9UAKGMoAIAAAAoYuidCiCUEVQAAAAARRS9UwGEIpYnBQAAAAAAIYOgAgAAAAAAhAyCCgAAAAAAEDIIKgAAAAAAQMggqAAAAAAAACGDVT+AQsjtdltfs9QYAAAAgHBCUAEUIna7XTabTS6Xy9pms9nkdrsJKwAAAACEBYIKoBBxOBxyu93yer2SzvWscLlc8nq9BBUAAAAAwgJBBVDIOBwOQgkAAAAAYYvJNAEAAAAAQMggqAAAAAAAACGDoR8AAAAA8P8uXD1NYgU1IBgIKgAAAAAUeTmtniaxghoQDAQVAAAAAIq8rKunSaygBgQLQQUAAAAAiNXTgFDBZJoAAAAAACBkEFQAAAAAAICQQVABAAAAAABCBkEFAAAAAAAIGQQVAAAAAAAgZLDqB8KGx+Oxlotyu91Brg0AAAAAID8QVCAseDweOZ1O+Xw+a5vNZpPdbg9irQAAAAAAgUZQgbDg9Xrl8/mUnJwsp9MpSbLb7axzDQAAAACFDEEFworT6VRSUlKwqwEAAAAAyCdMpgkAAAAAAEIGPSqAIuDCyUcZMgMAAAAglBFUAIWY3W6XzWaTy+WyttlsNrndbsIKAAAAACGJoAIoxBwOh9xut9+yri6XS16vl6ACAAAAQEgiqAAKOYfDQSgBAABwFRhGCxQsggoAAAAAyAHDaIHgIKgAAAAAigCPx+M3HBSXxjBaIDgIKgAAAIBCzuPxyOl0yufzWdtsNpvsdnsQaxUeGEYLFDyCCgAAAKCQ83q98vl8Sk5OltPplMRcCwBCF0EFAAAAUEQ4nU4lJSUFuxoAcFGRwa4AAAAAAADAeQQVAAAAAAAgZBBUAAAAAACAkEFQAQAAAAAAQgZBBQAAAAAACBkEFQAAAAAAIGQQVAAAAAAAgJBBUAEAAAAAAEJGdH698OTJkzVt2jRFRESoS5cuuvnmm/3279ixQx9//LG2b9+u2rVrq1+/frLb7flSBoA/t9ttfW232+VwOIJYGwAAAAD4n4D3qDDGqGvXrnryySdVq1YttWzZUkOGDNEPP/xgldm2bZuSkpLkdrt16623KjU1Vc2aNdOBAwcCXgbA/9jtdtlsNrlcLjVt2lRNmzaV0+mUx+MJdtUAAAAAQFI+9KgYOXKkxo0bp9WrV6tu3bqSpG7duunQoUNWmSFDhqhy5coaN26cIiMj9dhjj6lWrVp699139eabbwa0DMKXx+OR1+uV5N8DAFfO4XDI7Xb7fa4ul0ter5deFQAAAABCQsB7VHz66ae69957rZBCkiIiIlSuXDnr+bRp09ShQwdFRp57+9jYWN13332aNm1awMsgPHk8HjmdTuuuv8vlks1mY1hPADgcDiUlJSkpKUlOpzPY1QEAAAAAPwEPKlatWqUWLVro888/12OPPabnnntOCxcutPafPHlSe/fuzXb31uFwaMuWLQEtk5P09HQdPXrU74HQ4/V65fP5lJycrLS0NKWlpcntdnPXHwAAAAAKuYAGFcYY+Xw+DR8+XDNmzFCbNm0UHR2t2267TZ9//rmkc0GBJNlsNr9j4+PjderUqYCWycnQoUNVunRp61GlSpUrPV0UAKfTad39J6QAAAAAgMIvoHNUREREqHTp0kpMTNSYMWOs7WfPntXrr7+uHj16KD4+XlFRUX5zVkjSgQMHVKZMGUkKWJmcDBo0SM8++6z1/OjRo4QVAAAAAACEiIBPptm4cWNVrFjRb1vt2rW1d+9eGWMUHR0tp9Op1atX+5VZvXq1rrvuunOVClCZnMTExCgmJuaKzw8AAISeX3/9Venp6XI6ndbcVRc6cOCAtmzZoipVqqhChQo5vkagygAAgKsT8DkqunfvrtmzZ2v//v2SzvWmmDBhgm644QZFRERIklwul8aOHaudO3dKktatW6fp06fL5XJZrxOoMgAAoPBavny5GjVqpJtuukmPP/64mjRporVr1/qVeemll1SpUiV1795d1apVU58+fZSZmZkvZQAAwNULeI+K7t27a9GiRXI6nWrevLk2bNig2NhYTZw40Srz5z//WQsXLlTjxo3VpEkTLV26VF26dNFjjz0W8DIAAKBw2rNnj9q2bauuXbvqvffeU1RUlH799Vdt2bJFDRs2lCSNGzdOw4YN09y5c9WyZUu53W61aNFCjRs31tNPPx3QMgCKDrfbbX1tt9uZSw0IsAhjjMmPF960aZM2bNiga6+9Vo0bN1Z0dPZMZNWqVfJ4PKpVq1auyyQGqkxujh49qtKlS+vIkSMqVarUZR2L/LN8+XI1bdpUaWlpSkpKCnZ1AuPECSk+/tzXx49LJUoEtz4qpJ8zgMsWrtfCF154QV999ZW2b9+u4sWL51imffv2ioyMVEpKirXtiSee0Jo1a7Rs2bKAlrmUcP2cUTiE5TU/BNtOHo9HTqdTPp/P2maz2VidDsijvF4LA96j4rzatWurdu3aFy3TuHFjNW7cuEDKAACAwmXWrFm66667FBERoeXLl6ts2bKqWrWq3xwVK1asUO/evf2Oa968ub7++mtlZmYqMjIyYGWySk9Pt1Ypk8SS6ChwHo9HXq9Xkn8PAFw5h8Mht9vt97m6XC55vV6CCiCA8i2oAAAAyE+7d+9W9erVVb9+fZUoUUJ79uxR+fLl9e2336pRo0aSpIMHDyohIcHvuISEBJ05c0bHjh1T6dKlA1Ymq6FDh+r1118P8FkDeZPbnX+73R7EWhUODoeDUALIZwGfTBMAAKAgFCtWTFOnTtXYsWO1cuVK7dixQzVq1FDXrl39ypw6dcrvuJMnT0qSNVwkUGWyGjRokI4cOWI9duzYcaWnClw2r9crn8+n5ORkpaWlKS0tjeEJAMIGPSoAAEBYqlatmqpXr67rr79e0rnA4IknnlCHDh106NAhayjIrl27/I7btWuX7Ha74uLiJClgZbJiSXSEAqfTGT5zUgDA/6NHBQAACEt33XWX9uzZowvnBd+5c6diY2NVsmRJSVLbtm01depUZWRkWGUmTZqktm3bWs8DVQYAAAQGQQUAAAhL/fr109mzZ/X4449rxowZ+uSTT/Taa6/pueees1Ybe/7553X48GF17dpVKSkp6tWrl3755Re9/PLL1usEqgwAAAgMggoAABCWSpcurUWLFql8+fL6xz/+oblz5+qjjz7Sm2++aZWpVKmSlixZorJly+q9997T2bNntXDhQtWvXz/gZQAAQGAwRwUAAAhbiYmJeueddy5apkaNGvrkk08KpAwAALh69KgAAAAAAAAhg6ACAAAAAACEDIIKAAAAAAAQMpijAoDcbrf1td1ul8PhCGJtAAAAABRlBBVAEWa322Wz2eRyuaxtNptNbrebsAIAAABAUBBUAEWYw+GQ2+2W1+uVdK5nhcvlktfrJagAAAAAEBQEFUAR53A4CCUAAAAAhAwm0wQAAAAAACGDoAIAAAAAAIQMggoAAAAAABAyCCoAAAAAAEDIIKgAAAAAAAAhg6ACAAAAAACEDIIKAAAAAAAQMqKDXQEAAAAACGdut9v62m63y+FwBLE2QPgjqAAAAACAK2C322Wz2eRyuaxtNptNbrebsAK4CgQVAAAAQJjzeDzyer3W8wvv8CP/OBwOud1u67N3u91yuVzyer0EFcBVIKhAyLjwAsvFFQAAIG88Ho+cTqd8Pp/fdpvNJrvdHqRaFR0Oh4NQAggwggqEhJwusFxcAQAALs3r9crn8yk5OVlOp9PazlwJAMIVQQVCQk4XWC6uAAAAeed0OpWUlBTsagDAVSOoQEjhAgsAAAAARVtksCsAAAAAAABwHj0qAGSTdTJThuEAAAAAKCgEFQAsOa0FLrEeOAAAAICCQ1ABwJJ1LXCJ9cABAAAAFCyCCgB+WAscAAAAQDAxmSYAAAAAAAgZBBUAAAAAACBkEFQAAAAAAICQQVABAAAAAABCBkEFAAAAAAAIGQQVAAAAAAAgZBBUAAAAAACAkEFQAQAAAAAAQkZ0sCsAAAAAAIWJ2+22vrbb7XI4HEGsDRB+CCoAAAAAIADsdrtsNptcLpe1zWazye12E1YAl4GgAgAAAAACwOFwyO12y+v1SjrXs8Llcsnr9RJUAJeBoAIAAAAAAsThcBBKAFeJyTQBAAAAAEDIoEcFgsbj8fh1iwMAAAAAgKACQeHxeOR0OuXz+axtNptNdrs9iLUCAAAID9zwAVCYEVQgKLxer3w+n5KTk+V0OiWxdFOoY5ktAABCAzd8ABR2BBUIKqfTqaSkpGBXAxfBMlsAAIQWbvgAKOwIKgBcFMtsAQAQmrjhA6CwIqgAcEksswUAAACgoLA8KQAAAAAACBkEFQAAAAAAIGQQVAAAAAAAgJBBUAEAAAAAAEIGQQUAAAAAAAgZBBUAAAAAACBkEFQAAAAAAICQQVABAAAAAABCBkEFAAAAAAAIGQQVAAAAAAAgZBBUAAAAAACAkEFQAQAAAAAAQgZBBQAAAAAACBnRwa4AgPDkdrutr+12uxwORxBrAwAAAKCwIKgAcFnsdrtsNptcLpe1zWazye12E1YAAAAAuGoEFQAui8PhkNvtltfrlXSuZ4XL5ZLX6yWoAAAAAHDVCCoAXDaHw0EoAQBAAfJ4PH43CQCgMCOoAAAAAEKYx+OR0+mUz+ezttlsNtnt9iDWCpeDub2Ay0NQAQAAAIQwr9crn8+n5ORkOZ1OSfyzGy6Y2wu4MgQVAAAAQBhwOp1KSkoKdjVwGZjbC7gyBBUAAAAAkE+Y2wu4fJHBrgAAAAAAAMB5BBUAAAAAACBkEFQAAAAAAICQQVABAAAAAABCRr5Oprl582aNGDFCjRo10qOPPuq37/jx4xo7dqy2b9+u2rVrq1OnTipevHi+lAEAAAAAAOEh33pUnD59Wo888ohGjBihCRMm+O07cOCAmjVrpo8++kinTp3SX//6V910003y+XwBL4PQ4fF4tHz5ci1fvlxutzvY1QEAAAAAhKB861ExcOBANWrUSImJidn2vfXWW8rIyNC8efNks9n0wgsvqE6dOvrwww81YMCAgJZBaPB4PHI6nX4hks1mk91uD2KtAAAAAAChJl96VEybNk0pKSl6//33c9w/fvx4Pfzww7LZbJKkhIQE3XfffRo/fnzAyyA0eL1e+Xw+JScnKy0tTWlpaXK73awpDQAAAADwE/AeFXv27NEf/vAHjR8/XvHx8dn2p6ena/v27apVq5bf9lq1aiklJSWgZXKSnp6u9PR06/nRo0cv7wRxVZxOp5KSkoJdDQAAAABAiApoUJGZmalHH31Uffr0UcuWLXMsc77rf6lSpfy2ly5dWidOnAhomZwMHTpUr7/+el5PCQAAhIE9e/bo73//u6pXr64//elPfvtOnjypr7/+Wps2bVLVqlX16KOPqnTp0vlSBgAAXL2ADv2YMWOGFixYoOPHj2vgwIEaOHCg1q9fr7Vr12rgwIE6dOiQSpQooYiICB05csTv2MOHD1s9MAJVJieDBg3SkSNHrMeOHTsCceoAACBIMjMz5XK59PXXX2vs2LF++44fP64bb7xRH3zwgWJjY/XVV18pKSlJ+/fvD3gZAAAQGAENKmrVqqXXX39dCQkJKlOmjMqUKaNixYqpWLFiKlOmjCIjI1W8eHHVqFFDGzZs8Dt2w4YNcjqdkhSwMjmJiYlRqVKl/B4AACB8DR06VCVLltT999+fbd/777+v3bt3a+7cuXr99df1008/KTIyUm+++WbAywAAgMAIaFBRu3ZtqyfF+UfNmjVVt25dDRw40Ooe+fDDD2vs2LHW/BB79uzRlClT9PDDD1uvFagyAAqG2+22lp/1eDzBrg6AImLhwoX65JNP9Omnn+a4f9KkSbrvvvusNkhcXJwefvhhTZw4MeBlAABAYOTLqh+XMnDgQJUtW1YtWrRQnz59dMMNNygpKUl9+vQJeBkA+ctut8tms8nlcqlp06Zq2rSpnE4nYQWAfHfo0CF17dpVI0aMyHW5640bN6pmzZp+22rWrCmPx2NNrh2oMlmlp6fr6NGjfg8AAHBpAV/1I6vu3bsrOtr/bUqXLq0lS5YoJSVFHo9H9913n+6++25FRUUFvAyA/OVwOOR2u+X1eiWd61nhcrnk9XpZfhZAvvrDH/6gBx54QO3atcu1jM/nU8mSJf22nR/2efLkScXExASsTFZM4A0AwJXJ96Di97//fY7bixcvro4dO1702ECVAZC/HA4HoQSAArV8+XJNmDBBvXv31jPPPCNJWrx4sQ4dOqRnnnlGf/7zn1W1alWVLFlShw8f9jv20KFDioiIUIkSJSQpYGWyGjRokJ599lnr+dGjR1WlSpUrP2kAAIqIfA8qAAAAAq1ixYr65z//6betZMmS8vl8qlatmtXDoX79+lq/fr1fufXr16t27doqVqxYQMtkFRMTk2NPCwAAcHEEFQAAIOxUqFDB6klx3tq1axUVFeW3vVOnTnr55Ze1a9cuVapUSQcPHtSYMWPUo0ePgJcBgLxyu93W13a7nZ6pQBYEFQAAoNDq06ePpkyZohYtWuiOO+7Q/PnzVblyZb344osBLwMEksfj8Zv/CYXDhZOQn2ez2eR2uwkrgAsQVAAAgELhkUce0aFDh/y2FStWTNOnT9esWbO0adMmderUSXfddZffRN+BKgMEisfjkdPplM/ns7bZbLZcV7dB+GASciBvuLoCAIBCoU2bNjluj4iIUJs2bXLdH8gyQCB4vV75fD4lJyfL6XRKYnhAYcIk5MClEVQAAAAAIcjpdCopKSnY1QCAAhcZ7AoAAAAAAACcR1ABAAAAAABCBkEFAAAAAAAIGQQVAAAAAAAgZDCZJoB8ceGa78xUDgAAACCvCCoABJTdbpfNZpPL5bK22Ww2ud1uwgoAAAAAl0RQASCgHA6H3G63vF6vpHM9K1wul7xeL0EFAAAAgEsiqAAQcA6Hg1ACAAAAwBVhMk0AAAAAABAyCCoAAAAAAEDIIKgAAAAAAAAhg6ACAAAAAACEDCbTRL7xeDx+Kz8AAAAAyO7CtrLdbmdSchR5BBXIFx6PR06nUz6fz9pms9lkt9uDWCsAAAAgdNjtdtlsNrlcLmubzWaT2+0mrECRRlCBfOH1euXz+ZScnCyn0ymJdBgAAAC4kMPhkNvt9uuF7HK55PV6aTejSCOoQL5yOp1KSkoKdjUAAACAkORwOAglgCyYTBMAAAAAAIQMggoAAAAAABAyGPoBAAAABBmrpQHA/xBUAAAAAEHEamkA4I+gAgAAAAgiVksDAH8EFQAKxIXdWGl8AQCQHaulAcA5BBUA8pXdbpfNZpPL5bK22Ww2ud1uwgoAAAAA2RBUAMhXDodDbrfbb4Iwl8slr9dLUAEAAAAgG4IKAPnO4XAQSgAAAADIk8hgVwAAAAAAAOA8ggoAAAAAABAyCCoAAAAAAEDIIKgAAAAAAAAhg6ACAAAAAACEDIIKAAAAAAAQMggqAAAAAABAyIgOdgUAAAAAAP/jdrutr+12uxwORxBrAxQ8ggoAAAAACAF2u102m00ul8vaZrPZ5Ha7CStQpBBUAAAAAEAIcDgccrvd8nq9ks71rHC5XPJ6vQQVKFIIKgAAAAAgRDgcDkIJFHlMpgkAAAAAAEIGQQUAAAAAAAgZDP0AAAAACpjH4/GbhwAA8D8EFQAAAEAB8ng8cjqd8vl81jabzSa73R7EWgFA6CCoABAUrA8OACiqvF6vfD6fkpOT5XQ6JXEtBIALEVQgYOjCiLxgfXAAAM5xOp1KSkoKdjUAIOQQVCAg6MKIvGJ9cAAAAAAXQ1CBgKALIy4H64MDAAAAyA1BBQKKLowAAAAAgKsRGewKAAAAAAAAnEdQAQAAAAAAQgZBBQAAAAAACBkEFQAAAAAAIGQQVAAAAAAAgJBBUAEAAAAAAEIGQQUAAAAAAAgZBBUAAAAAACBkRAe7AgAAAACA3Lndbutru90uh8MRxNoA+Y+gAkBI4AIMAADgz263y2azyeVyWdtsNpvcbjdtJRRqBBUAgooLMAAAQM4cDofcbre8Xq+kczd2XC6XvF4v7SQUagQVAIKKCzAAAEDuHA4HbSIUOQQVAIKOCzAAoLDzeDx+oTwAIHcEFQAAAEA+8ng8cjqd8vl81jabzSa73R7EWgFA6CKoAAAAAPKR1+uVz+dTcnKynE6nJCaOBoCLIagAAAAACoDT6VRSUlKwqwEAIS8y2BUAAAAAAAA4jx4VuGJMCoX8dOHPFN1jAQAAgKKDoAJXhEmhkF/sdrtsNptcLpe1zWazye12E1YAAAAARQBBBa4Ik0IhvzgcDrndbr/eOi6XS16vl58vAAAAoAggqMBVYVIo5AeHw0EoAQAAABRRTKYJAAAAAABCBkEFAAAAAAAIGQz9AAAAAIAwknXFPeaKQ2FDUAEAAAAAYSCn1dEkVkhD4UNQAQAAAABhIOvqaBIrpKFwypeg4uDBg/r555915swZNWnSRFWqVMlWxhijRYsWafv27apdu7aaNWuWb2UAAACAguTxePyW2gYChdXRUBQEPKh48cUX9c0336hBgwaKiorS7Nmz9cILL+i1116zypw6dUr33Xef1q5dq+bNm2v+/Plq06aNvvnmG0VFRQW0DAAAAFCQPB6PnE6nfD6ftc1ms8lutwexVgAQPgIeVDRs2FCbNm1SbGysJGnGjBlq166d7rzzTt14442SpGHDhmnt2rVauXKlEhMTtWnTJjVp0kRffPGFevbsGdAyAAAAQEHyer3y+XxKTk6W0+mUxGSHAHA5Ar486WOPPWaFFJJ01113qXjx4lq7dq217ZtvvlHnzp2VmJgoSapdu7buueceffPNNwEvAwAACqf9+/frlVdeUcuWLdWwYUO5XK4cu9gvWbJE7dq1U82aNXX77bdr5syZ+VYGuJDT6VRSUpKSkpIIKQDgMgQ8qMjqp59+0unTp9W4cWNJ0pkzZ7RhwwY1aNDAr1zDhg2tMCNQZXKSnp6uo0eP+j0AAED46datm+Li4vTee+9p9OjRioyM1I033qitW7daZTZs2KDbb79dDRs21MSJE3X77berffv2WrRoUcDLAACAwMjXVT/27dunnj17qmPHjmrRooUk6fjx48rMzFSZMmX8ypYrV05HjhwJaJmcDB06VK+//vrVnRiAAnfhXVK6zwKQpClTpvjNSTVy5Ehdc801GjdunJ5//nlJ0jvvvKPatWtr2LBhkqRGjRppzpw5Gjp0qCZPnhzQMgAAIDDyrUfFwYMHdeedd6pKlSr68ssvre1xcXGSpBMnTviVP3bsmLUvUGVyMmjQIB05csR67Nix40pOD0ABuXC98KZNm6pp06ZyOp3yeDzBrhqAIMs6cXZmZqYyMjJUrFgxa9vcuXPVtm1bv3Lt2rXTvHnzAl4GAAAERr70qDh06JDatGmj+Ph4TZs2TSVKlLD2xcbGqmLFitq2bZvfMdu2bVPNmjUDWiYnMTExiomJufKTA1Cgsq4XzlrhAHLzt7/9TadPn9aDDz5obdu5c6cqVKjgVy4xMVFHjhzRiRMnVKJEiYCVySo9PV3p6enWc4abAgCQNwHvUXH48GG1adNGNptN06dPV3x8fLYy9957r8aNG6ezZ89KOtcrIiUlRffee2/AywAIfw6Hw5qM7Pzs6QBwobFjx+qNN97Qp59+qqpVq1rbMzMzFR3tf1/mfI+LjIyMgJbJaujQoSpdurT1qFKlypWeHgAARUrAg4q7775bGzdu1COPPKJx48bpiy++0BdffKGVK1daZV555RUdOHBA7du317vvvqu2bduqZMmSeuaZZwJeBgAAFG4TJkzQY489po8++khdu3b121e+fHmrR9Z5Xq9XsbGxKlmyZEDLZMVwUwAArkzAg4qGDRuqY8eOWrZsmebMmWM9LhyiUblyZa1YsUI333yz1q9frw4dOmjp0qV+E2MGqgwAACi8Jk6cqEceeUTDhw9Xr169su1v3ry5FixY4LctNTVVzZo1U0REREDLZBUTE6NSpUr5PQAAwKUFfI6KESNG5KlcYmKiXnrppQIpAwAACp+UlBQrpOjdu3eOZZ5++mndeeed+v777/X73/9es2bN0qRJk/TFF18EvAwAAAiMfFv1AwAAID/1799fGRkZeuONN1S5cmXrMWjQIKvM7bffro8++kh9+vRRfHy8HnroIQ0ePFiPPPJIwMsAAIDAyJdVPwAgv7ndbr/ndrudVUCAImbhwoXWhNoXyjqRd+/evdWzZ08dPHhQZcuWzTYpZiDLAECwXNg2ol2EcMcVFkBYsdvtstlscrlcftttNpvcbjcXZaAIybpc6MVERUWpfPnyBVIGAApSTm0j2kUIdwQVyDOPx2PNeJ71bjZQUBwOh9xut9/s+263Wy6XS16vlwsyAAAoUrK2jWgXoTAgqECeeDweOZ1O+Xw+a5vNZpPdbg9irVBUORwOLrwAAAD/j7YRChuCCuSJ1+uVz+dTcnKynE6nJMa+AQAAAAACj6ACl8XpdCopKSnY1QAAAAgpDJEFgMAhqAAAAACuAkNkASCwCCoAAACAq8AQWQAILIIKAAAAIAAYIgsAgREZ7AoAAAAAAACcR48KAIXGhZOX0eUWAAAUZbSLEM4IKgCEPbvdLpvNJpfLZW2z2Wxyu91clAEAQJFCuwiFAUEFgLDncDjkdrv9loVzuVzyer1ckAEAQJFCuwiFAUEFgELB4XBw8QUAABDtIoQ/JtMEAAAAAAAhg6ACAAAAAACEDIIKAAAAAAAQMpijAgAAALhMHo/Hb7JCAEDgEFQgRxdefCUuwAhPrB8OAMgPHo9HTqdTPp/P2maz2WS324NYKwAoPAgqkE1OF1+JCzDCB+uHAwDyk9frlc/nU3JyspxOpyQCcYQ+buAgnBBUIJucLr4Sf9AQPlg/HABQEJxOp5KSkoJdDeCiuIGDcERQgVxx8UU4Y/1wAAAAbuAgPBFUAAAAAEAhxg0chBuWJwUAAAAAACGDoAIAAAAAAIQMggoAAAAAABAyCCoAAAAAAEDIYDJNAEUG64cDAAAAoY+gAkChx/rhAAAAQPggqABQ6LF+OAAAABA+CCoAFAmsHw4AuFIej8cKuyX/oYRAuGJILEIZQQUAAACQC4/HI6fTKZ/P57fdZrPJbrcHqVbAlWNILMIBQQUk+d8p4C4BigruJAAALsXr9crn8yk5OVlOp9PaznUD4YohsQgHBBXI8U4BdwlQmHEnAQBwuZxOp5KSkoJdDSAgGBKLUEdQgRzvFHCXAIUZdxIAAACA0EVQAQt3ClCUcCcBAAAACE2Rwa4AAAAAAADAeQQVAAAAAAAgZDD0AwAAALgAq6EBQHARVAAAAAD/j9XQACD4CCoA4P9deNeMlW8AoGhiNTQUVbSDEEoIKgAUeXa7XTabTS6Xy9pms9nkdru5SANAEcVqaCgqaAchFBFUACjyHA6H3G6333hkl8slr9fLBRoAABRqtIMQiggqAEDnLtJZL8Z0gQQAAEVBTu0gIJgIKgAgC7pAAgAAAMFDUAEAWdAFEgAAAAgeggoAyAFdIAEAQFHGEFgEE0FFEeXxePzuFgMAAAAAQ2ARCggqiiCPxyOn0ymfz2dts9lsstvtQawVAABAcHADB/gfhsAiFBBUFEFer1c+n0/JyclyOp2S6M4F5AVdIAGg8OEGDpAdQ2ARbAQVRZjT6VRSUlKwqwGEPLpAAkDhxQ0cAAg9BBUAcAl0gQSAwo8bOAAQOggqACAP6AIJAAAAFIzIYFcAAAAAAADgPHpUAAAAoEhhlQ/g8jGpOAoSQUURwQUZAACAVT6Ay8Wk4ggGgooigAsyAADAOazyAVweJhVHMBBUFAFckIH8QRdIAAhfrPIB5B2TiqOgEVQUIVyQgcDIrQvk+PHjVb58easMF3QAAADg8hFUAMBlytoFcv/+/erQoYPatWtnlWHsJgCEDubqAgKPnqXITwQVAHAFsnaBzGnsZmpqqjXcSuIiDgDBwFxdQGAxuSYKAkEFAATAhcFFThdwiYs4AAQDc3UBgcXkmigIBBUAEGBZL+BSzr0saCgDQMFhri4gcJhcE/mNoAIA8kHWCzjdJAGg4DAnBVDwmLMCgURQAQAFgG6SAFAwmJMCKFjcjEF+IKgopLiTAIQeukkCQP5jTgqgYHEzBvmBoKIQ4k4CAAAo6piTAig4Od2MyXqzlMAQl4OgohDiTgIQPhjPCQAAChNWP0MgEFQUYtxJAEJXbuM5x48fr/Lly+d6DBd3AMiOIa9A6LjY6mcMB0FeEVQAQBBkvYjv379fHTp0ULt27XI9hjsRAJAdQ16B0MO8XLhaBBUAECRZL+JZ7z5ciDsRAJAzhrwC4YMhr8grgopCgi6PQPjj7gMA5E1O7R6GvAKhKy9DXgkucCGCikKALo9A0cGdCABFHe0eIPzkZcgrQ1xxIYKKQoAuj0Dhx50IADiHdg8Qni425JUhrsiKoKIQocsjUHjl9U4EwQWAwuDCoR1ZMdQDKBxyGvJKz1GcR1ARppiTAih6LnYngi6UAAqLnIZ2ZMVQD6BwoecosioUQcXevXu1Y8cO1ahRQwkJCcGuTr5jbCYAKW9dKFNTU+kaDQTQgQMHtGXLFlWpUkUVKlQIdnUKhay9J9xud7ahHVnx9wwoXJjDAlmFdVCRmZmpp556Sl9++aVq1qypzZs367nnntNbb70V7KoFXNYeFIzNBJDVhcFFXu5M5AV/W4D/eemll/TOO+9YbY7HH39cH3/8sSIjI4NdtbCVW+8Jm82m1q1b8/cHKEIu9wZMVrRZCpewDio+/vhjjRkzRitXrlS9evW0aNEi3XLLLWratKk6dOgQ7OpdlQuDifOJYtYeFFzAAeQmL3cm8uJS4QaNAhQV48aN07BhwzR37ly1bNlSbrdbLVq0UOPGjfX0008Hu3oh7VLzTeTUe4K/LQAudQMmK4aKFC4RxhgT7EpcqaZNmyopKUkjRoywtt19992KiorSlClT8vQaR48eVenSpXXkyBGVKlUqv6p6SXkJJvjFC3MnTkjx8ee+Pn5cKlEiuPVBkXOxfxZyktPfoqwut5cGf7tCU6hcC0NZ+/btFRkZqZSUFGvbE088oTVr1mjZsmV5eo3C+jlf7G9LXv+O0J0bOaLthAtc7t+aK7nZkvU9aLcEXl6vhWHboyIjI0Nr1qxRr169/LY3b95cn376aa7HpaenKz093Xp+5MgRSec+sEDZu3ev9u7dm+fyXq9XLpdLJ0+etLbFxcVp3Lhx1rwTCQkJqlKlit9xgawzCsCJE//7+uhRKSMjeHVBkVSmTBmVKVMmz+Vr1aqln3/+WQcOHMhx//m/XZfTSyMuLk7JycnMqXOVKlSoEND5Ec5fT8L43kW+W7FihXr37u23rXnz5vr666+VmZmZ4/CPUGxzBFpObZissrZpskpISFCZMmVo1yA72k64wMXaMVnbLHlpo2Rtk+T2P1lRb7cErc1hwtShQ4eMJDN27Fi/7f/6179MbGxsrscNHjzYSOLBgwcPHjx4ZHns2LEjvy/fYSsmJsa8//77ftvGjBljJJnDhw/neAxtDh48ePDgwSPnx6XaHGHbo6JYsWKSpFOnTvltP3nypIoXL57rcYMGDdKzzz5rPc/MzNTBgweVkJCgiIiI/KlsmDl69KiqVKmiHTt2FKquqaGAzzb/8NnmLz7f/BMKn60xRseOHVPFihWD8v7hoFixYjm2OSTl2u4oiDZHKPz8FLSieM4S5815F35F8ZylonfeeW1zhG1QUaJECSUkJGjXrl1+23ft2nXRcUQxMTGKiYnx23Y5XaGLklKlShWJX5Zg4LPNP3y2+YvPN/8E+7MtXbp00N47HFStWjXHNofdbldcXFyOxxRkmyPYPz/BUBTPWeK8i5qieN5F8ZylonXeeWlzhPV6Wm3btvWb1CojI0NTp05V27Ztg1grAABQ2LRt21ZTp05VxgVj5CdNmkSbAwCAfBDWQcUrr7yiVatWqU+fPkpJSVGXLl107NgxPffcc8GuGgAAKESef/55HT58WF27dlVKSop69eqlX375RS+//HKwqwYAQKET1kFF/fr1tWjRIp0+fVrvvfeeypUrp8WLF6tSpUrBrlpYi4mJ0eDBg7N1V8XV47PNP3y2+YvPN//w2YaHSpUqacmSJSpbtqzee+89nT17VgsXLlT9+vWDWq+i+PNTFM9Z4rw578KvKJ6zVHTP+1IijGEtMgAAAAAAEBrCukcFAAAAAAAoXAgqAAAAAABAyCCoAAAAAAAAISM62BVAwcvIyNDatWsVFRWl+vXrKzLy4nnVwoULlZmZ6betevXqTFqaC6/Xq/Xr18vpdCohISFPx+zZs0c7d+5UzZo1Va5cuXyuYfgyxmj58uWKjo5W48aNL1o2MzNTCxcuzLa9Tp06uuaaa/KrimErIyNDGzduVHR0tKpXr67o6LxdHnbs2KHffvtNderUKTJrf18Jr9erHTt2qFq1aipbtuxFy546dUrLli3Ltr1hw4YqU6ZMPtUQoWzTpk06duyY6tevr9jY2FzLnT59Wj///HOO+8Ltun369GmtW7dOcXFxqlevXp6OOXr0qLZu3arixYurZs2aKl68eD7XMvAOHjyozZs3q3Llyrr22mvzdMyePXu0a9cu1ahRI2zbEB6PRx6PR0lJSbLZbHk6ZuvWrTp48KDq1aunEiVK5HMN88eqVat05swZNWvWLF+PCSWnTp3S8uXLde2116p69ep5OmbXrl3yer2qWbOm4uPj87mG+WPPnj3avHmzGjVqpNKlS+fpmK1bt+r48eOqXr162J73VTEoUpYuXWocDoepXLmySUxMNLVq1TLr1q276DExMTGmfv36plWrVtbj22+/LaAah4+1a9eaRx991FSoUMFIMt99990ljzl79qzp3r27iY2NNfXr1zcxMTHm9ddfL4DahpfMzEwzbNgwU7t2bVOmTBnTtGnTSx5z7NgxI8lcd911fj+706ZNK4Aah5c333zTJCYmGqfTaapVq2YqVapkUlJSLnrMyZMnTYcOHUxcXJypV6+eiYuLM++//34B1Th8/Pzzz+bWW281drvdNGnSxMTFxZnHH3/cpKen53rMpk2bjCTTrFkzv5/dxYsXF2DNEQp+++0307JlS1OmTBlTq1YtU7ZsWTNp0qRcy+/bt8/vZ6ZVq1amYcOGRpL58ssvC7DmV+fHH3805cuXN9WrVzcJCQmmcePGxuPxXPSYv/zlL8Zms5nGjRubGjVqmGuuucZ8//33BVTjwHj99detNldMTIzp3r27OXv2bK7l9+7da+644w5TpkwZ06RJExMfH28GDhxYgDW+evPmzTPt27c3CQkJRpJZs2bNJY85cuSIadOmjYmPjzd16tQx8fHxZtSoUQVQ28D5+OOPTYMGDUzZsmVN1apV8+2YUOL1es3zzz9vKlasaOLi4kz//v0vecz06dNN48aNTcWKFc11111nbDZb2P2ML1261HTs2NGUL1/eSDKzZ8++5DGTJk0y9erVM3Xq1DH169c3NpvNvPLKK/lf2RBDUFGEpKenm6pVq5qePXsaY4zJyMgwDz30kGnQoIHJzMzM9biYmJhL/tMCY8aOHWtGjRplDh8+nOegYtiwYaZcuXLm119/NcYYM2fOHBMVFWWmTp2a39UNK+np6ebZZ581GzZsMP3797+soGLRokUFUMPwdfbsWfPSSy+ZAwcOWNsGDx5sbDab2bNnT67HDRw40FSuXNns3r3bGGPMhAkTjCT+mc7i66+/NnPmzLGeb9682ZQvX/6iDY7zQcXWrVsLoIYIZQ8++KBp1qyZOXHihDHGmKFDh17ydzOrV1991ZQsWdIcO3Ysv6oZUIcOHTJly5Y1f/nLX4wx5/7+t27d2tx22225HjNv3jwjycycOdPaNmjQIBMXF3fRUDCUTJkyxURHR5t58+YZY4zZuHGjKVOmjHnnnXdyPeahhx4ySUlJ5siRI8YYY7Zu3Wrsdrv56quvCqTOgfDhhx+alJQUk5aWlueg4g9/+IOpW7euOXjwoDHGmBEjRpjo6GizYcOG/K5uwDz33HNm9erVZsiQIXkOHa7kmFCyYsUK8/e//93s37/fNG3aNE9BxYcffmhWrVplPV+0aJGJiYkJq+B15MiRZuzYsWbLli15Dir+/e9/m+3bt1vPZ82ale1vXFFAUFGETJs2LVvjd9myZZf8Zy4mJsaMGDHCLF261Ozfv78AahreTp48meegon79+qZfv35+22699VbTsWPH/Kpe2LvcoOK7774zy5YtM4cOHcr/yhUSe/fuNZIuGpglJiaa1157zW9bw4YNTZ8+ffK7emHP5XJd9J+u80HFTz/9ZJYvX26OHj1agLVDqNi/f7+JjIw0o0ePtradPHnSlCxZ0rz77rt5eo2MjAzjcDhM796986mWgff555+b4sWLW/98G3Pun/iLhXfjxo0zkvzCmJSUFCPJL4QNZR06dDBt2rTx2/bkk0+aBg0a5HpMqVKlzD//+U+/bd26dTOtW7fOlzrmpxUrVuQpqDh16pSx2Wzmgw8+sLZlZmaaihUrmpdeeim/qxlwVxI6hGtQcaG8BhU5admypenVq1dgK1QAduzYkeegIqvTp0+b6Oho88UXXwS+YiGMyTSLkBUrVighIUHVqlWztiUlJSkqKkorVqy46LEvvPCCevbsqSpVqqh9+/bau3dvPte28Dt16pTcbreaNm3qt7158+aX/H4g75566il1795diYmJ6ty5sw4fPhzsKoW8pUuXSpJq1qyZ4/7du3frt99+42f3CmRmZmr58uWqVavWJcu6XC65XC4lJCSod+/eOnnyZAHUEKFi9erVyszM9Ps9i42NVaNGjfL8ezZz5kx5PB716tUrv6oZcCtWrFDt2rX95rxp3ry5tS8n7du3180336zHH39cP/zwg77//nsNGjRIgwYNCps5G1asWJHj31S326309PQcjylXrpx27drlt23Xrl1KS0vLt3oG24YNG+Tz+fw+q4iICDVr1ozrTxFw7NgxbdiwIU/X0HB3+PBhzZ8/X9OmTVOXLl3UuHFjdezYMdjVKlBMplmEHDx4MNvkjhERESpXrpwOHjyY63EffPCBevToocjISO3atUt33323Hn/8cc2YMSO/q1yoHT58WMaYbN+ThISEi34/kDfR0dEaNWqUHnvsMUnSli1b1KZNG/Xt21dff/11kGsXurxer/74xz+qU6dOqlu3bo5lzv988rN7+d544w1t27ZN48ePz7VMiRIlNHnyZN13332SpDVr1ui2225TyZIl9c477xRUVRFkgfg9++yzz9SkSZOwmnQvp7bK+bAht/OOiYnRn/70J/Xr108bN27U8ePHlZCQoEceeSTf6xsoOZ13QkKCMjMzdfjwYSUmJmY75plnntHAgQNVoUIFXXfddZoxY4bS0tLk8/mUnp6umJiYgqp+gbnY74Xb7Q5GlVCAnnrqKcXGxqpnz57Brkq+27BhgwYOHCiv16t9+/Zp2LBhRW5CTXpUFCHFihXTqVOnsm0/efLkRWfG/sMf/mCtDFKpUiW9+uqr+u9//8s/JFepWLFikpTte3Kp7wfyJjY21gopJKlGjRp64YUX9P333+vs2bNBrFnoOnLkiNq1a6cKFSro008/zbUcP7tX5pNPPtHQoUP17bff5hoCSdK1115rhRSS1KhRI/Xr10+jR48uiGoiRFzt79mBAwc0adKksOpNIeXcVjn/PLfznjFjhh555BF9++23WrNmjbZu3aoHH3xQt9xyi7xeb77XORByOu/zvahyO+/+/fvr22+/1fLly/X3v/9dNptNr776qqKiogrt32KuP0XXgAEDNGXKFE2ePDnPq+qFsxYtWmj+/Plav369xo8fr6eeekpjxowJdrUKFEFFEVK1alXt27dPZ86csbYdPXpUx48fl8PhyPPrnE/1d+/eHfA6FiXlypVTyZIlc+y2eTnfD+RdYmKiTp8+HTYN14J09OhR3XnnnYqKitIPP/ygkiVL5lq2SpUqVg+rC/Gzm7sRI0aof//+Gj16tO6///7LPj4xMVF79uzJtlQ0Cq+qVatK0hX/nn311VeKiorSo48+mi/1yy9Vq1bN8Zwl5XreU6ZMUf369XXrrbda2/r27avDhw9r3rx5+VbXQMrtvEuWLHnRJY0ffPBBffPNN/rxxx81ZMgQrVmzRvXr11dERER+Vzkorvb3AuFp0KBB+s9//qMZM2aEVQ+xQLn11lvVvHlzTZ8+PdhVKVAEFUVImzZtlJ6erpkzZ1rbJk2apOjoaN12223WtgULFlgXgBMnTmR7nf/+97+Ki4tTjRo18r/Shcz27dutNe4jIiJ0xx13aPLkydb+M2fOaNq0aWrbtm2wqhi2MjIyNH/+fO3bt09S7j+75cuX1zXXXFPQ1Qtp50MK6dxnlNP63r/++qs1/tdms+nGG2/0+9k9ceKEfvzxR352c/Dpp5+qX79++vbbb/XQQw9l23/y5EnNnz/fmj8lt5/d+vXrW73bUPg1atRIiYmJfr9nGzdulNvt9vs9W7NmjdavX5/t+M8++0ydOnXK8fc5lLVt21a7du3S8uXLrW2TJk1SqVKl1KJFC0nS2bNnNX/+fO3fv1+SVL58ef32229+N2J27Nhh7QsHbdu21fTp0/16/E2aNElt2rSxnu/bt0/z589XRkaGJGWbt2bfvn0aN26cunXrVjCVLiBut1tr1qyRJFWuXFn16tXz+73Yt2+fFi1aVOiuPxe2x4uKQ4cOaf78+X49Zv7yl7/oo48+0owZM6y/AYXNhW2szMzMbD2Gzpw5I4/HUyR6kvgJ9myeKFi9evUyFStWNF999ZX5/PPPTbly5cwLL7zgVyYqKsr84x//MMYYk5ycbO69917z5ZdfmunTp5sBAwaYYsWKmWHDhgWj+iHt4MGDJjU11VpCaMiQISY1NdVs3rzZKvPiiy+axMRE6/mKFStMXFyc6du3r5k8ebJ58MEHTYUKFcxvv/0WjFMIacuXLzepqanm4YcfNnXr1jWpqakmNTXVZGRkGGPOLWknyYwcOdIYY8y//vUv07lzZ/P111+bqVOnmr59+5ro6GhrP845ffq0ufHGG80111xjpkyZYn2uqampZu/evVa5nj17+s0+P2fOHFOsWDEzcOBAM2nSJNOmTRtTs2bNsFkCsaB8++23JiIiwjz77LN+n21aWppVxu12G0lm+vTpxhhjXnrpJdOjRw8zZswYM3nyZNOtWzdTrFgxM2XKlGCdBoLks88+M8WLFzf//Oc/zffff28aNWpkbrnlFr8lxVu1apVtpajFixcbSWb+/PkFXeWAuO+++0ydOnXMmDFjzIcffmji4uL8VrfYv3+/kWQtw7l161ZTunRp88ADD5ipU6ea0aNHmwYNGphmzZqZ06dPB+s0LsvevXtNhQoVTIcOHczkyZPNU089ZeLi4szKlSutMiNHjjSSrFWsJk+ebLp06WImT55svvnmG9OgQQNz2223hc05G2PMzp07TWpqqnVuo0aNMqmpqdbS18YY88ADD5hbbrnFej5hwgQTFRVlhgwZYiZMmGBuuOEG07hx47A679WrV5vU1FTTq1cvU6FCBevacOrUKavMhe3xvB4Tys6ePWvVuW7duubhhx82qampZsWKFVaZ86v1bNq0yRhjzFtvvWUiIiLMsGHD/K6h69atC9JZXL7ffvvNpKammvHjxxtJ5v333zepqal+y49e2Mby+XymQYMG5t133zUzZswwY8aMMXfccYex2+1my5YtwTqNoIgwxpigpSQocGfPntXHH3+sadOmKTIyUg888IDfHBTSue5FTz75pDUJ1U8//aTk5GTt2rVL1atX1+OPP64bbrghWKcQshYuXKgXXngh2/aOHTvqz3/+s6RzY9SnT5+uSZMmWftXrFih9957Tzt37lTdunX14osvWl0b8T/dunXTli1bsm2fOXOm4uLidPz4cbVr104vvfSS7r77bklSSkqKvvvuO+3fv181atRQ79691bhx44Kuekg7cuSI2rdvn+O+gQMH6t5775Ukvf3221q3bp1GjRpl7V+wYIE+/PBD/fbbb2rUqJE1qRv+569//aumTZuWbXvVqlWtSV09Ho+6du2qd955Ry1atJAxRmPGjNHkyZN15MgR1alTR08//bRq165d0NVHCJg4caJGjRqlY8eOqVWrVhowYIBKlChh7X/66aeVkJCgIUOGWNuGDx+uOXPmaNy4ccGo8lU7deqU3n33Xc2ePVtxcXF65JFH1KVLF2v/+b9br776qtUbbNu2bXr//ffldrsVGxur5s2bq1+/fhcdxhZqtm/frrffflsbNmxQ5cqV9cwzz+j666+39k+fPl1vvfWWfvjhB2tSvXHjxik5OVkZGRlq166devXqZc3jEA5Gjx6tDz74INv2Z599Vh06dJB07o76sWPHNHz4cGv/zJkzNWLECB08eFDNmjXTiy++eNEhMqHmj3/8Y46rlHz//ffWdTRrezwvx4Sy8+20rGrXrq2RI0dK+l9besyYMapUqZKefvpprV69OtsxrVq10ttvv53vdQ6E87+3WfXo0UM9evSQlL2NtXv3bn3wwQdasWKFSpQooeuvv15PPvlkketRQVABAAAAAABCBoNdAQAAAABAyCCoAAAAAAAAIYOgAgAAAAAAhAyCCgAAAAAAEDIIKgAAAAAAQMggqAAAAAAAACGDoAIAAAAAAIQMggoABWL06NHat29fsKsRdvbt26fRo0cHuxoAAISNSZMmyePxBLsaYefYsWMaPXq0Tp06FeyqAAQVQH5bsGCBRo8erfXr12fbt379eo0ePVoLFiwIQs0KVpcuXfTLL79csty0adO0e/fuAqhR/tm1a5fGjRsXkNf65Zdf1KVLl4C8FgCgcFu+fLlGjx6tlStXZtu3bds2jR49WrNmzSr4ihWwPn36aOHChZcs99NPP2nz5s0FUKP8c+DAAY0ePVpnz5696tfatWuXunTposOHD199xYCrRFAB5LN3331XXbp00bPPPpttX//+/dWlSxe9++67QahZ6Dl8+LAefPBBZWZmBrsqV2Xp0qV6/PHHg10NAEAR8/nnn6tLly7q0aNHtn2DBw9Wly5dNHjw4CDULPScPXtWHTt2DPt/yjdt2qQuXbrQCwKFTnSwKwAUBTfeeKNmzZqlnTt3qnLlypIkj8ejuXPn6oYbbshWPjMzU0uXLtXevXtVs2ZNNWzY0G//4sWLtW3bNkVERMhut6tx48ay2+3ZXsftdmvTpk1yOBy67rrrFBl5LpucOnWq6tWrp5o1a1plU1NTFR8fr+uvv16S9OOPP6pixYq65pprtGTJEpUqVUqtW7eWJB0/flwLFy5URkaGmjRpomuvvTbbe2/evFlr165VjRo1stU/N9OnT1eDBg2szygnF3vv6dOnq0yZMn6f6erVq7Vp0yZ16NBBO3bsUFpamu6//36tXr1aHo9HSUlJqlKlSrb32b17t5YtW6bSpUsrKSlJJUuWzFZm7969WrZsmcqWLavmzZurWLFi8nq9mj9/vs6ePWsN2ahXr56aNGmS59ddunSp9u3bl+fPDQCA85o0aaJff/1VK1assK7pR48e1ffff6+bb75ZGRkZ2Y5ZsWKFtm/frqpVq6pJkyaKiIiw9q1cudLqFVq2bFk1btxYFSpUyPYamzdv1i+//KIKFSro+uuvV3T0uX8zZs2apcTERL9r2tKlS3X27Fnrej1//nzZbDbVqFFDixcvVlRUlNq2bStJOnXqlBYuXKiTJ0+qUaNGcjgc2d57586dSktLU5UqVdS4ceM8fU6pqamKi4tTUlJSrmUu9t5z5syRMUa33XabtW3jxo1asWKFOnTooEOHDumnn35Sp06dtG7dOm3ZskUNGjRQrVq1sr2P1+vVkiVLFBsbq6SkJJUtWzZbmYMHD1plbrjhBsXGxur48eP68ccfJUnjxo1TTEyMqlWrppYtW+b5dVevXq3t27erXr16efrcgIJCUAEUgGuvvVZ33XWXvvjiC7388suSzt31aN++fbYGw44dO3T//ffr5MmTqlOnjlasWKHGjRtbFyBJSktLU2pqqqRz//iuXLlSn376qTp16mS9zuOPP66pU6eqVatW2r17t2w2myZPnqzSpUurf//+ev755/2Cirffflu1atWyGjUvv/yyYmJitH37djVq1EitW7dW69atlZKSou7du6t+/fqKj4/XokWL9Oqrr/r1GHnvvfc0cOBAtWrVSl6v96LBw4WmTJmie++9N9f9l3rvQ4cOqWPHjlq2bJnq16+vvXv3qk2bNnrmmWcUERGhhQsXqlevXrr++ut16tQpxcXFacmSJfrPf/6jbt26We/z2muv6V//+pdatmyp48ePa+PGjfr22291++23W2UGDx6sf/zjH2rWrJkk6eTJk5oyZYoOHz6sJUuWKCMjQxMnTpQk3XfffWrSpMklXzczM1OdO3fWzJkzdcMNN+iXX35R/fr18/TZAQAgSSVKlNAjjzyizz77TB988IEk6ZtvvlGTJk3kdDq1du1aq+zBgwf10EMPyePx6LrrrtPatWtVqVIlTZo0yfqndt26dUpJSZEk7d+/X0uWLNE//vEPPfXUU9brDBgwQCNGjFDr1q114MABnTlzRhMnTlSlSpU0ePBgtWnTxi+o+Pe//63jx49bQcWwYcO0Z88eeb1e1a1bV9dff73atm2r+fPn6+GHH5bD4VD58uW1cOFCPfnkk/rrX/9qvdbXX3+tnj17qnnz5kpPT1exYsV08uTJS35OU6ZMUfv27f1CmQtd6r3Pnj2ru+++W7Nnz9ZNN92kY8eO6Z577tE999yjzp07W0M3R40apZ07d6p8+fJasGCB3nrrLT333HPW+3zwwQd6+eWX9bvf/U6ZmZlasWKFRowYoY4dO1plhg8froEDB+q6665TXFyc9u/frwkTJqhkyZKaO3eupHNtpOjoaLVu3VotW7bM0+v27dtXX375pVq1aqVff/1VdevWveTnBhQYAyBfdezY0XTs2NFMnDjR1KhRw2RmZpqMjAzjcDjM1KlTzQMPPGA6duxolb/55ptNv379TGZmpjHGGJ/PZ66//nozZMiQXN9jzJgxply5csbn8xljjFm/fr2RZLZv326VWbhwodmzZ48xxpiaNWuajz/+2O812rdvb/r37289b9GihbHb7Wb37t3Wth07dpj4+HgzY8YMa9uqVatMbGysWbVqlTHGmG3btpnixYub7777zirTvXt3I8nMnj0713M4e/asKVeunFm8eHGO+/Py3sYY8+ijj5rrrrvOnDx50rRr187cfPPNJiMjwxhjzLfffmskmRdffNEqP3z4cFOqVCnj9XqNMcaMHz/eVKhQwezYscMq88knn5iKFSua9PR0Y8y5zzs6OtqkpqZaZVauXGm2bt1qjDFmwoQJpkSJEn71z8vrfv311yY+Pt5s2bLFGGPM0aNHTcOGDQ1/qgEAedG3b1/TqlUrs3jxYlO2bFlz8uRJY4wxTZs2NZ999pnp06ePadWqlVW+U6dOpnPnzubs2bPGGGNOnz5t2rZta/r27Zvre8ydO9fExsaavXv3GmOMOXz4sImIiDA///yzVWb16tVm06ZNxhhjWrVqZQYPHuz3Gj179jSdO3e2nj/wwAMmLi7ObNy40dp29OhRY7fbzVdffWVt27ZtmyldurSZOXOmMcaYgwcPmtKlS5vhw4dbZV588UUjyXz77bcX/azq1KljJk6cmOO+vLy3McY899xzpmrVqubw4cOmW7dupkGDBtZnPnv2bCPJdO/e3WrTfffdd6ZYsWJm8+bNxhhjFixYYEqVKmV++eUX6zUnTpxoSpUqZQ4cOGCMMWbOnDkmIiLCTJgwwSrz66+/mrVr1xpjjFm0aJGRZI4dO2btz8vr/vTTTyY6OtqsWLHCGGNMenq6ufnmm40kq70IBBNzVAAFpH379vL5fJozZ47++9//KiMjQ3fddZdfmS1btmjevHmqWbOmxo0bp++++04pKSmqWbOmZs+e7Vd23759mjVrlsaMGaP09HQdPHhQv/76qyQpNjZWERERWr16tVX+hhtuyLGr5sU8/PDDfkMrxo4dq5IlS+ro0aP67rvv9N1332n9+vVKTEy0Ev0JEyaoYsWK+v3vf28d98ILL1zyvRYuXKjo6Gj97ne/y3F/Xt5bkj766CMdO3ZMzZs31+LFi/XVV19ZQ17OGzBggPX1k08+KUn64YcfJEkjR45Uw4YNtXjxYn333XcaO3asihcvrt27d2vDhg2SpC+//FIPPvigbrrpJut1GjdurGrVquV6fnl53bFjx+r3v/+9qlevLkkqWbKknn766Ut+dgAAXKhFixaqWLGiJkyYoFWrVmnDhg1+vS6lcys8jBs3TvXq1dOECRP03XffacKECapatWq2NsehQ4c0Z84cjR07Vrt371ZkZKTVxoiKilKxYsW0evVqGWMkSY0aNcpxiMPFtGvXTrVr17aep6Sk6MSJE4qJibGu+0uWLFG1atWs+s2YMUOZmZnWtVySXnzxxUu+16+//qrt27erTZs2Oe7Py3tL0l//+leVK1dOrVu31pgxY/TNN98oNjbW77Wef/55q9fG73//e1WpUkUTJkyQJH3xxReqU6eO1q1bZ7UNTp06pVOnTiktLU2SNGrUKLVq1UoPPvig9Zo1a9ZUgwYNcj2/vLzumDFj1KZNG2toavHixfXMM89c8rMDCgpDP4ACEh0drW7duunzzz/XyZMn1b17d0VFRfmV2bZtm6Rz4yaLFStmbY+KivIbQ/nOO+9o8ODB1jjR82XPL/9ZtWpVDR8+XD179lRMTIxuu+02de/e3W8cZV5knXti27ZtysjI0Pfff++3vWXLlkpMTJR0bu6NrP+wn//H+2KmTJmie+65J1uocDnvLUmlSpXSgAED9PTTT+uNN97INpY1Pj5eCQkJ1vPo6GhVqVJF27dvt97HGJPtfTp37mw1NDwej+6///5LnlPW+ufldc8PvTkvL58dAABZ9ezZU59//rnq1q2rzp07Kz4+3m//jh07lJGRoaVLl2ZbmaxVq1bW11988YX+9Kc/yel0qlKlSipevLgyMzOtNkd8fLxGjhypgQMH6uWXX9att96qRx999KJDOXOSU5sjOjo62ypa9erVU9WqVSWdu25WqlTJmg9DOjePRunSpS/6XikpKbrttttUokSJHPfn5b2lc//cv/LKK+rQoYOefvppXXfdddleK6c20YVtjkOHDmVrGzz00EOKi4uzzrFOnToXPZ+c6p+X172S9hpQUAgqgALUo0cPXX/99crMzNTbb7+dbX+pUqUknZv/IKeLnXRuZYwBAwZoxowZ1kRTXq9XY8aMse5kSOfGHT711FNas2aNJk6cqDvvvFMTJ05U+/btFRkZmW1ljZxmi846brNUqVIqUaKENUlkThISEnTo0CG/bVmf52TKlCl64403ct2fl/eWzo23feutt1S/fn199NFHeuqpp/wmGj1x4oROnz6t4sWL+9XvfJlSpUqpTp06+vzzz3N9jzJlyujAgQOXPKes9b/U617pZwcAQFaPPfaYBg0apCVLlmjGjBnZ9p9vc/Tv31933nlnjq+RmZmpP/7xj/roo4/02GOPWdvi4uL82hxdu3ZV165d5Xa7NWXKFHXu3FnDhw9Xjx49rqrNkZmZqW+++SbXmxg5XTfPnj2r48eP51j+vClTpqhDhw657s/Le0vn5qd6+eWXVb9+fX311Vd67rnnVKNGDb8yhw4d8gtEcmpzXKxtczVtjvxorwEFhaEfQAGqW7eu+vbtq2eeecZvIsvzzveQ+OSTT7Lt2717t6RzvSaMMX4THmVNzA8cOKBTp04pMjJSjRs31uDBg5WUlKQlS5ZIkipVqmQNE5GkI0eOaMWKFZesf7t27bR169ZsDZ4TJ07oyJEjkqSbbrpJa9as8VuXfPz48Rd93W3btunXX3/NtaGU1/eWpN69e6tSpUpKS0tT1apV1bNnT7/yxhhNnjzZer5kyRLt2bNHN954o/U+48ePt+4Unbdr1y7r6zvvvFMTJkzQsWPHrG2nT5/W0aNHJZ27u5Senu7XiMvL6950002aMmWK31rol/rsAADIid1u1yuvvKKuXbvmuMJY5cqV1aBBgxzbHOevTcePH9fx48f92hyTJ0/W6dOnrefHjh2zrodOp1MDBgxQmzZttHjxYknZ2xynT5/WwoULL1n/O++8UydPnlRycrLf9tOnT8vr9Uo6t6ravn37rPeSpEmTJuW4ssl5R48eVWpq6kV7fOTlvaVzwzrOnDmjJUuW6Oabb9Zjjz2W7b3PT6wtnVsZZdWqVVaPlXbt2mnWrFnatGmT3zH79u3TmTNnrLr8+OOP2rNnj7U/MzNTBw8elCSrp8yF4U9eXvemm27Sjz/+qBMnTlj7aXMglNCjAihg//jHP3LdV6xYMX322Wfq2LGj9u/frzvvvFNer1cpKSnq2rWr+vXrZy1X+uijj+qJJ57Q+vXr9eWXX/q9ztatW/Xoo4/q4YcfVu3atbV69WqtWbPGmv27W7du6tevnxISElSuXDmNHDkyT3W/6aab9Kc//UkdOnRQv379VLduXW3cuFHjx4+3VhS55ZZbdPfdd+vOO+/UM888I6/Xe9FeBNK5Lpg333xzjkt1Xs57f/755/rhhx+0cuVKxcbG6uuvv1aTJk30n//8R71795b0vzGYGzduVFxcnP7+97+rW7du1mzkzz77rKZOnarf/e53evrpp1WmTBmlpaVp0aJFWrNmjSTpz3/+syZOnKgWLVqoV69eks7NqD5q1CiVKlVKjRo1UlRUlAYOHKjrr79e9erVy9Pr9u/fXyNGjNBdd92lzp07a+HChZozZ06evjcAAGT10ksvXXT/p59+qnbt2unuu+/WAw88oOPHj+uHH35Qy5Yt9eabb6pUqVK67bbb1KdPH/Xt21c7d+7Uhx9+aK1CJp3r1dm2bVt16NBB9evX1+bNmzVz5kzrn97HHntMDz74oGrUqCGHw6FvvvnGL+jPTZ06dTRkyBD16tVLy5YtU+PGjbVt2zZ9//33+vTTT2W32+V0OtWjRw899NBDGjBggNLT07PVL6sZM2aobt26fkM4ruS9p02bphEjRmjBggWKj4/X559/ruuuu05vvfWWXn31Veu1/va3v2nfvn265ppr9N5776lNmza64447JEndu3fXhAkT1KpVK/3xj3/Utddeq9WrV2v69OlatWqVihUrpu7du2vMmDFq2bKl+vbtq7i4OI0ZM0ZDhw5V69atVaNGDZUtW1Yvv/yybrnlFlWvXj1Pr/vEE0/oX//6l26//XY98cQTWrduncaMGXPJ7wtQUOhRAeSzm266yW/Sxaxat27tt/+ee+7R2rVr1aBBA82fP1/Hjx/Xe++9p379+kk6N1/F7Nmzdccdd2ju3LkqVqyYFi9erM6dO1uTZTZr1kyzZs2SzWbT3LlzFRsbq7S0NGuiyieeeELJycnasWOHtm7dquHDh+vFF1/0mwejbdu2OU7U9K9//UuTJk1Senq65s+fr4SEBKWmpvqtv/3999+rX79+SktLU2RkpBYtWqTOnTvrmmuuyfEzuNSypHl57/T0dC1atEgjR460JvCqWbOmRo4cqSVLllh3GkqWLKk5c+bo1KlTWrVqlV599VV9+umn1nvYbDbNmzdPQ4YM0ebNm7Vy5Uq1aNFCy5Yts8rEx8drwYIF6t+/v9asWaO9e/dq5MiRcjqdkqTExETNnDlTJ0+e1KRJk7Ru3bo8vW6ZMmX0888/q3nz5vr555+VlJSkGTNmqHPnzpf8bAAAaNq0aa4TRErS7373O7/9LVu21C+//KLWrVtr0aJF2r9/v15++WW9+eabVplJkyapU6dOSk1NtSYF7969uzW/QfXq1bV48WKVL19e8+bN0+nTp5Wamqp27dpJOteumTp1qg4dOqT169fr1Vdf1VtvvWX1ZJTOtYWaNm2arb5/+ctfNGfOHBUvXlypqamKi4vTDz/84DeHxr///W+99tprWr16tXw+n3766Sd179491yAir22Oi713ZmamZsyYoQ8//NBqW11zzTVKTk7Wpk2b/IZQLFy4UMWKFVNaWpqeeuopvx4W0dHRSklJ0UcffaTdu3fr559/Vr169bRixQrZbDZJ525i/fDDDxoyZIg2btyorVu36p///Kdat24t6Vy7ZdasWYqLi1NKSorS0tLy9LrFixfX/Pnzdffdd+vnn3+Ww+HQvHnz1LlzZ2seCyCYIsyFfZMBoICdOHFCCQkJWrduXY7DYQJp9OjR6tevn1+3TQAAUDRkZmaqQoUKVm+D/DRnzhzddtttOnPmjN9knwDyhh4VAIJq69at6tmzZ76HFAAAoGjbvXu32rVrp5YtWwa7KgAugXgPQFA1bNhQH374YYG8V9WqVS86yzcAACi8KleurFGjRhXIe11zzTXq3LnzRVcNAZA7hn4AAAAAAICQQcQHAAAAAABCBkEFAAAAAAAIGQQVAAAAAAAgZBBUAAAAAACAkEFQAQAAAAAAQgZBBQAAAAAACBkEFQAAAAAAIGQQVAAAAAAAgJBBUAEAAAAAAELG/wGbQtIWFDklSwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 1280x1280 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from matplotlib import pyplot as plt\n",
    "\n",
    "fig, axes = plt.subplots(\n",
    "    nrows=2,\n",
    "    ncols=2,\n",
    "    figsize=(12.8, 12.8),\n",
    "    )\n",
    "\n",
    "axes[0, 0].set_title('Bias')\n",
    "axes[0, 0].hist(np.mean(bias_vars / true_bias_var, axis=0).flatten(),\n",
    "        bins=100,\n",
    "        histtype='step',\n",
    "        color='k',\n",
    "        )\n",
    "\n",
    "axes[0, 1].set_title('Dark')\n",
    "axes[0, 1].hist(np.mean(dark_vars / true_dark_var, axis=0).flatten(),\n",
    "        bins=100,\n",
    "        histtype='step',\n",
    "        color='k',\n",
    "        )\n",
    "\n",
    "axes[1, 0].set_title('Flat')\n",
    "axes[1, 0].hist(np.mean(flat_vars / effective_flat_var, axis=0).flatten(),\n",
    "        bins=100,\n",
    "        histtype='step',\n",
    "        color='k',\n",
    "        )\n",
    "\n",
    "axes[1, 1].set_title('Total')\n",
    "axes[1, 1].hist((empirical_err / measured_err).flatten(),\n",
    "        bins=100,\n",
    "        histtype='step',\n",
    "        color='k',\n",
    "        )\n",
    "\n",
    "for ax in axes.flatten():\n",
    "    ax.axvline(1, c='r', ls='-')\n",
    "    ax.set_xlabel('Measured / expected')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ba1508f",
   "metadata": {},
   "source": [
    "Clearly, `phoptic` is correctly propagating the errors from each calibration. However, this doesn't prove that the calibrations have actually been applied correctly. To verify this, we can check that the mean flux of the calibrated images is close to the expected value:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "c591828b",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:38.254370Z",
     "iopub.status.busy": "2026-06-16T12:21:38.254295Z",
     "iopub.status.idle": "2026-06-16T12:21:38.272036Z",
     "shell.execute_reply": "2026-06-16T12:21:38.271609Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9999779237287508\n"
     ]
    }
   ],
   "source": [
    "print(np.mean(calibrated_images / science_flux))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "f52323de",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:38.272829Z",
     "iopub.status.busy": "2026-06-16T12:21:38.272755Z",
     "iopub.status.idle": "2026-06-16T12:21:38.288975Z",
     "shell.execute_reply": "2026-06-16T12:21:38.288645Z"
    },
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [],
   "source": [
    "assert np.isclose(np.mean(calibrated_images / science_flux), 1., rtol=0.01, atol=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "03ec43d8",
   "metadata": {},
   "source": [
    "Another reassuring result! Comparing the calibrated fluxes graphically:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "d75e9d7a",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-06-16T12:21:38.290158Z",
     "iopub.status.busy": "2026-06-16T12:21:38.290087Z",
     "iopub.status.idle": "2026-06-16T12:21:38.460570Z",
     "shell.execute_reply": "2026-06-16T12:21:38.460209Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAHGCAYAAACmS4sdAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlcelbwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAR99JREFUeJzt3Xl8jOf+//F3KJEhojJ2xlbLKEqilmqrB6VUUV0cbbpSWj2l2tNFq6eni+p6qpyuqIP0VFutrRxR20FxWvvSaYuSUbFkEhLNRERy/f7wdf8MSUxI3Flez8cjj4e55r7u+cwtk7xz39d9XSHGGCMAAACblLG7AAAAULoRRgAAgK0IIwAAwFaEEQAAYCvCCAAAsBVhBAAA2IowAgAAbEUYAQAAtrrM7gIAFL6kpCTNmDFDvXv3VtOmTXNtS0hI0Jdffql+/fqpYcOGdpZ80TIyMvThhx/q2muvVbt27c67fXx8vL7//nslJibqzLkg+/TpoyuuuKIwSwVKPc6MAEVUZmamvv/+e02dOlWffPKJ/vOf/+j333+/oH0dOHBAo0aN0saNG/Ns++233zRq1Cjt2LHjousPxuHDhzV+/Hjt2rWrwPednp6uUaNGacmSJefddvHixXK73frqq6+0Z88e7d27Vz/88INGjRqlzZs3F3htAAJxZgQogqZOnarnnntO5cqVU5cuXVSxYkV99tlnWrNmjXr06KFp06bJ6XRe1Gs4nU6NHDlSzZo1K6Cq88/r9WrUqFGqW7eurWcfxo0bp44dO2r27NlW2/r16/X555/bVhNQmhBGgCLmtdde0/PPP69XX31Vzz77rMqWLWs9t3XrVt111106ePDgRYeRmjVravz48RdZbcnw66+/qmfPnnaXAZRahBGgCPnll1/0wgsv6K677tLzzz9/zvOtW7fWf//7X2VnZ0uSjh8/ro8++sh6PjQ0VA0aNNCf/vQnVahQIc/XymnMyJkyMjK0cOFCHT58WG3atFGHDh0Cnj9zfEndunUVFxenPXv2aODAgapcufJ569q1a5d15mHBggXWJahrrrlG7du3t/r6/X4tX75cXq9XERER6tatm2rUqHFOvbt379aKFStUvnx59erVS5dddv4fbzt27NB3332n5ORkeTweK5w9+OCDOW7/66+/auHChbr77rtVrVo1q/30/8P111+vqKgoSdJnn32mlJQUDRs2LCBQbtmyRcuXL1ePHj3UokWL89YIlAaMGQGKkE8//VTZ2dl67LHHct0mMjLS+kWYlZWlvXv3Wl/r1q3TI488ooYNG2r79u15vlZOY0ZOO3TokP70pz/pq6++0sKFC3Xdddfp7rvvVlZWlrXN6fEly5cv15/+9CfNmDFDsbGx8nq9QdWVlpamAwcOSJISExOtbVNSUqzXWLhwoRo0aKC//vWvWr9+vaZOnaqGDRvq008/Daj37bffVrNmzfSvf/1Ly5YtU/fu3bV27drzHG3p2LFj2rt3r7Kzs61/7927VydPnsxx+40bN2rUqFHav39/QPsff/yhUaNGadmyZVZbtWrV9NhjjwWEysTERPXp00exsbFq3LjxeesDSg0DoMj405/+ZEJCQsyJEycueB8nTpwwXbp0MZ06dbLatm3bZiSZzz//PM+2VatWGUmmRYsWZteuXVb7ggULjCTz5ptvnrOt2+02O3fuNMYYk5mZaRITE4Ou68cffzSSzFdffXXO9jt27DAVKlQw999/vzl58qTV/u6775qyZcuaTZs2GWOMWblypZFkXnjhBWub1NRU06dPHyPJjBs37rzHLCIiwgwePDigLafaPv/8cyPJeu3TEhMTjSTz1ltvBbS/+OKLJiQkxMydO9dkZWWZ7t27m8svv9zs2bPnvDUBpUmxvEyzY8cOHThwQF26dFG5cuUueD9//PGHtm3bplq1aqlBgwYFVyBwgY4cOaKKFSvm+/t606ZN2rJli1JTU5Wdna3w8HAtWrRImZmZF/QZ6d27d8Bf7r1791aXLl300Ucf6amnngrYtlevXtbg08suuyxgLMvF1DVhwgRlZmZq/PjxAZc5HnvsMY0dO1ZTpkzRxIkTNWnSJFWuXFmjR4+2tgkPD9fgwYP17bff5vu9F6S//e1vWrt2re677z7ddtttWrp0qebPn8/PG+AsxSqMzJkzR2+++aZ2796tw4cPKzEx8YIH8b322mt67bXX1Lx5cx0/flytWrXS9OnTLyrcABcrPDxcfr9fWVlZAb+Ac3PkyBH17dtXmzdvVteuXVW3bl2VK1dOqampOnnypJKTk3McX3E+V111VY5t//3vf5Wenq6wsDCrvXXr1oVS14YNG1S1alXNmDFDkmSMsb7CwsL0888/S5K2b9+uJk2aBNSU23u41MqUKaPPPvtMrVu31pQpUzR69GjdfPPNdpcFFDnFKozs3btXb731lo4dO6ZevXrlul12dra2b9+u9PR0ud1uVa5cOeD5iRMnauzYsVqyZIk6deokSZo1a5bS09MJI7BVmzZttGrVKm3fvj2oX6Zjx47V+vXrtXXrVjVp0sRqf/7557Vy5cqAybvyI6cgdHpA6OnBs6dFRkYWSl2ZmZkKCQnJcQ6SAQMGqFGjRlY9OX1ugxnAml+nj8vZxyA9PT3XPnv37lVycrKkU2NxAJyrWIWRxx9/XJK0aNGiXLfZuHGj7rzzToWEhKhq1ar6+eef9dJLL1l9s7Ky9Oqrr+ovf/mLFUQk6fbbby/M0oGg3HPPPZo4caKmTJmiCRMm5LhNenq6Tpw4oYiICG3atEnNmzcP+IUvnbo8cjFOn3U4k8fjUa1atVSxYsXz9g+2rjJlch9D36JFC+3Zs0fvvPNOnmeJmjRpouXLl+vkyZMBAcTj8Zy3zvyqXbu2pFOTtZ0pp+MlScnJybr99tvVrFkz3X333XrmmWd07bXX6oEHHijw2oDirETdTZOWlqY+ffpo1KhR2rlzp/73v/9p1apVGjNmjP73v/9JOjXe5PDhw+rdu7f279+v77//3hrRD9jt6quv1vDhw/XBBx9o5syZ5zx/4MAB3Xzzzdq3b58kqUGDBoqPj1dqaqq1zcqVK7Vu3bqLqmPWrFnWX/PSqRARFxene++9N6j+wdZVvXp1SaduMz7bo48+qj/++EOvvPLKOc8dPXrUOmNyzz33KCkpSVOmTLGeP3nypCZNmhRUrflx5ZVXqmLFipo3b955X8sYo3vuuUdHjhzRrFmz9PTTT2vgwIF69NFHtXXr1gKvDSjOitWZkfOZN2+eUlNT1aJFCy1fvlzSqR8ITZo0UVxcnDp06GAFj5kzZ2revHlyuVzaunWr+vXrp2nTpnGZBrabOHGiIiMj9cADD+if//ynunbtqooVK8rj8ejrr7+Wy+VS1apVJUlPP/20vv76a1133XW64447tG/fPm3evFkjR47U3//+9wuu4aGHHlLv3r3VtWtX+f1+/etf/1LHjh31t7/9Laj+wdZVp04dXX311XrjjTfk8/lUsWJFa56Rzp07a8qUKRo+fLiWLl1qDVj/+eeftWbNGn3yySe64oor1LdvXz388MMaPny41qxZo8aNG2vJkiUaPHiwZs2adcHHICdVqlTR6NGjNWbMGKWlpalRo0ZatmyZ/vKXv+irr74K2Hbs2LFauHChvvnmG+sM0eTJk3X11Vfrjjvu0Pr16xUeHl6g9QHFVYk6M7Jz506FhIRo7Nix1tdrr72myMhIRURESJIVNnbt2qXffvtNa9eu1bZt27Rw4UK9++67dpYPSDp16eLll1+W1+vVsGHDJEkpKSnq0KGDFi9erB07dliXC5o1a6ZffvlFgwcPVlpamq655hotX75cN9xwg0aOHGldUslp6vec2urUqaORI0eqR48emjNnjpxOpypUqKCPPvpI//3vf+VwOM7Z9vTYjTMFW1dISIji4uL01FNP6ejRo+fMM3L//fcrPj5e999/v7KyslSmTBkNGDBAHo8nYMbUDz/8UEuWLFG9evUkSe+//77uvPNOjRw5UldfffV5j/kjjzyiG2+8MaCtRo0aGjly5DmXmp5//nktXrxYderU0WWXXaYpU6aod+/eGjlypKKjoyVJqampOnLkiD7++GPdeuutVt9KlSrp66+/Vq9evfT999+fty6gtAgxFzrCzUaLFi1Sr169zrmb5p133tHYsWMDTi+fbc+ePWrUqJGmTp2q+++/32rv37+/jDGaO3duYZYOAADOUqLOjHTt2lVHjhwJuJ4ryZpdUZIaNmyo5s2bKz4+PmCb+Pj4C7oFEgAAXJxiNWZk586dio+Pt5b0XrlypSpXrqzWrVurevXqatu2rUaOHGmNWm/VqpX27Nmj6dOn67333tN1110n6dQZlEGDBiksLEzNmzfX7NmztXPnTv373/+28d0BAFA6FavLNJ988om+/PLLc9pffPFFK2hI0ty5c/XNN98oOTlZTZo00YMPPqiWLVsG9Pn+++81efJkJSUlqUmTJnrssceYFREAABsUqzACAABKnhI1ZgQAABQ/hBEAAGCrYjGANTs7WwkJCQoPD1dISIjd5QAAgCAYY3Ts2DHVrl07z+UfikUYSUhIsCYzAgAAxcu+fftUt27dXJ8vFmHk9JTJ+/btO2cFXgAlUFqa9H+zzCohQQpicT4ARU9qaqrq1at33qUPikUYOX1ppnLlyoQRoDQ4c5XeypUJI0Axd74hFgxgBQAAtiKMAAAAWxFGAACArQgjAADAVoQRAABgK8IIAACwFWEEAADYijACAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWhBEAAGCrYrFqL4DSZd++far3f//evHmzssPC5HQ65XK5bK0LQOEgjACwndfrlc/nkyQlJiYq5tZblfh/z3W+9lr5JTkcDnk8HgIJUAIRRgDYyuv1yu12y+/3W23OsDDr39+vXq0de/cqJiZGPp+PMAKUQIQRALby+Xzy+/2KjY2V2+2WJFVzOKT/+3ebNm2U/X/hxOPxWP24bAOUHIQRAEWC2+1WVFTUqQdpaQHPOZ1OORwOxcTEWG1ctgFKDsIIgCLP5XLJ4/FY40o8Hg+XbYAShDACoFhwuVwED6CEIowAuKTOvHNGChwHAqB0IowAuGRyunNGOjX+w+l02lQVALsRRgBcMjndOSNxZwxQ2hFGAFxyAXfOXISzL/EQaoDiiTACoNjJ6VZfidt9geKKMAKg2Dn7Vl+J232B4owwAqBY4lZfoOQoY3cBAACgdCOMAAAAW3GZBkChOnOSMyY4A5ATwgiAQpPTJGdMcAbgbIQRAIUmp0nOmAsEwNkIIwAKXUFNcgagZGIAKwAAsBVnRgCUKGcOkuWSEFA8EEYAlAg5TRHP9PBA8UAYAVAinD1FPNPDA8UHYQRAicEU8UDxxABWAABgK8IIAACwFWEEAADYijACAABsxQBWAAXmzEXxJBbGAxAcwgiAApHTongSC+MBOD/CCIACkdOieBKzoAI4P8IIgALFongA8osBrAAAwFacGQFQop09iJbLRkDRQxgBUCLltHCexOJ5QFFEGAFQIp29cJ7E4nlAUXVBYSQ1NVW//vqrGjVqpKpVqwbVJzk5Wfv371f9+vVVuXLlC3lZAMgXFs4Diod8DWDdvXu3hg0bpmbNmunqq6/WwoULz9vHGKORI0eqVq1a6t+/v6pXr65XXnnlggsGAAAlS77CyKZNm9S2bVvt3Lkz6D6TJk3S1KlTtX79eu3evVtxcXF6+eWXNW/evHwXCwAASp58Xaa5/fbb8/0CkyZN0u23365WrVpJkrp06aKuXbtq8uTJ6tu3b773BwAASpZCnWckKytLW7Zs0dVXXx3Q3qFDB23cuLEwXxoAABQThXo3zbFjx5SZmanIyMiAdqfTqaSkpFz7ZWRkKCMjw3qcmppaaDUCAAB7FeqZkXLlyklSQLCQpPT0dOu5nIwbN04RERHWV7169QqzTAAAYKNCDSMVK1ZU1apVlZCQENCekJCQ5+12o0ePVkpKivW1b9++wiwTAADYqMDDSEJCgjZv3mw97tatmxYsWGA9zs7O1sKFC9WtW7dc9xEaGqrKlSsHfAEoerxerzZu3KiNGzeeM+06AAQrX2NGjh07pl9++cV6vGfPHq1fv17VqlVT/fr1JUmffPKJxo8fr6NHj0qSxowZo44dO2rEiBG65ZZbNH36dCUnJ+vJJ58suHcB4JLzer1yu93y+/1Wm8PhkNPptLEqAMVRvs6MeDwePfzww3r44YcVHR2tuXPn6uGHH9aMGTOsbWrXrq22bdtaj1u3bq1Vq1YpMTFRf//731W2bFmtWbOGWRGBYs7n88nv9ys2NlYbNmzQhg0bWPMFwAXJ15mR9u3ba/369XluM3ToUA0dOjSgLTo6Wp9//nn+qwNQ5LndbkVFRdldBoBirFAHsAIAAJwPYQQAANiKMAIAAGxVqDOwAkBRdOZtyE6nk0G3gM0IIwBKDafTKYfDoZiYGKvN4XBwFxBgM8IIgFLD5XLJ4/HI5/NJOnWGJCYmRj6fjzAC2IgwAqBUcblcBA+giGEAKwAAsBVhBAAA2IowAgAAbEUYAQAAtiKMAAAAWxFGAACArQgjAADAVoQRAABgKyY9AxAUr9drzVwqBa7vAgAXgzAC4Ly8Xq/cbrf8fn9Au8PhkNPptKkqACUFYQTAefl8Pvn9fsXGxsrtdlvtrHgLoCAQRgAEze12Kyoqyu4yAJQwDGAFAAC2IowAAABbcZkGQKl39p1BjIUBLi3CCIBSy+l0yuFwKCYmJqDd4XDI4/EQSIBLhDACoNRyuVzyeDznzJ8SExMjn89HGAEuEcIIgFLN5XIROgCbMYAVAADYijACAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWhBEAAGArwggAALAVYQQAANiKMAIAAGxFGAEAALZibRoAOfJ6vdYCch6Px+ZqAJRkhBEA5/B6vXK73fL7/Vabw+GQ0+m0sSoAJRVhBMA5fD6f/H6/YmNj5Xa7JUlOp5PVbQEUCsIIgFy53W5FRUXZXQaAEo4BrAAAwFaEEQAAYCvCCAAAsBVjRgAgB2fezszgXaBwEUYA4AxOp1MOh0MxMTFWm8PhkMfjIZAAhYQwAgBncLlc8ng8ARO+xcTEyOfzEUaAQkIYAYCzuFwuggdwCTGAFQAA2IowAgAAbEUYAQAAtiKMAAAAW+U7jEyfPl1NmzZVaGioWrZsqfnz5+e5fXp6ukaNGiWXy6UKFSqoUaNGeuGFF5SVlXXBRQMAgJIjX2Fk0aJFGjx4sMaMGaODBw/qwQcf1IABA7Rhw4Zc+zz77LP66quvNHv2bB09elSffvqpJkyYoLfeeuuiiwcAAMVfvsLI22+/rVtuuUX33nuvLr/8cj3xxBNq06aNxo8fn2ufjRs36qabblJ0dLQqVKigG264Qddee22eAQYAAJQeQYcRY4zWrVunG264IaC9W7duWrNmTa79/vznPysuLk6bNm1SRkaGVq5cqXXr1mnQoEEXXDQAACg5gp707NixY0pLS1O1atUC2qtXr66DBw/m2u/RRx/V7t27FRUVJUkqU6aMXnvtNQ0YMCDXPhkZGcrIyLAep6amBlsmAAAoZi76bhpjjEJCQnJ9/plnntHXX3+tdevWKS0tTd99953efPNNTZgwIdc+48aNU0REhPVVr169iy0TAAAUUUGHkfDwcFWsWFGJiYkB7YmJiapRo0aOfbKzszVx4kQ98cQT6tChgxwOh7p27aqhQ4fq3XffzfW1Ro8erZSUFOtr3759wZYJAACKmaAv04SEhKhjx45avny5RowYYbUvW7ZM11xzTa59ypYte86ZE2OMypYtm+trhYaGKjQ0NNjSAFwkr9drLQwnnVocDgAulXwtlPfXv/5Vt9xyi6ZPn65bbrlF//rXv7Rp0yZ9+OGH1jZ///vfNX78eB09elQhISHq37+/3nnnHXXo0EGtW7fWunXrNGnSJA0ePLjA3wyA/PN6vXK73fL7/QHtDodDTqfTpqoAlCb5CiM33XSTpkyZoldeeUUPPfSQmjRpoq+//lrR0dG59vnggw/0t7/9TQMHDtShQ4dUp04dDR8+XGPGjLno4gFcPJ/PJ7/fr9jYWLndbqvd6XSyci2ASyLEGGPsLuJ8UlNTFRERoZSUFFWuXNnucoASZePGjYqOjtaGDRusu95sl5YmVap06t9//CFVrGhbKUXy+ADFRLC/v1mbBgAA2IowAgAAbEUYAQAAtsrXAFYAKK3Ovt2ZAb5AwSGMAEAenE6nHA6HYmJiAtodDoc8Hg+BBCgAhBEAyIPL5ZLH4zlnUriYmBj5fD7CCFAACCMAcB4ul4vQARQiBrACAABbEUYAAICtCCMAAMBWhBEAAGArwggAALAVYQQAANiKMAIAAGxFGAEAALYijAAAAFsRRgAAgK0IIwAAwFaEEQAAYCsWygNKIa/Xa61C6/F4bK4GQGlHGAFKGa/XK7fbLb/fb7U5HA45nU4bqwJQmhFGgFLG5/PJ7/crNjZWbrdbkuR0OuVyuWyuDEBpRRgBSim3262oqCi7ywAABrACAAB7EUYAAICtCCMAAMBWjBkBgAt05m3RDAIGLhxhBADyyel0yuFwKCYmxmpzOBzyeDwEEuACEEYAIJ9cLpc8Hk/AxHExMTHy+XyEEeACEEYA4AK4XC6CB1BAGMAKAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWhBEAAGArwggAALAVYQQAANiKMAIAAGxFGAEAALYijAAAAFsRRgAAgK0IIwAAwFas2guUcF6v11rqXjq13D0AFCWEEaAE83q9crvd8vv9Ae0Oh0NOp9OmqgAgEGEEKMF8Pp/8fr9iY2PldrutdqfTKZfLZWNlAPD/EUaAUsDtdisqKsruMgAgRwxgBQAAtiKMAAAAW3GZBgAKyNl3KjE2BwgOYQQALpLT6ZTD4VBMTExAu8PhkMfjIZAA55HvyzQ//PCD7rzzTnXo0EExMTFBzVmQkZGhf/zjH+revbtuvPFGTZ069YKKBYCiyOVyyePxaMOGDdZXbGys/H5/wBwvAHKWrzMjW7ZsUZcuXfTII4/okUce0bRp09S5c2dt2rRJ9evXz7HPiRMn1L17d6WkpOjll19WtWrVNH36dIWHh+v2228vkDcBAHZzuVycAQEuUL7CyKuvvqoOHTroH//4hySpS5cuat68ud555x1NmDAhxz4TJkzQli1btHPnTtWoUUOS1LlzZ504ceIiSwcAACVBvi7TLFu2TL179/7/ncuUUe/evbV06dJc+3z22Wfq37+/FUROK1++fD5LBQAAJVHQYSQtLU3JycmqVatWQHvt2rW1b9++XPt5PB5deeWVeu6559SpUyf169dPM2bMkDEm1z4ZGRlKTU0N+AIAACVT0GEkMzNTkhQaGhrQHhoaaj13NmOMTpw4obFjx0qS3n77bfXp00ePPvqoXn311Vxfa9y4cYqIiLC+6tWrF2yZAACgmAk6jISHh6t8+fJKSkoKaE9KSlJkZGSOfUJCQlS1alVFRUXptddeU+fOnfXQQw/piSee0Pvvv5/ra40ePVopKSnWV15nXgAAQPEWdBgpW7asrrrqKv3www8B7WvXrlV0dHSu/a6++mpVr149oK169eo6duxYrpdqQkNDVbly5YAvAABQMuVrAOvQoUM1a9Ysbd68WdKpAa3Lly/XQw89ZG3z0Ucf6brrrrMeP/LII1q0aJF++uknSdLRo0f16aefqkePHgoJCSmAtwAAAIqzfN3aO3jwYHk8HnXo0EE1a9bU4cOH9corr6hPnz7WNgcPHtS2bdusx3379tWYMWPUsWNH1axZUwkJCbr++uv18ccfF9y7AAAAxVa+wkhISIjeeecdvfDCCzpw4IDq1aunSpUqBWzz8MMP67bbbgtoe/rpp/WXv/xFXq9XtWrVUkRExMVXDgAASoQLWpumSpUqqlKlSo7P1axZUzVr1jyn3eFwqHnz5hfycgAAoATL99o0AAAABYkwAgAAbHVBl2kAFF1er9daKTaYVbUBwG6EEaAE8Xq9crvd8vv9VpvD4ZDT6bSxKgDIG2EEKEF8Pp/8fr9iY2PldrslSU6nk6XtARRphBGgBHK73YqKirK7DAAICgNYAQCArQgjAADAVlymAYBCdOYdTYzfAXJGGAGAQuB0OuVwOBQTE2O1ORwOeTweAglwFsIIABQCl8slj8cTMOdLTEyMfD4fYQQ4C2EEAAqJy+UieABBYAArAACwFWEEAADYijACAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWhBEAAGArwggAALAVYQQAANiKMAIAAGxFGAEAALYijAAAAFuxai9QjHm9XmuJeunUMvUAUNwQRoBiyuv1yu12y+/3B7Q7HA45nU6bqgKA/COMAMWUz+eT3+9XbGys3G631e50OuVyuWysDADyhzACFHNut1tRUVF2lwEAF4wwAgCX0NnjejiTBRBGAOCScDqdcjgciomJCWh3OBzyeDwEEpRqhBEAuARcLpc8Hs85dz/FxMTI5/MRRlCqEUYA4BJxuVyEDiAHTHoGAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWhBEAAGArwggAALAVYQQAANiKMAIAAGxFGAEAALYijAAAAFsRRgAAgK0IIwAAwFas2gsUI16v11qC3uPx2FwNABQMwghQTHi9Xrndbvn9fqvN4XDI6XTaWBUAXDzCCFBM+Hw++f1+xcbGyu12S5KcTqdcLpfNlQHAxSGMAMWM2+1WVFSU3WUAQIHJdxjx+Xz69NNPFR8fryZNmmjw4MEKDw8Pqu/27ds1fvx4tWvXTg8//HC+iwWAkujM8T+c7UJplK+7aRISEtSmTRvFxcWpQYMG+vzzz9W+fXulpqaet6/f79ef//xnLViwQEuWLLngggGgpHA6nXI4HIqJiVF0dLSio6Pldrvl9XrtLg24pPIVRl555RVVqVJFixYt0lNPPaUlS5boyJEjGj9+/Hn7jhgxQt27d1d0dPSF1goAJYrL5ZLH49GGDRu0YcMGxcbGyu/3W3dMAaVFvsLI/Pnzddttt6lcuXKSpPDwcN1yyy2aP39+nv2++OIL/e9//9Prr79+4ZUCQAnkcrkUFRWlqKgoa2AyUNoEPWbk+PHj2r9/vxo0aBDQ3rBhQ82aNSvXfr/99ptGjBih7777ThUqVAjqtTIyMpSRkWE9DuYyEAAAKJ6CPjOSnp4uSecMVg0PD7eeO1tmZqYGDRqkZ599Vq1btw66qHHjxikiIsL6qlevXtB9AQBA8RJ0GKlUqZLKlCmjI0eOBLQnJyercuXKOfaJi4vT5s2btX37dg0ZMkRDhgzR1q1btWHDBg0ZMkRJSUk59hs9erRSUlKsr3379uXjLQEAgOIk6Ms05cqVU9OmTfXTTz8FtO/YsUMtW7bMsU+rVq30/vvvB7T9+OOPKleunDp27KjQ0NAc+4WGhub6HAAAKFnyNc/IoEGD9NFHH+m5555TtWrVtHv3bi1cuFDvvvuutc28efO0ZMkSTZgwQfXr19eQIUMC9jFnzhxVqFDhnHYAAFA65etumqeeekrNmjVT27Ztddttt6lTp07q2bOnHnzwQWubjRs3avr06QVeKAAAKJnydWYkLCxMS5cu1apVq+T1evXUU0+pY8eOAdv07dtXjRs3znUfI0aMUNmyZS+sWgAAUOLkezr4MmXKqEuXLrk+f/p++dz06NEjvy8JAABKsHxdpgEAAChohBEAAGArwggAALBVvseMALg0vF5vwIJpZy4zDwAlCWEEKIK8Xq/cbrf8fn9Au8PhkNPptKkqACgchBGgCPL5fPL7/YqNjQ1YydXpdMrlctlYGQAUPMIIUIS53e48b5VHyXT2JTlCKEo6wggAFBFOp1MOh0MxMTEB7Q6HQx6Ph0CCEoswAgBFhMvlksfjOWfgckxMjHw+H2EEJRZhBACKEJfLRehAqcM8IwAAwFaEEQAAYCvCCAAAsBVhBAAA2IowAgAAbEUYAQAAtiKMAAAAWxFGAACArQgjAADAVoQRAABgK8IIAACwFWvTAEWE1+u1Fkg7ewl5ACjJCCNAEeD1euV2u+X3+602h8Mhp9NpY1UAcGkQRoAiwOfzye/3KzY2Vm63W5LkdDpZvRWWM8+W8b2BkoYwAhQhbrdbUVFRdpeBIsTpdMrhcCgmJsZqczgc8ng8BBKUGIQRACjCXC6XPB5PwHiimJgY+Xw+wghKDMIIABRxLpeL4IESjVt7AQCArQgjAADAVoQRAABgK8IIAACwFWEEAADYijACAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBWTAcP2MDr9VprjUiBK7ICQGlDGAEuMa/XK7fbLb/fH9DucDjkdDptqgoA7EMYAS4xn88nv9+v2NhYud1uq93pdLIYGoJ29tk0vn9QnBFGAJu43W5FRUXZXQaKGafTKYfDoZiYmIB2h8Mhj8dDIEGxRBgBgGLE5XLJ4/GcM+YoJiZGPp+PMIJiiTACAMWMy+UidKBE4dZeAABgK8IIAACwFWEEAADYijACAABsRRgBAAC2IowAAABbEUYAAICt8j3PSFZWlpYuXar4+Hg1adJEXbp0UUhISJ59EhIStGbNGmVmZio6OlpNmza94IIBAEDJkq8zI2lpabr++us1bNgwrVixQnfddZduueUWZWZm5trnL3/5izp37qyZM2dqzpw5atu2rf76179edOEAAKBkyNeZkTfeeEPx8fHaunWrqlatqvj4eLVs2VKTJ0/WI488kmOfrl27avz48brsslMvtXLlSnXp0kV9+/bV9ddff/HvACgGvF6vNX332QucAQXlzO8tFs5DcZKvMPLFF19o4MCBqlq1qiSpfv366tOnj7744otcw8iAAQMCHl9//fUqX768fv75Z8IISgWv1yu32y2/32+1ORwOOZ1OG6tCSZLT4nksnIfiJOgwkpmZqZ07dwYseS6dWnl06dKlQb/gokWLdOLECUVHR+e6TUZGhjIyMqzHqampQe8fKGp8Pp/8fr9iY2Otzw9/taIgnb14HgvnobgJOoykpaXJGKMqVaoEtF9++eVBh4X9+/dr8ODBiomJyTOMjBs3Ti+99FKwpQHFgtvtVlRUlN1loIRi8TwUZ0EPYK1QoYIk6dixYwHtqampcjgc5+1/6NAhde/eXVdeeaUmT56c57ajR49WSkqK9bVv375gywQAAMVM0GdGKlSooHr16mnPnj0B7b/99puaNGmSZ9/Dhw+ra9euqlu3rubOnavQ0NA8tw8NDT3vNgAAoGTI1629ffv21VdffWWN50hJSdH8+fPVt29fa5sffvhBkyZNsh4nJiaqa9euql27tubNm6ewsLACKh0AAJQE+bqb5oUXXtC3336r7t27q0ePHpo9e7Zq1KihESNGWNssXLhQ48eP10MPPSRJ6tmzp+Lj43XfffdpypQp1nbt27dX+/btC+htAACA4ipfYaRGjRratGmTpk+fLq/Xq6FDh+ree+8NGDPSvn17DR061HrcpUsXXXPNNYqPjw/YV6NGjS6ydAAAUBLkezr4yy+/XCNHjsz1+d69e6t3797W43fffffCKgMAAKUCC+UBAABbEUYAAICt8n2ZBkDezlyHRmItGtjn7O89Zv5FUUUYAQpQTuvQSKxFg0srp7VqJNarQdFFGAEKUE7r0Ej8RYpL6+y1aiTWq0HRRhgBCgHr0MBurFWD4oQBrAAAwFaEEQAAYCvCCAAAsBVhBAAA2IowAgAAbEUYAQAAtuLWXuAinTnjKrOtoqg783uU+W9QVBBGgIuQ04yrzLaKoiinWVmZkRVFBWEEuAg5zbjKX5sois6elZUZWVGUEEaAAsCMqygOmJUVRRUDWAEAgK0IIwAAwFaEEQAAYCvCCAAAsBUDWAGgFDt7bhzuBoMdCCNAPpw5wZnEJGcovnKad0Ri7hHYgzACBCmnCc4kJjlD8XT2vCMSc4/APoQRIEg5TXAmcVobxRfzjqCoIIwA+cQEZwBQsLibBgAA2IowAgAAbEUYAQAAtmLMCJCHM2/l5TZelBZnfq8zQBuXAmEEyEVOt/JyGy9KspzmHmHeEVwKhBEgFzndystfiSjJzp57hHlHcKkQRoDz4FZelCbMPQI7MIAVAADYijMjAIA8sZgeChthBPg/LIIHBGIxPVwqhBFALIIH5ITF9HCpEEYAsQgekBsGtOJSIIwAZ+DOGQC49AgjAIB8Y5ZWFCTCCEotpnoH8o9ZWlEYCCMolZjqHbgwzNKKwkAYQanEVO/AhWNQKwoaYQSlGgNWgYLBxGi4GIQRlApMaAYUDiZGQ0EgjKDEY0IzoPAwMRoKAmEEJR4TmgGFK7cxJNz+i2ARRlAi5XTbLuNDgEuD23+RX4QRlDjctgvYK7fbf1etWsXZSeSIMIJiL6fBqdy2C9jrzEs3DHLF+RBGUKzlNTj1uuuu44ccUATkNcj1zLMl/NFQel1QGDl48KD27dunRo0aKTIystD6AGcL5iyIxA81oKg5e5Ar40pwpnyFkezsbD3yyCOaNm2aGjdurN27d+vJJ5/U2LFjC7QPcNqZ4SMxMVEDBgzgLAhQAjCuBGfKVxj58MMP9cUXX2jz5s1q3ry51q5dqy5duig6OloDBgwosD4onc4+65FT+HA4HFq0aJGqVatmtfHDCiiegh1X8s033wR85s/Gz4DiL8QYY4LdODo6WlFRUZo0aZLV1qtXL5UtW1bffvttgfU5W2pqqiIiIpSSkqLKlSsHWy6KkLODxtnyOutx5g8ifuiUEmlpUqVKp/79xx9SxYr21oNLIpg/SHKSU2DhZ0XREOzv76DPjGRlZWnbtm166KGHAtrbt2+vyZMnF1gfScrIyFBGRob1OCUlRdKpN1XQDh48qIMHDxb4fvH/+Xw+xcTEKD09Pc/twsLC9PXXXwfcghsZGal69eoFbFcY3wcoYtLS/v+/U1OlrCz7asElU6VKFVWpUsV6fMUVV+iHH35QUlJSrn1O/3y56aabAtrDwsIUGxvLLf1BqlmzpmrWrFng+z398/p85z2CDiPHjh1TZmbmOYNPIyMjlZycXGB9JGncuHF66aWXzmk/+5cSSpb09HTddtttdpeBoqZ2bbsrQDHEz5Oi5dixY4qIiMj1+aDDSLly5SRJx48fD2hPT09X+fLlC6yPJI0ePVpPPPGE9Tg7O1vJycmKjIxUSEhIsCUXK6mpqapXr5727dtXqi9FcRxO4ThwDE7jOJzCcSiex8AYo2PHjqn2ef6oCDqMVKxYUZGRkdq/f39A+/79+3O9LnchfSQpNDRUoaGhAW1nnrorySpXrlxsvskKE8fhFI4Dx+A0jsMpHIfidwzyOiNyWpn87PDGG2/U/PnzrcdZWVlasGCBbrzxRqvN6/Vq3bp1+eoDAABKr3yFkRdeeEFbtmzRsGHDNH/+fA0aNEjHjh3Tk08+aW3z6aefBgwkCqYPAAAovfIVRlq0aKG1a9fqxIkTGj9+vKpWrap169apTp061jYul0udOnXKVx+cujT14osvnnN5qrThOJzCceAYnMZxOIXjULKPQb7mGQEAACho+TozAgAAUNAIIwAAwFaEEQAAYKt8LZSH/Nm3b58OHTqkpk2bBnVPuDFGe/bs0ZEjR1SvXj1Vr1494Pndu3frwIEDAW0VK1ZU27ZtC7TugpSVlaUdO3YoJCREV155pcqUOX/+9fv92rVrlypUqKCmTZsW2H7tdPToUe3atUs1a9ZU3bp1g+pz6NAh7du3T/Xr1z9nkbDk5GT99NNP5/Rp3759nhMK2s3r9crr9SoqKkoOhyOoPnv27FFycrKaN2+uirmsURPMNkXJpk2bJCnoz+7Jkye1fft2lS9fXm63+5zJH1evXn1OnyuuuKJQpvcuKH/88Ye2bNmi+vXrB/2ZCKZPUlKSfvvtN9WrV69Iv//TfvvtNyUkJOTrs5tXn+3bt+vo0aMBbZGRkQErIRdJBgUuPT3dDBgwwISFhZnmzZubsLAwM2HChDz7bNmyxTRv3tzUqFHDREVFGYfDYQYMGGD8fr+1zbBhw0xkZKTp3Lmz9XXPPfcU9tu5YBs3bjT169c3derUMTVr1jSNGzc227Zty7PPG2+8YSpVqmRatWplateubdq2bWvi4+Mver92euutt0yFChWM2+02YWFhZuDAgSYjIyPX7ZOSkkzv3r1N5cqVTZs2bUx4eLgZMWKEyc7OtraZPXu2KVOmTMD3QufOnc3hw4cvxVvKt5UrV5qbb77ZREZGGklB/X+lpKSY7t27m0qVKpmmTZuaSpUqmenTp+d7m6Lk3XffNc2aNTNVqlQxV155ZVB9Vq9ebWrVqmVcLpepVq2aadGihdm1a5f1fGZmppFkWrZsGfC9MGfOnMJ6Gxdl3759Zvjw4aZmzZqmXLlyZty4cQXW529/+5sJDQ01LVq0MKGhoWbw4MEmKyuroN9CgYiLizPdu3c3VatWNZLMvn37CqRPt27dTN26dQO+F0aPHl0Yb6FAEUYKwbPPPmvq1q1rEhISjDGnfnFIMuvWrcu1zzXXXGN69uxpMjMzjTHGxMfHm4iICPPWW29Z2wwbNszcdttthVt8ATlx4oRp1KiRue+++4wxxmRnZ5s77rjDNG/ePNcfDnFxcaZMmTImLi7OGGPMyZMnzX333Wc6d+58Ufu104oVK0xISIj1nuLj40316tXNSy+9lGufe+65x7Ro0cIkJSUZY4zZv3+/qV27tvnwww+tbWbPnm0qVqxYuMUXoPfff9/Mnz/fbNiwIegwMmTIENOsWTOTnJxsjDFm0qRJ5rLLLjO//PJLvrYpSkaOHGl++ukn88wzzwQVRv744w9Ts2ZNM2LECGPMqc9Ez549zdVXX21tczqMLF++vLDKLlArVqww//znP01KSoqpU6dOUGEkmD5z5swx5cqVM99//70xxpiff/7ZREREmPfee6/A30NBeOedd0xcXJxZtmxZ0GEkmD7dunUzzzzzTGGUXKgII4WgRo0a5u9//3tAW8uWLc2wYcNy7dOkSRMzZsyYgLYrr7zSPP3009bjYcOGmZtvvtls2LDB7Nq1q0j+8j1t8eLFRlLAX3CbN282ksyqVaty7PP000+bpk2bBrStWbPGSDI7duy44P3a6d577zUdO3YMaPvrX/9q6tevn2uf2rVrm5dffjmgbfjw4aZt27bW49Nh5KeffjJbt2416enpBVp3Ydm0aVNQYeT48ePG4XCYf/7zn1Zbdna2qV27tnn++eeD3qaoCjaMfPnll6ZMmTLm0KFDVtuKFSsCjuHpMPLZZ5+Z9evXWyG2OAg2jATTp2/fvuamm24KaBsyZIi56qqrLqbEQrd8+fKgw0gwfbp162aGDx9ufvzxR+P1egPOqBZlRftCezGUkJCgQ4cOKTo6OqC9ffv21nXinLz88suaOnWqpkyZoiVLlujpp5+W3+/X8OHDA7aLi4vT/fffr06dOql+/fpauHBhobyPi7Vp0yZFRESocePGVttVV12l8uXL53ocqlatKp/Pp4yMDKvt9LpGGzZsuOD92mnTpk05fi/Ex8fryJEjOfapWrVqjus5bdu2TZmZmVZbWlqabrnlFt166626/PLL9corrxT8G7DJL7/8Ir/fH3DsQkJC1K5dO+v/OZhtirtNmzadM36sffv21nNnGjFihB544AHVqlVLt912m5KSki5prXbL7bO2ffv2gM9NafDpp59qyJAhat26tVq1aqUff/zR7pLOizBSwJKTkyWdGjB0psjISOu5nHTt2lXt27fX6NGj9fTTT+uTTz7R0KFDAxYU7Nmzp/bv36+tW7fqwIEDGjRokO644w7t3r27cN7MRTi9yvLZ8joO99xzjyTpzjvv1MKFCzVt2jSNGTNGl112mdXnQvZrp5zqPf04t3off/xxffrppxo3bpwWL16s5557TitXrtTJkyeVmpoqSapXr57Wr1+vXbt2adeuXfrmm2+sQFsSBPM5utDPWnGS0/dPWFiYwsLCrPcYEhKiKVOmKDExUVu3btWvv/6qbdu2aejQoXaUbJvcPmtZWVnW56Y0ePDBB5WYmKjNmzcrISFBrVq1Uv/+/ZWSkmJ3aXkijBSwcuXKSZKOHz8e0J6enp7rSGljjHr16iVJ+v3337Vx40Zt27ZN7733nsaOHWttd+utt1p/IZUtW1bjxo1TuXLlNG/evMJ4KxelXLly5xwDKe/jULt2bW3atEkNGzbUe++9p0WLFumzzz5Tdna2wsLCLni/dsqp3vT0dEnKtd7Bgwdrzpw52rFjh958801lZWXptddekyTrOERHRwf8FdirVy/169dPM2fOLIy3cckF8zm6kM9acZPT948xRidOnLDeY9myZfXggw9ad9jUr19fzz33nObOnZvjZ6WkupDPWkl01113qVKlSpJO/bwYP368EhIStHLlSpsryxthpIDVq1dPZcqUyfE0+5lnOc6UkJCgjRs36qGHHrI+NPXq1VPfvn3zDBply5ZVZGTkOa9VFNSvX18+n08nTpyw2tLS0pSSkpLrcZBOrW00fvx4xcXF6fPPP1dWVpays7PVsmXLi9qvXerXr5/j90K5cuXyvO2wd+/eio2N1ZIlS/TGG2/op59+UsOGDfO8HbZGjRpF8nvhQtSvX1+S8vwcBbNNcVe/fn0dOHBA5oxVOw4cOKCsrKw832ONGjWUlZWlgwcPXooyi4TcPmtVqlRReHi4TVXZLzIyUmXLli3yPxsIIwXM4XDommuuCQgRaWlpWrJkiW688UarbdeuXdY136pVq6pMmTL6/fffA/a1b98+a34JY4yV8k/buXOn4uPjrV/URUm3bt2UmZmpRYsWWW3z5s1TmTJl1LVrV6ttzZo1Ae/77Pf44YcfqmnTpurYsWO+9ltU3HjjjVq8eHHAOJi5c+fqhhtusP6yT0pK0urVq62AdfYxOHr0qGbOnKl7773XaktLSwvY5uTJk1qxYkWR/F4Ilsfj0bZt2yRJdevWVfPmzQM+R4cPH9batWutz1Ew2xQ3xhitXr3aChE33nijjhw5olWrVlnbzJ07VxUqVNB1110n6dzvBUlavHixqlSpEvT8HUXNoUOHtHr1amVnZwfd58Ybb9TChQuVlZVltc2dO7fYfi9I0rZt2+TxeILePiMjI+D9S9LSpUuVlZVV9H822Dp8toRasWKFKVeunHn22WfN3LlzTffu3U3jxo3NsWPHrG0GDx4cMJp+2LBhxul0mg8//NDExcWZJ5980oSEhJiFCxcaY4zJyMgwbrfbvPPOO2bRokXmk08+MQ0aNDBXX321OX78+CV/j8F45JFHTM2aNc306dPN1KlTTWRkpHniiScCtgkNDQ0YGX/ttdeaSZMmmbi4OPPQQw+ZSpUqmbVr1+Z7v0VFcnKycblcpnfv3mbevHnmiSeeMOXLlw94T1999VXAyPhly5aZAQMGmDlz5pgvvvjCtGnTxnTq1Cngjpk777zTPPXUU2bu3Llm1qxZpnv37qZKlSpm+/btl/w9BuP33383q1atMlOnTjWSzPTp082qVaus29+NMaZfv36mS5cu1uPZs2ebsmXLmldeecXMnj3bdOrUyVx11VXmxIkT+dqmKNm8ebNZtWqVufvuu03Dhg3NqlWrzKpVq6xb+tPT042kgNu4Bw0aZBo0aGD+/e9/m08++cSEh4cH3G318ccfmwEDBpgZM2aYhQsXmpEjR5rLLrvMfPTRR5f8/QXD7/db77tatWrm4YcfNqtWrQr43p00aZKRZP3MDKZPQkKCqV69urn99tvNvHnzzLBhw4zD4SiycxDFx8ebVatWmQkTJhhJ5ptvvjGrVq0KmCuoS5cupl+/fkH32bVrl4mKijIffPCBiYuLM++++66JjIw0t95666V+e/nGqr2F5Pvvv9f777+vQ4cOqVWrVnr22WcDTsu/8cYb2rFjh6ZPny7p1Iyi06dPV1xcnJKSklS/fn0NGTLEOiMgnTo9O3HiRG3atEmXX365rrvuOg0ZMsT6C7uoycrK0ocffqgFCxYoJCREffv21dChQwNmS+3WrZsefPBB3X333ZJOzdD5+uuva+fOnWrRooVGjBgRcOdMsPstSvbv32/9f9eqVUuPPfaYOnToYD2/YsUKjRkzRrNnz7bOhH377beaOnWqMjIydOONN+rhhx8OWDY8IyNDkyZN0rJly6zLWCNHjjxnptaiYubMmfrnP/95TvsTTzyhAQMGSJKee+45HTt2TBMnTrSe/+677zRp0iQlJyerXbt2euaZZ3T55ZcH7COYbYqKIUOG6Oeffz6nfcGCBYqIiNCJEyfUtWtXPfXUU+rXr58k6cSJE5o4caIWL16s8uXL67bbbtP9998f0P8///mPvvjiCx08eFCNGjXS4MGDz7mzpKjwer266667zmnv2LGj3n77bUmnvv9ff/11fffddwoLCwuqj3RqJt4333xTv/zyi1wul0aNGqWrrrqq8N7MRZg0aZKmTZt2Tvvf/vY39ejRQ5L02GOPKTw83BozFkyfX3/9VR988IE8Ho9q1qypXr16aeDAgefM2lvUEEYAAICtiuafkgAAoNQgjAAAAFsRRgAAgK0IIwAAwFaEEQAAYCvCCAAAsBVhBAAA2IowAqDQ/fLLL1q8eLHdZRRLq1at0ubNm+0uAyhUhBGgEKWlpWnmzJn68ssvldP8gnPmzNHMmTOVlJRkQ3WXzvz58/Xcc8+dd7uEhAQtWLDgElRUuJYuXart27cXyL7eeOMNxcbGFsi+gKLqMrsLAEqyxMREDRo0SJLkdDoDFvPbvn27BgwYIGOM1q5dq8jISLvKLDI+/vhj/frrr7r55pvtLuWivPjii+revXvRX5wMKCI4MwJcAtdff72mTJkS0DZ58mRdf/31OW6flJSkBQsWaMmSJUpOTg547vjx45o5c6Zmzpypr7/+Whs3bsxxddP09HStWLFC//nPf3T48GGrff/+/fr6668Dtj127JhmzpxprRjs8/k0c+ZMZWVlaf369Zo1a1bAPn7++WfNmTNH//vf/5SZmXnOa588eVLLly/Xd999p8TExPMcnf/v22+/VZ8+ffLcJrfXPnDggGbOnKkjR45YbcYYzZo1yzpLsXjxYnk8Hvl8Pi1ZskTLli3Lsf7s7GytW7dOc+fO1Y4dO3KsIzs7Wz/++KPmzp2r+Ph4q3316tXy+Xzavn279f90/PjxoPebkpKib7/9Vt9//738fn+exwIoMWxcpA8o8fbs2WMkmWnTppmwsDBz5MgRY4wxx48fN5GRkWb69OlGUsAqvh988IGJiIgw3bt3N926dTNVqlQxX3zxhfX80aNHzcCBA83AgQPNrbfeaurVq2c6duxojh49am2zbds2U716ddOuXTvTp08f06BBA2sV19mzZ5uKFSsG1OnxeAJWDl61apWRZG6++WZz1VVXmTvvvNPs2LHDHD9+3Nx5552mRo0apk+fPqZly5amRYsWZvfu3da+kpOTTZs2bUydOnXMTTfdZGrXrm26dOlioqOj8zxWCQkJ5rLLLjNJSUk5Pn++187MzDTt27c3AwYMsPq88cYbJjIy0uzfv98YY0x0dLS5/vrrrdrq1q1rWrdubXw+n9Vn7969pnXr1qZ58+bmlltuMXXq1DF9+/Y1GRkZAf+vV111laldu7bp3bu3adSokXn99deNMca8/fbbxul0miuvvNL6fzp69GhQ+/3xxx9NZGSkadWqlenatatp0qSJcbvd5sknn8zz2AHFHWEEKESnw8jatWvNtddea95//31jjDGff/65cbvdZufOnQFhZN26dSY8PDxg2fNvv/3WhIeHBywtfqaMjAzTuXNnM2bMGKttyJAh5s477wzY5ttvvzXG5C+MjBw5MmC75557znTq1MmkpaVZbQ899JDp2bOn9XjUqFGmZcuWJjU11RhjzG+//WbCw8PPG0Y++eQTc+211+b6fDCvvXPnTlOpUiUzadIks379elO+fHkze/Zs6/no6GhTqVIl89tvvxljjElNTTUtW7Y0I0aMsLa55pprzOOPP26ys7ONMcakpaWZ1q1bm3HjxlnbtG/f3vTo0cP4/X5jjDFZWVlm/vz51vOdO3c2L774YkD9wew3Ojra3H///dY2X331lZFEGEGJx5gR4BIZPHiwJk6cqOHDh2vy5MkaPHjwOdtMmzZNV1xxhX7++Wd5PB6ZU38wKDMzU+vXr1evXr2sbTdv3qz4+Hilp6erTp06+uGHH6znwsLC5PF4lJSUpMjISJUvX/6CxmGMGDEi4PHUqVPVt29fLVy40KqtevXqmjFjhrKzs1WmTBl9+eWXGj16tMLDwyVJDRs21B133KEtW7bk+Vrnu0QTzGtfccUVmjBhgh577DHVqFFD9913n/r37x+wn9tvv10NGzaUJIWHh2v48OF6+eWX9d577+nXX3/VmjVr9Oc//1lff/219TqNGzfW8uXL9eyzz8rj8eiHH37QunXrFBYWJkkqU6ZMnrUHs99du3Zpw4YNmjZtmrXc++23365GjRrledyAkoAwAlwid9xxh0aMGKHZs2dr9erV+ve//63U1NSAbfbu3auUlBTNmjUroL1fv36qWLGipFPjSXr06KGEhAS1bdtWlStX1q+//hqw/XPPPachQ4aobt26ioqK0k033aRHH31UVatWzVfNtWrVsv594sQJHThwQB6PR0ePHj2nvuPHj6t8+fJKSEhQgwYNAp5v2LBhnmEkIyNDS5cu1WuvvZbj88G8tsPhkCQ98MADGjt2rH7//Xe9+eab5+wrp9oOHjyojIwM7d27V5K0cuVKlS1b1tqmfPnyatq0qSTJ6/VKkvU4GPnZb071ASUdYQS4RCpWrKg///nPeuCBB3TzzTerevXq54SRypUrq3Hjxpo5c2au+3nvvfckSfv27dNll536CD/77LNatGiRtU3NmjX17bffKiUlRStXrtSbb76pzz//XD/99JPKlClzzoDX0wMsz3b6L3RJKleunCpUqKCBAwdq+PDhudYXERERMIhU0jmPz7Zs2TJVq1ZNV155ZY7PB/vakvT+++/L5/MpMjJSr7/+ul5//fU8azly5IjCw8MVGhqqypUrS5JeeukltWjRIsf9V6lSRdKpUHj55ZfnWctpwez39N1UR44csYJnTvUCJRF30wCX0KOPPqqbbrpJo0aNyvH5m266SStWrNDPP/8c0O7z+XTixAlJ0sGDB9W4cWMriJw8eVJz584N2H7//v2STgWDW265Ra+99po8Ho9SU1NVp04dpaenKyEhwdp++fLl5609JCREPXv21OTJk5WVlZXj60nStddeqzlz5liPT548qXnz5uW57/Ndogn2tX/66Sc99dRT+uCDDzRjxgy98847+u9//3vOa508edJ6/M0336hz586SpKioKFWvXl0fffTROTWcPl5t27ZVtWrVNH369IDnz7xrqFKlSgEBL5j9NmvWTE6nM+DY7d69+7yXt4CSgDMjwCV01VVX5XnW45577tE333yj6667TiNGjFCtWrW0bds2LVy4UBs3blT58uXVr18/9e/fXy+88ILq1q2rGTNm6ODBg6pfv761nzFjxig5OVndunVThQoV9Mknn6hHjx6qXLmy2rRpo1atWunOO+/UAw88oF9//VWff/55UPX/4x//UJcuXXTNNdcoJiZGxhitXLlS5cuX17///W9J0ssvv6xrrrlG9913nzp37qwvv/xSR48eVURERK77XbBggT7++OOLeu0TJ07o7rvv1oABA3TXXXdJkh5//HHdc8892rp1q3VGIyUlRT179tTAgQO1Zs0azZ8/X6tXr5Z06rLJ5MmTdccdd+jQoUPq3r27fD6f5s2bp3vvvVePPPKIypcvr48++kiDBg3S/v371bFjR23btk2JiYnWcWzXrp1mzpypZs2aKSwsTP379z/vfitUqKBXXnlFjz/+uA4fPqzq1atrwoQJ1lkVoCQLMSaHaSEBFIjExEQ99thjevXVV3XFFVec8/yhQ4c0cuTIgOezs7M1e/ZsLV26VJmZmWrTpo3uvfdea0CoJMXFxWnOnDnKzMzUddddpypVqmjt2rUBlyTmz5+vJUuWKCMjQ+3atVNMTIwqVKggSUpOTtaECRP0+++/68orr1SfPn30wgsv6IMPPlDVqlX1yy+/6MUXX9SMGTNUrly5gJqPHj2qadOmafPmzapSpYq6dOmifv36BVzS2bJli3UWo3PnzoqIiNDatWs1duzYc47Btm3b1KlTJyUlJSk0NDTP45nXa8fFxemzzz7TxIkTreBz4sQJPfbYY7rhhhs0aNAgtWvXTrfffrvq16+v1atXq0yZMho8eLDatGkT8Dq7du3SjBkztHfvXtWrV0/9+/dXu3btArbZvn27ZsyYocTERLVr105DhgxR+fLlJUl+v18fffSRtm3bpvT0dH388ceKiIgIar/z5s3T3Llzdfnll6t///5as2aN6tSpo7vvvjvPYwMUZ4QRALYaN26cfvjhB82ePbvQX+t0GHn22WcL/bUABI8xIwBsFRISoqFDh9pdBgAbMWYEgK0u5VmKnj175no3CwD7cJkGAADYiss0AADAVoQRAABgK8IIAACwFWEEAADYijACAABsRRgBAAC2IowAAABbEUYAAICtCCMAAMBW/w+ij4YvTJvd+gAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots()\n",
    "\n",
    "ax.hist(\n",
    "    calibrated_images.flatten() / science_flux,\n",
    "    bins=100,\n",
    "    histtype='step',\n",
    "    color='k',\n",
    "    )\n",
    "\n",
    "ax.axvline(1, c='r')\n",
    "ax.set_xlabel('Measured / expected')\n",
    "ax.set_title('Calibrated flux')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dbabc8bd",
   "metadata": {},
   "source": [
    "We can see that, unsurprisingly, the calibrated fluxes follow a Gaussian distribution centred on the expected value.\n",
    "\n",
    "That concludes this demonstration of `phoptic`'s calibration error propagation."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "opticam7",
   "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.14.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}