diff --git a/CHANGES.md b/CHANGES.md index 297625b..081b6d7 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -2,7 +2,9 @@ ## Version 0.3.0 (unreleased) -- No changes yet. +- The ``get_limmag`` function now accepts ``synphot.SourceSpectrum`` + instance, rather than a spectral model class, for more consistency with the + other methods and for more flexibility. ## Version 0.2.0 (2021-03-03) diff --git a/dorado/sensitivity/__init__.py b/dorado/sensitivity/__init__.py index 2d81b3f..f066fb0 100644 --- a/dorado/sensitivity/__init__.py +++ b/dorado/sensitivity/__init__.py @@ -10,7 +10,7 @@ from astropy import units as u import numpy as np from synphot.exceptions import SynphotError -from synphot import Observation, SourceSpectrum +from synphot import Observation from . import backgrounds from . import bandpasses @@ -75,13 +75,13 @@ def _amp_for_signal_to_noise_oir_ccd( return 0.5 * snr2 / signal * (1 + np.sqrt(1 + 4 * noise2 / snr2)) -def get_limmag(model, *, snr, exptime, coord, time, night): +def get_limmag(source_spectrum, *, snr, exptime, coord, time, night): """Get the limiting magnitude for a given SNR. Parameters ---------- - source_model : synphot.Model - The spectral model of the source. + source_spectrum : synphot.SourceSpectrum + The spectrum of the source. snr : float The desired SNR. exptime : astropy.units.Quantity @@ -99,11 +99,13 @@ def get_limmag(model, *, snr, exptime, coord, time, night): astropy.units.Quantity The AB magnitude of the source """ + mag0 = Observation(source_spectrum, bandpasses.NUV_D).effstim( + u.ABmag, area=constants.AREA) + result = _amp_for_signal_to_noise_oir_ccd( snr, exptime, - constants.APERTURE_CORRECTION * _get_count_rate( - SourceSpectrum(model, amplitude=0*u.ABmag)), + constants.APERTURE_CORRECTION * _get_count_rate(source_spectrum), ( _get_count_rate(backgrounds.get_zodiacal_light(coord, time)) + _get_count_rate(backgrounds.get_airglow(night)) @@ -113,7 +115,7 @@ def get_limmag(model, *, snr, exptime, coord, time, night): constants.NPIX ).to_value(u.dimensionless_unscaled) - return -2.5 * np.log10(result) * u.ABmag + return -2.5 * np.log10(result) * u.mag + mag0 def _exptime_for_signal_to_noise_oir_ccd( diff --git a/dorado/sensitivity/tests/test_sensitivity.py b/dorado/sensitivity/tests/test_sensitivity.py index 83ce122..d913a5e 100644 --- a/dorado/sensitivity/tests/test_sensitivity.py +++ b/dorado/sensitivity/tests/test_sensitivity.py @@ -31,7 +31,8 @@ def test_round_trip_snr_limmag(snr, exptime, ra, dec, time, night): coord=SkyCoord(ra * u.deg, dec * u.deg), time=time, night=night) - limmag = get_limmag(ConstFlux1D, snr=snr, **kwargs) + limmag = get_limmag( + SourceSpectrum(ConstFlux1D, amplitude=0*u.ABmag), snr=snr, **kwargs) snr_2 = get_snr(SourceSpectrum(ConstFlux1D, amplitude=limmag), **kwargs) assert snr_2 == approx(snr) diff --git a/example.ipynb b/example.ipynb index 2eabb79..e64b0e9 100644 --- a/example.ipynb +++ b/example.ipynb @@ -48,7 +48,7 @@ "output_type": "display_data", "data": { "text/plain": "
", - "image/svg+xml": "\n\n\n\n \n \n \n \n 2021-03-01T17:59:58.862923\n image/svg+xml\n \n \n Matplotlib v3.3.3, https://matplotlib.org/\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n\n", + "image/svg+xml": "\n\n\n\n \n \n \n \n 2021-03-18T19:44:39.674934\n image/svg+xml\n \n \n Matplotlib v3.3.3, https://matplotlib.org/\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n\n", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEYCAYAAACtEtpmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA3vElEQVR4nO3dd5hcZdn48e+9vdfsJpuyaUt6QltKqAnNUAIISFGsQF5EFFReAVFQfiigvgoookGQTkCkhSpEQglJSCOk9152s5tsyfZy//44Z5Zh2c22mTmzs/fnuuZi5plznvPMYbL3PF1UFWOMMeZQorwugDHGmPBnwcIYY0yHLFgYY4zpkAULY4wxHbJgYYwxpkMWLIwxxnTIgoUxxpgOWbAwxhjTIQsWxvQyIrJKRKZ09/2uHmcMWLAwvZCIbBWRYhFJ9ku7WkTmus9VRApanfMrEXnKff6WiNzZRr4XiMheEYkJ8kfoEVUdr6pzoeVenNHe+53Np728jPGxYGF6q2jghm6e+zhwpYhIq/RvAk+ramOPSmZMBLJgYXqr3wM3iUhGN859GcgGTvYliEgmcB7wRFczE5GbRWSXiFSKyDoROd1NHygi/xaRfSKyRUR+1Oq8rSJyk4h8JiLlIvKciCQcKk+/884QkSeBfGC2iBwUkZ+1ev9mEXmh1TXvF5EH/I9zn7fO62YR+Xercx8Qkfu7en9MZLBgYXqrxcBc4KaunqiqNcDzwLf8ki8F1qrq8q7kJSKjgeuBY1Q1FfgKsFVEooDZwHJgEHA6cKOIfKVVFpcC04DhwCTgO+3l2cbn+CawHZiuqimq+rtWh8wCzhGRVLes0e71nukoL5ygOc0XjN2mucvpRjA1kcGChenNbgd+KCI53Tj3ceAS3y95nMDxuO9NERnm1gj+KyIfiMj/iUhSG/k0AfHAOBGJVdWtqroJOAbIUdU7VbVeVTcDD+P8wfX3gKruVtX9OMHliEPk2SWqug1YCnzVTToNqFbVBZ04dw/wAfA1N2kaUKKqS7paDhMZLFiYXktVVwKvAbe0eqsJiG2VFgs0+J37EVACXCgiI4Fj+fIv7vdV9TTgVKAa+HUbZdgI3Aj8CigWkVkiMhAYCgwUkTLfA/g50L9VFnv9nlcDKYfIszueAa5wn3+9jc94KI8DV7rPrwSe7GYZTASwYGF6uzuAa3Caeny2A8NaHTcc2NYq7QmcGsWVwNuqWtTWBdTZ9OX/Aee38/4zqnoSToBQ4F5gB7BFVTP8Hqmqek5nPlQ7ebZ5aAdZ/QuYIiKDcWoYhwoWrfN6GZgkIhNw+nOe7qjcJnJZsDC9mvsr/DnAv/P4OeAXIjJYRKLcTtzpwAutTn8COAMn2DzOIahqPRDXOl1ERovIaSISD9QCNUAz8AlQ6XYUJ4pItIhMEJFjOvpMh8izLUXAiEOUex9O384/cYLXmkNc+gt5qWotzj17BvhEVbd3VHYTuSxYmEhwJ5Dc6vXHwEfAAeB3wDfcZqsWqrrVPS4ZePVQF3D/cNe18VY8cA9Ok9ZeIBe4VVWbcH6NHwFscd//B5Deic/TZp7tHHs3TmAsE5H2OvufwQmKHTVBtZXX48BErAmqzxPbVtWYLxORYcAfVPUS9/WvgARVbd0/EtFEJB9YCwxQ1Qqvy2O8E9YzVY3x2Kki8h7OBMCFOKOv+gx3+O9PgFkWKIzVLIwxX+IupVKEMyhgmqru8LhIxmMWLIwxxnQo7Juh3F84fwXqgbmqasP3jDEmxDypWYjIozgjRYpVdYJf+jTgfpw24n+o6j0i8k2gTFVni8hzqnpZR/n369dPhw0bFqTSG2NMZFqyZEmJqra5IoJXNYvHgL/gt86Mu27Ng8CZwE5gkYi8CgwGVriHNXUm82HDhrF48eJAltcYYyKeiLSeuNrCk3kWqvoBsL9V8rHARlXd7E6AmgVcgBM4BrvHtFteEZkhIotFZPG+ffuCUWxjjOmzwmlS3iCcJRJ8drppLwIXi8hDOAuttUlVZ6pqoaoW5uR0Z105Y4wx7Qn7Dm5VrQK+25ljRWQ6ML2goKDDY40xxnReONUsdgFD/F4PdtM6TVVnq+qM9PTOrKhgjDGms8IpWCwCDhOR4SISh7Pu/yHX62lNRKaLyMzy8vKgFNAYY/oqT4KFiDwLzAdGi8hOEbnK3ff4euBtYA3wvKqu6kq+VrMwxpjg8KTPQlWvaCf9DeCN7uZrfRbGGBMc4dQM1WNe1yy2lFTx4Hsb+XRHmSfXN8aYYImoYOGlNXsqOOf+D/n92+u4+KGPeW9dsddFMsaYgImoYOFlB/ddr68mITaK//z4FEb1T+Wm55dTWdvQ8YnGGNMLRFSw8KoZak95DfM2lnL1ySMY1T+Vey+eSGlVPf/4cEtIy2GMMcESUcHCK3PXOcuLnDmuPwCTBmdwxthcnlywjdqGTi1nZYwxYS2igoVXzVDzNpaQl57AYbkpLWlXnTSC/VX1vLSsS/MKjTEmLEVUsPCqGWpj8UHG5qUhIi1px4/IYvzANB79aAu2wZQxpreLqGDhheZmZUtJFSNzkr+QLiJ878ThbCg+yAcbSjwqnTHGBEZEBQsvmqF2ldVQ19jMiJyUL703/fCB5KTG88hH1tFtjOndIipYeNEMtWnfQQBGthEs4mKi+M4Jw/hg/T4+2dJ6+w5jjOk9IipYeGHHgRoAhmYntfn+d08cRl56Ar+evYqGpuZQFs0YYwLGgkUP7auoJUqgX0p8m+8nxcVw+3njWLW7gjteXUVzs3V2G2N6HwsWPVRUUUd2SjzRUdLuMWdPzOPaU0fyzMLt/GjWMuoabe6FMaZ3Cfud8rrCi1VniytryU1tu1bh7+Zpo8lIiuWeN9dSU9/Ew98qJOoQAcYYY8JJRNUsvOjgLq6s61SwEBGuPXUkt583jjlri3lh6c4QlM4YYwIjooKFF5xgkdDp479zwjAmDkrnEVs3yhjTi1iw6IGmZqX0YB25aR3XLHyiooRLCwezrqiSNXsqglg6Y4wJHAsWPVBaVUezQk4nmqH8nT0xD4D/rrU9L4wxvYMFix4or3b2q8hMiuvSef1S4hmRk8yy7QeCUSxjjAk4CxY9UFbjBIuMpNgun3t0fiZLth2wRQaNMb1CRAWLUK8NVebWLDISu1azADgyP5MD1Q3s2F8T6GIZY0zARVSwCPXQ2bLqeqB7NYvRA1IBWF9UGdAyGWNMMERUsAi1crcZKr0bwaLA3ShpfbEFC2NM+LNg0QNl1Q1ERwmp8V2fCJ+eGMuAtAQ2FB0MQsmMMSawLFj0QFlNPWkJMV/YIa8rDuufwgarWRhjegELFj1QVt1ARheHzfob3i+ZbSXVNiLKGBP2LFj0QHlNA+mJXe+v8BmanUxlXSMH3FFVxhgTrixY9EB5TUO3RkL5DM1yNkzaVloVqCIZY0xQhH2wEJERIvKIiLzgdVlaK6tuIKNHNQtfsKgOVJGMMSYoghosRORRESkWkZWt0qeJyDoR2SgitxwqD1XdrKpXBbOc3VVWXd+jPoshWUmIWLAwxoS/YG9+9BjwF+AJX4KIRAMPAmcCO4FFIvIqEA3c3er876lqWK6219SsVNQ29qjPIiE2mgFpCWzbb81QxpjwFtRgoaofiMiwVsnHAhtVdTOAiMwCLlDVu4HzunstEZkBzADIz8/vbjadVtGDdaH8Dc1OspqFMSbsedFnMQjY4fd6p5vWJhHJFpG/AUeKyK3tHaeqM1W1UFULc3JyAlfadvRkEUF/Q7OS2b7fgoUxJryF/R7cqloKXNuZY0O5B3fLulDdWETQ3+DMRPZV1lHb0ERCbHQgimaMMQHnRc1iFzDE7/VgN63HQrmQYFkP1oXyNygzEYDdZbb6rDEmfHkRLBYBh4nIcBGJAy4HXg1ExqFcory8ZXnyHgaLDCdY7LJgYYwJY8EeOvssMB8YLSI7ReQqVW0ErgfeBtYAz6vqqkBcL6Q1i5blyXvYDOVOzNt5wIKFMSZ8BXs01BXtpL8BvBHo64W0z8JthkpL6Nkt7J8aT3SUsMuChTEmjIX9DO6uCGXNorymgdT4GGKie3YLY6KjGJCWwM4DNiLKGBO+IipYhFJ5TQNpPeyv8BmcmWh9FsaYsBZRwSKUHdwVNY0BCxaDMhOtGcoYE9YiKliEshmqoqaB9MTAdPkMzkhkb0UtDU3NAcnPGGMCLaKCRSiV1zSQlhCoZqgkmhX2ltcGJD9jjAm0iAoWIW2Gqu3Zxkf+fBPzdlgntzEmTEVUsAj1aKiABQvfxDzrtzDGhKmIChah0tDUTHV9U8A6uPMyEhCxiXnGmPAVUcEiVM1QvuXJA1WziI+JJjc13obPGmPCVkQFi1A1Q5X7Zm8HaDQUOE1RNjHPGBOuIipYhEp5gGsW4IyIspqFMSZcWbDohoraRiCwwWJQZiJ7ymppataA5WmMMYFiwaIbWpqhAjTPApxmqMZmpajC5loYY8JPRAWL3trBDZ/PtbCmKGNMOIqoYBH6Du4A9lnYXAtjTBiLqGARKhU1DcTFRAV0z2yrWRhjwpkFi24I5FIfPklxMWQmxdrEPGNMWLJg0Q2BXOrD3yDb18IYE6YsWHSDs+Js4HekHZSRyC6bmGeMCUMWLLqhoqYxODWLDGdinqrNtTDGhJeIChahGjobzGao2oZm9lfVBzxvY4zpiYgKFqEcOhvIYbM+LUuVW7+FMSbMRFSwCIXmZqUyCKOhAAZn2lwLY0x4smDRRQfrG2nWwM7e9rGahTEmXFmw6KLy6sCvC+WTkRRLUly0BQtjTNixYNFFFbWBX+rDR0Tc4bMWLIwx4cWCRRcFY+MjfzYxzxgTjixYdFEwVpz1NyjDgoUxJvz0imAhIheKyMMi8pyInOVlWSpqAr/xkb9BmYmUVTdQVdcYlPyNMaY7gh4sRORRESkWkZWt0qeJyDoR2SgitxwqD1V9WVWvAa4FLgtmeTsSjOXJ/dmIKGNMOApFzeIxYJp/gohEAw8CZwPjgCtEZJyITBSR11o9cv1O/YV7nmcqahuIEkiJC06fhc21MMaEo+D8xfOjqh+IyLBWyccCG1V1M4CIzAIuUNW7gfNa5yEiAtwDvKmqS9u6jojMAGYA5OfnB+4DtOKbvR0VJUHJf1BGEgA7rWZhjAkjXvVZDAJ2+L3e6aa154fAGcAlInJtWweo6kxVLVTVwpycnMCVtBVnxdngNEEB5KbGExstVrMwxoSVoNcsAkFVHwAe6Og4EZkOTC8oKAhaWcqqG8hICl6wiIoS8tJtRJQxJrx4VbPYBQzxez3YTeuRUCwkWFZdT0ZSXNDyB9vXwhgTfrwKFouAw0RkuIjEAZcDr/Y001AsUX6guoHMINYswCbmGWPCTyiGzj4LzAdGi8hOEblKVRuB64G3gTXA86q6qqfXCkXN4kB1PZkhqFkUV9ZR39gc1OsYY0xnhWI01BXtpL8BvBHIawW7z6KxqZnK2sag9lmAU7NQhT3lNQzNTg7qtYwxpjN6xQzuzgp2zcI3IS/YNYvBGTbXwhgTXiIqWATbAXd58lDULMDmWhhjwkdEBYtgd3CXVTt7Ywe7ZpGXnoiI1SyMMeEjooJFsJuhfDWLYAeLuJgoclPjbUSUMSZsRFSwCLYDbs0i2M1QgG2CZIwJKxEVLELWDJUc3JoFwKDMJKtZGGPCRkQFi1A0Q8VGC8lx0UHJ39/gzET2lNfQ1KxBv5YxxnSkU8FCREaKSLz7fIqI/EhEMoJasjBUVl1PemIcziK4wZWflURDk7LbahfGmDDQ2ZrFv4EmESkAZuKs6/RM0ErVTcFuhjpQFfylPnyGZjlLlW/fb2tEGWO819lg0ewu0fFV4M+q+r9AXvCK1T3Bb4YK/lIfPvnZTrDYVmrBwhjjvc4GiwYRuQL4NvCamxaan9hhJNjLk/vLS08kNlrYtr8qJNczxphD6Wyw+C4wGfiNqm4RkeHAk8ErVngKZc0iOkoYkpnEDmuGMsaEgU4tJKiqq4EfAYhIJpCqqvcGs2DhRlWdmkVy6CpU+dlJ1gxljAkLnR0NNVdE0kQkC1gKPCwifwxu0boumB3c1fVN1Dc1h6xmAU4n9/bSalRt+KwxxludbYZKV9UK4CLgCVU9DmdP7LASzA7uAy3rQoWyZpFMZV1jyzIjxhjjlc4GixgRyQMu5fMO7j6l9KATLPqlxIfsmr7hs9tKrZPbGOOtzgaLO3F2tdukqotEZASwIXjFCj8lB+uAEAeLbJtrYYwJD53t4P4X8C+/15uBi4NVqHDkq1lkp4Suz2JIls21MMaEh852cI8SkTkistJ9PUlEfhHcooWXkqrQ1ywSYqMZkJZgwcIY47nONkM9DNwKNACo6mfA5cEqVDgqqawnJT6GhNjgLyLoLz87ie02Mc8Y47HOBoskVf2kVVpjoAvTU8EcOltaVRfSJiifoVk218IY473OBosSERkJKICIXALsCVqpuimYQ2dLDtaRHYJ9LFobmp1EcWUdNfVNIb+2Mcb4dDZY/AD4OzBGRHYBNwLfD1ahwlHpwfqQ9lf4DOuXDMCWEmuKMsZ4p1PBQlU3q+oZQA4wRlVPUtWtQS1ZmCk5WE+2B8FiZE4KABv3HQz5tY0xxueQQ2dF5CftpAOgqmG35EcwNDUr+6vq6OdBn8XwfslECWwstmBhjPFOR/MsUg/xXp9ZsKisup5mDe2wWZ+E2GiGZCWxyWoWxhgPHTJYqOqvAUTkRFWd5/+eiJwYzIKFk9Kq0E/I81eQk8Imq1kYYzzU2Q7uP3cyLSKVVDoT8rKTQ1+zACjITWFzSRVNzX2mMmeMCTMd9VlMBk4Aclr1X6QBIZmdJiJjgRuAfsAcVX0oFNf1V+LWLHJSvalZjMxNob6xmR37q1tGRxljTCh1VLOIA1Jwgkqq36MCuKSjzEXkUREp9i0T4pc+TUTWichGEbnlUHmo6hpVvRZnxVtPmr68rlm0jIiypihjjEc66rN4H3hfRB5T1W3dyP8x4C/AE74EEYkGHgTOBHYCi0TkVZyayt2tzv+eqhaLyPk48zo82cq1tKqO6CghPdGbbccLcj8fPnsG/T0pgzGmb+uoGeo+Vb0R+IuIfKnBXFXPP9T5qvqBiAxrlXwssNFduRYRmQVcoKp3A+e1k8+rwKsi8jrwTDtlnQHMAMjPzz9UsbqsqKKO3NR4oqIkoPl2VnpiLP3T4llfVOnJ9Y0xpqOhs75f8n8I4DUHATv8Xu8EjmvvYBGZgrNDXzzwRnvHqepMYCZAYWFhQHuCiypqyU1LCGSWXTZ6QBrr9lqwMMZ4o6NmqCXuf98PTXHaLMNcYG5njhWR6cD0goKCgJahqKKW4R53LI8dkMo/55XS2NRMTHRnB7EZY0xgdHY/ixNF5B0RWS8im0Vki4hs7uY1dwFD/F4PdtN6LFgLCe4tr6W/xzWLMXmp1Dc12xpRxhhPdGqnPOAR4MfAEqCny58uAg4TkeE4QeJy4Os9zBMITs2ipr6JitpG74PFgDQA1u6t5LD+h5pYb4wxgdfZ9oxyVX1TVYtVtdT36OgkEXkWmA+MFpGdInKVqjYC1+Ps6b0GeF5VV3X7E/gJRs2iqKIWwPNgMTInhZgoYe3eCk/LYYzpmzpbs3hPRH4PvAjU+RJVdemhTlLVK9pJf4NDdFZ3VzBqFr5gMcDjYBEXE8XInBTW7rFObmNM6HU2WPhGKx3t/ldwFhI8LeAl6gFVnQ3MLiwsvCZQee71BYt0bybk+RuTl8rirQe8LoYxpg/q7BLlr7n/VWAf8JGqbglmwcJFcYVTkfJ66CzA6AGpvPLpbipqG0hL8GaCoDGmb+qoz8K3vEeK+0gFCoE3ReTyIJety4KxB/feilqS4qJJje9sJSx4xrqd3DbfwhgTap1aorw1EckC3gVmBaNQ3RWsZqj+aQktGz55adxAJ1is3FXOMcOyPC6NMaYv6dbsLlXdj9NvEfGKK2rpn+Z9fwU4I7JyU+P5bGfgak7GGNMZ3QoWIjIVCLue1mA1Q3k9bNbfpMEZfLazzOtiGGP6mI46uFfw5e1Ts4DdwLeCVajuCnQzlKpSVFHn+bBZf4cPTmfO2iIqaxtItU5uY0yIdNRr23oVWAVKVbVPrDlRVt1AfWNzWIyE8pk4OB1VWLGrnBNG9vO6OMaYPqKjDu7u7GHhmUBPyiuqDI8Jef4mDc4A4LOdFiyMMaETUcuXBnq5j73l4TMhzycrOY4hWYnWb2GMCamIChaB5lvqIzc1fGoWAEflZ7Jo6wFUA7pthzHGtMuCxSEUtczeDp+aBcDkEdnsq6xj074+0XVkjAkDFiwOYW9FLVnJccTHRHtdlC+YPDIbgPmbSjwuiTGmr4ioYBHoeRbFYTbHwic/K4lBGYnM39zhKvHGGBMQERUsAt7BHUazt/2JCMePyGb+plKam63fwhgTfBEVLAIt3Cbk+TthZDYHqhtYV2SLChpjgs+CRTsampopOVgXVhPy/Pn6LeZttH4LY0zwWbBox77KOlTDa0Kev4EZiRTkpvCf1UVeF8UY0wdYsGhHURjtkNee8w8fyKKt+9ldVuN1UYwxEc6CRTvCdUKev/MPH4gqvPbZbq+LYoyJcBEVLAI5dNY3IW9AevgGi2H9kjkqP4OnF26nyUZFGWOCKKKCRSCHzu6tqCU2WshKigtAyYLnmpNHsK20mrdX7fW6KMaYCBZRwSKQiipqyU1NICoqvDcEPGv8AIZlJ/H39zfZWlHGmKCxYNGOoorasFsTqi3RUcI1p4xg+c5y5m+yGd3GmOCwYNGOveW1YTtstrWLjxrMwPQE7nxtNY1NzV4XxxgTgSxYtKO4oi4s14VqS0JsNLdPH8/avZX8+b8bvS6OMSYCdbStap9UVddIZV1jrwkWAF8Z35+LjhrE/XM2ECXC96eMZHPJQVbuqqC8poGE2Chio6Ooa2hiSFYS4wems6e8hoTYaEb1T/W6+MaYMGfBog29YUJeayLCPRdNAoU/vbueP727vtPn3nfZEVx45KAgls4Y09v1imAhIsnA+8CvVPW1YF9vrxss+ofxhLy2xMVE8X+XHs65k/JYvqOMgRmJHD00k/SkWIor6mhsVgakJbC1tIol2w4wKCORWYu2c9O/llNysI6zJ+YxKCPR649hjAlDQQ0WIvIocB5QrKoT/NKnAfcD0cA/VPWeDrK6GXg+aAVtpdidkNc/jCfktUdEOH1sf04f2/8L6f4z0QekJ3D8CGchwtPH5nL144u56/U13PX6Go4YksG9F09i9ABrmjLGfC7YHdyPAdP8E0QkGngQOBsYB1whIuNEZKKIvNbqkSsiZwKrgeIgl7VFS82iF/VZdFdqQiyzZhzPuz85hZ+fM4adB6q58MF5vPLpLpu3YYxpEdSahap+ICLDWiUfC2xU1c0AIjILuEBV78aphXyBiEwBknECS42IvKGqXxofKiIzgBkA+fn5PSr33vJaUuJjSInvFa10PSYiFOSmUpCbyoVHDOK6p5dyw6xPWbhlP3eeP56YaBs0Z0xf58VfgUHADr/XO920Nqnqbap6I/AM8HBbgcI9bqaqFqpqYU5OTo8KWFzZOybkBUNuWgLPXHM83zlhGM8s3M61Ty2lvtHmbhjT1/Wan4yq+lhHnduBWkiwN03IC4a4mCjumD6OW84ew7trivjBM0ttsp8xfZwXwWIXMMTv9WA3rccCtZBgOG+nGioiwrWnjuTmaWN4Z3URNz73qfVhGNOHeREsFgGHichwEYkDLgdeDUTGgahZNDer2wzVt4OFz/enjOSnZ47itc/28NyiHR2fYIyJSEENFiLyLDAfGC0iO0XkKlVtBK4H3gbWAM+r6qpAXC8QNYsD1fU0NCkD+mifRVuum1rASQX9uP2VVSzfUeZ1cYwxHghqsFDVK1Q1T1VjVXWwqj7ipr+hqqNUdaSq/iZQ1wtEzaIvDZvtrOgo4c9XHElOajzXPrWEkoN1XhfJGBNivaaDuzMCUbPwLfXRGyfkBVNmchx//+bR7K+q587Zq70ujjEmxCIqWATC3nJ3O1WrWXzJhEHpfOfEYcz+bDfriyq9Lo4xJoQiKlgEqhlKBHJSrc+iLf9zykjSEmL5xUsrabZ9v43pMyIqWASkGaq8ln4p8cTarOU2ZSXHcdu5Y/lk635m2egoY/oM+4vYyt6Kvj0hrzO+dvRgjh+RxT1vrmFfpXV2m76hpr7J6yJ4KqIWPxKR6cD0goKCbuext7yWIVlJgStUBBIRfvPViZx934f87IXlPPytQls/ykQkVeWzneW89tluHv94G/1S4hg3MJ342CjOnZjH0OwkBmUkkpEU53VRgy6igoWqzgZmFxYWXtPdPPZW1HLM8MwAlioyjcxJ4dZzxvDr2at5euF2vn3CMK+LZEzAzFlTxL+X7mTexlLKaxoA6JcSx2H9U5m3sYSahiZe/2xPy/GHD07nuBHZnHxYPyaPyI7IH08RFSx6qrahifKaBmuG6qTvnDCMN1fs5c7XVnP00EwmDOrZMivGeGl/VT0vL9vFn/+7gQPVToA4sSCbKaNyOX1sLiNyUgCntrF0exmVtQ3MWVPMS8t2sXxnOct3ljPzg80kx0Vz+tj+zDhlRET9m5BIWu/Hrxnqmg0bNnT5/K0lVUz5w1x+f8kkvlY4pOMTDHvKa5j6h7mMzUvjie8dS2pCrNdFMqZLXvl0F0/O38bibQcAGJuXxnmT8vj6sflkJne+eamoopb/rC7indVFfLRhHwmx0Xzz+KFMP3wg4wemISLB+ggBIyJLVLWwzfciKVj4FBYW6uLFi7t83oLNpVw+cwFPXXUcJx3WLwgli0yzl+/mh88u46xx/bnrwgm2rpYJe/ur6tm07yD3vbueeRtLSU+M5ewJAzh7Yh4njMzu8WjIoopaZjyxmOU7nWH8507M41fnjw/7IfmHChbWDOXHN3t7QHp4/w8NN9MPH8i20irun7OBc//8Ec/NOL6lym5MONlVVsMdr6zi3TVFLWnfPXEYt549lriYwPUz9E9L4JXrT2LzvoM8/OFmnv1kB6+v2MMZY/tz78UTyU7pfX9jLFj42Vtu60J11/WnHcaZ4wZwyUMfc/eba/n7lUcTFRX+1W7TN+wuq2HmB5t5Yv5WmhUuOGIgJxb049RROUH99z4iJ4W7L3KatZ9asI0Xl+7inAfK+MHUAi4tHEJCbHTQrh1oFiz87K2oJTku2trdu2n0gFSunTKS37+9jpv//Rm/vWiiTW40nlqybT9/eHs9C7aUogpH5mfw6/PHM2lwRkjLcVR+JkflZ/LVIwdxw6xPuf2VVXy4oaRX/aiKqGDR03kWMVHC6AGpgS1UH3PdlJHUNTbzwJwN7C6vYeY3C0nuI3uZm/CxpaSK215awcebSgEoyE3hvsuO8Lyj+eTDcvj4ltP42/ubuO/dDVz9xGJuO3csI3tBs611cJugeH7xDn72wmdce+pIbjl7jNfFMX1EeU0D/5y3hUc+3EJlXSPJcdHcdu44Lj9mSFj9gldV7nh1FU/M30Z+VhJv3nByWPyosg5uE3KXFg5hweZS/vb+JuJjovjxmaO8LpKJUOXVDdw3Zz3/nLe1JS0vPYHvnDiM66YUkBgXfv0CIsKdF0zgvEkDuWzmfC58cB73X34k4wameV20dlmwMEHzu4snESXC/XM2cMqofhw9NMvrIpkIU1xRy/XPLuOTLfsBONedH3FiQe8Y+n7s8CxmnDKCv7+/mXMe+JDlt59FelJ49pla76MJmpjoKH51/nj6pcTx/aeWsrH4oNdFMhFkV1kNx/52Dp9s2c8tZ49h02/P4cGvH9VrAoXPT84cxRXH5gNw878/o7q+0eMStc2ChQmqlPgYnr76eJpVuXzmAjYW26ZJpufeX7+PM//4PgAnFfRjxskjiA6jPomuiI+J5u6LJnL1ScN5a9Vezn3gI2obwm+FWwsWJuhGD0hl1ozJqCq/eNk2TTLdV1PfxPQ/f8S3H/2E6vom/vqNo3jq6uPCqvO6u249Zyw3nTWKLSVV/PGd9V4X50siKlgEYqc8ExwFuSlcN7WABZv38+PnP6WhqdnrIplepK6xiY82lDD29rdYscv59/3Jz0/nnIl5HpcscKKjhOtPO4xTRuUw84PN3DhrGRW1DV4Xq0VEBYtA7JRngud7Jw7jf78ymlc+3c2vXl3ldXFML1Fd38iJ97zHlY8sbEl784aTI3YNsoe+cRRDshJ5+dPd/G3uJq+L0yKigoUJbyLCD6YWcOXx+Ty/eAfz3QlTxrSnqVn53mOLKDno7Mj4h68dzopfncXYvPAdYtpTyfExfPC/U/nK+P488tGWls/uNQsWJuRuOH0UQ7OT+fY/P2GO34Juxvj7y383MPLnb7Bg837OGJvLlrvP4ZKjB/eJ5Xh8P6zqGpu58MF5LYucesmChQm5nNR4nv+fyYzun8oNsz6lMozaZU14eGHJTv7wH6eTd+roHP525dG9Yj+IQJo0OIOvH5fPzgM1HPfbOZ4HDAsWxhNZyXHcdeEEquobufKRT9hXGR5VbeOtxz/eylf/Oo87XlkJOMuH//O7x0bkNqWd4b9UztcfXuDpwJC++X/AhIXDh2TwtyuPZt3eCi5+6OOwHFtuQuehuZu449VVlFU3MGVMLi9ddwJ3TB/vdbE8lZYQy5o7p3HRkYPYtK+K37y+xrOy2EKCxnPvrS3mu48tYsYpI7j17DF9rrmhL1q4uZQdB2oYlp3E3HX7KK6s5eVlu4mPjeLf3z+BUf1t9Wd/zc3K4Xf+h8raRmZffxITBwdnxGevXkhQRKYA/w9YBcxS1blelscE3pTROVx+zBBmfrCZ6vpG7rpwotdFMkFS29DEhxtKuOaJL/6YG5CWwJi8VP546REU5Ib/ct2hFhUlzJpxPOc+8BE3zFrGnJ+eGvIfVUENFiLyKHAeUKyqE/zSpwH3A9HAP1T1nkNko8BBIAHYGcTiGo+ICHdfNJGYaOHphds5Z0IeJ/Sy9X1M2x79aAuTBqczbmAa59z/IVtLqwFIjI3mR6cfxr1vreW44Vk89z+TPS5p+Bs/MJ0rj8/nqQXbeWrhdr55/NCQXj+ozVAicgrOH/onfMFCRKKB9cCZOH/8FwFX4ASOu1tl8T2gRFWbRaQ/8EdV/UZH17VmqN6p5GAdl89cwJaSKv502RGcf/hAr4tkumFveS2b9x1k+/5qbnlxBQDHDc9i4Zb9fOeEYYwbmMb4gWmMy0uj5GA9mUmxfbYDu6uq6xs54tfvUN/UzLJfnklmclxA8/esGUpVPxCRYa2SjwU2qupmt3CzgAtU9W6cWkh7DgDt7nIuIjOAGQD5+fk9KbbxSL+UeF667gS+/vBC7n5jDRMHpTO8X7LXxTJdsK20iql/mEvr5b8qahv56pGDuO3csV/Yajcntd1/0qYNSXEx/OK8sdz+yip+9/Za7r5oUsiuHfQObjdYvOZXs7gEmKaqV7uvvwkcp6rXt3P+RcBXgAzgoc70WVjNondbvHU/Vz2+mIamZn53ySTOm2Q1jHBUVl3Pqt0VbN9f3fJYs6eCzfuq+NuVRzNmQCoD0hNIiA2/zYd6s+ZmZeztb1HX2Mz6u84mLiZwtbJe3cGtqi8CL3bm2J7uwW3CQ+GwLN668WSufXIJt720kiPzMxmUkeh1sfqkusYmNhVXkRIfQ352Uku6qnLhg/Na+iBiooTBmYkMyUriq0cMYtqEAV4VOeJFRQn3XXYE3396KXe9vpo7L5jQ8UkB4EWw2AUM8Xs92E3rMVWdDcwuLCy8JhD5Ge/kpSdy14UTuXzmfKbd9wG/v2QS0yZEzgqj4aqusYnH5m1l5e4K1u6pYEtJFY3NSr+UOB777rFsKali076DbCg6yNbSan4wdSRXHJtPXnpir91Pojc6c1x/MpNieWL+Ni4tHMKEQcFfPNWLZqgYnA7u03GCxCLg66ra42VI/WoW12zYsKGn2ZkwsK20ih8+u4y1eyp57Ucn2fj7IJu9fDc/fHYZgzISGZuXyugBqdTUN/PovC0tx4jA4MxERvdP5Z6LJ9EvxfodvPDkgm388uWVnHxYP5686riA5HmoZqhgj4Z6FpgC9AOKgDtU9REROQe4D2cE1KOq+ptAXtf6LCLLxuJKpv95HjUNTVx05CB+/7XD7VdsAKgqW0qqWLa9jBW7ylm5q5xVuyuIiRaW/fLMlhFK9Y3NzFq0nazkOEbmpDC8X7L1Q4SJW19cwbOfbOfdn5xCQW7Pf0h5FixCzWoWkWtfZR0PvreRxz7eyg+mjuS6KQUkx4d9l1tYKa9pYPmOMpZtL2PZjgN8uqOMsmpnEcekuGjGD0xj/MB0zhrX3+a59BLzN5VyxcML+OmZo/jh6Yf1OL8+Eyx8rGYRmRqamrnu6aW8s7qIjKRY/nTpEUwdk+t1scJWycE6Fm7ez4LNpSzcUsr6ooOA04w0KjeVI/Mz3EcmI3NSrLbWS53yu/fYvr86ICOjevVoKGN8YqOjePhbhSzbfoCf/ms5N//7M34wtYDzJuWRbe3mlB6sY+EWJzgs2Px5cEiOi6ZwWBbnHz6QI/MzmTQ4vU/sCdFXHDc8ywkWRZVB7eiOqJqFNUP1HUu2HeAXL69kzZ4KYqKEv37jKM4a37eGazY0NbNk2wHeX7+P99ftY/WeCsBpUioclsXkEdkcPyKLCYPSvzARzkSWdXsr+cp9H3DFsfncfVHP1lWzZigTsdbureC6p5ZSVFHLlNG5nDW+P+cfPjAiVq6trm9k/qZS3ltXzNJtZdz11Qnkpsa3BIePN5VysK6RmCjh6KGZnDIqh8kjs5lowaFPqW9sZtQv3iQ/K4kPfja1R3lZM5SJWGMGpPH3bx7NIx9tYc7aYl5fsYcNRQe54rh8BqYn9LqgsaWkivfWFvPeumIWbtlPfWMzSXHRNDYpV/5jIdX1zp4fgzISOf+IgZw6KocTRmZbs1IfFhcTxfVTC/jr3I00NDUH7YdCRNUsrBmqb2tuVi6fuYBPtu4HID0xltvOHculhUM6ODO4Sg7WsX5vJct3lvPZzjK2llZz01mjmDo6l11lNawvquTDDSXMXVfcMiN6ZE4yU0bnMnV0LscMz2Tm+5tZvO0Ap4zK4dRROYzMSe51gdAEzwtLdnLTv5YzcVA6s394UrfzsWYo02fUNzazek8FK3aWce9b6zhY10i/lDgGZiS2NNfcdu44wAku9e4vse6OBGpqVtbsqWD+JqdT+bD+qVx98nAWbC7l402lfLyxpCUAAORnJbF9fzUp8TE0NjdT2+BskxkfE8UJI7OZOiaXKaNyv7C0hjEdOVjXyJTfz6Wsup71d51NVDe/zxYsTJ9UVFHL84t2sLu8hp0HathdVsOmfc46R/WNzdS7+xkfPiSDm6eNZm95LXvKa9m8r4rymgbuvmgiFbUNrN5dwYpd5azYWc620iruvGAC+6vqWbC5lLnr97G/qh6A1IQYKmsbW66fGh/DcSOyOG54NmPyUpk4KJ2MpDheXraLV5fvZkS/ZEbmpjCiXzKHD8mwiW6mR56cv5VfvrKKT247ndzUhG7lYcHCGOBAVT1/eW8jzarExUSREBPNjv3VvLjsi0uT5abGU1xZR5TQstR2XHQUYwemsW5vRUttIDs5jhML+nHamFwmj8ymqVn57RtrGDcwjRNG9mPCwDTbp8GEzJw1RVz1+GJeuu4EjszP7FYefSZYWJ+F6arGpmbeXVNMWmIMeemJDEhLIDEumr/O3UhJZT1j81IZm5fGqP6pxMVE8c7qIvaU1zB5RDYFuSnWb2DCxs4D1TyzcDuXHTOEodnd2wemzwQLH6tZGGNM1x0qWFgd2RhjTIcsWBhjjOlQRAULEZkuIjPLy8u9LooxxkSUiAoWqjpbVWekpwd/1yhjjOlLIipYGGOMCQ4LFsYYYzpkwcIYY0yHLFgYY4zpUEROyhORfcA2j4vRDyjxuAzhzO5P++zeHJrdn/b19N4MVdWctt6IyGARDkRkcXszIY3dn0Oxe3Nodn/aF8x7Y81QxhhjOmTBwhhjTIcsWATPTK8LEObs/rTP7s2h2f1pX9DujfVZGGOM6ZDVLIwxxnTIgoUxxpgOWbDoJBEZIiLvichqEVklIje46Vki8o6IbHD/m+mmi4g8ICIbReQzETnKL69vu8dvEJFve/WZgkFEokVkmYi85r4eLiIL3fvwnIjEuenx7uuN7vvD/PK41U1fJyJf8eijBJSIZIjICyKyVkTWiMhk++58TkR+7P67Wikiz4pIQl/+7ojIoyJSLCIr/dIC9n0RkaNFZIV7zgPSmS0fVdUenXgAecBR7vNUYD0wDvgdcIubfgtwr/v8HOBNQIDjgYVuehaw2f1vpvs80+vPF8D79BPgGeA19/XzwOXu878B33efXwf8zX1+OfCc+3wcsByIB4YDm4Borz9XAO7L48DV7vM4IMO+Oy33ZhCwBUj0+858py9/d4BTgKOAlX5pAfu+AJ+4x4p77tkdlsnrm9JbH8ArwJnAOiDPTcsD1rnP/w5c4Xf8Ovf9K4C/+6V/4bje/AAGA3OA04DX3C9iCRDjvj8ZeNt9/jYw2X0e4x4nwK3ArX55thzXWx9AuvvHUFql23dHW4LFDvePWoz73flKX//uAMNaBYuAfF/c99b6pX/huPYe1gzVDW6190hgIdBfVfe4b+0F+rvPff8AfHa6ae2lR4L7gJ8Bze7rbKBMVRvd1/6fteU+uO+Xu8dH4v0ZDuwD/uk20f1DRJKx7w4AqroL+AOwHdiD811Ygn13WgvU92WQ+7x1+iFZsOgiEUkB/g3cqKoV/u+pE6b75FhkETkPKFbVJV6XJQzF4DQpPKSqRwJVOM0ILfr4dycTuAAnqA4EkoFpnhYqzHnxfbFg0QUiEosTKJ5W1Rfd5CIRyXPfzwOK3fRdwBC/0we7ae2l93YnAueLyFZgFk5T1P1AhojEuMf4f9aW++C+nw6UEpn3ZyewU1UXuq9fwAke9t1xnAFsUdV9qtoAvIjzfbLvzhcF6vuyy33eOv2QLFh0kjta4BFgjar+0e+tVwHfKINv4/Rl+NK/5Y5UOB4od6uQbwNniUim+4vqLDetV1PVW1V1sKoOw+l0/K+qfgN4D7jEPaz1/fHdt0vc49VNv9wd8TIcOAynM67XUtW9wA4RGe0mnQ6sxr47PtuB40Ukyf135rs/ff6700pAvi/uexUicrx7v7/ll1f7vO7E6S0P4CScat9nwKfu4xycttI5wAbgXSDLPV6AB3FGZKwACv3y+h6w0X181+vPFoR7NYXPR0ONwPkHuxH4FxDvpie4rze674/wO/82976toxOjNHrDAzgCWOx+f17GGZ1i353PP9evgbXASuBJnBFNffa7AzyL03/TgFMzvSqQ3xeg0L3Xm4C/0GrwRVsPW+7DGGNMh6wZyhhjTIcsWBhjjOmQBQtjjDEdsmBhjDGmQxYsjDHGdMiChfGciPxJRG70e/22iPzD7/X/ichPAni9x0Tkko6P7HK+P/d7Psx/xdAOzrtRRL7l9zpGRPaJyD2BLqPf9ZKClPd5InJnMPI23rJgYcLBPOAEABGJAvoB4/3ePwH42INyddXPOz7ki9wZyN/DWanX50ycVY2/1qmlo7vuRqDNYCEi0T3M+3VgerCCkfGOBQsTDj7GWVUUnCCxEqh0Z57GA2OBpSJyu4gsEmfPg5nujNUxItIyS9f9Rb/CfX60iLwvIkvc2kpe6wu3d4yIzBWRe0XkExFZLyInu+lJIvK8OPuavCTOfgqFbi0gUUQ+FZGn3eyjReRhcfZp+I+IJLbx2U8DlurnC+aBswro/Tgzm333BRHZKiK/FpGl4uxFMMZNzxFnf4NV4ixSuE1E+olIsoi8LiLL3Xt2mYj8CGf9pfdE5D33/INu7W05MFlEfuIev9JX43Pv61q3VrZeRJ4WkTNEZJ44eyUcCy1rFs0FzuvU/3nTa1iwMJ5T1d1Ao4jk49Qi5uOs6DsZZ6bpClWtB/6iqseo6gQgEThPVdcCce7yDgCXAc+Js47Xn4FLVPVo4FHgN/7X7cQxMap6LM4v8TvctOuAA6o6DvglcLT7GW4BalT1CHWWOQFnuYkHVXU8UAZc3MbHPxFnhVVfmRJw1kqajTOL94pWx5eo6lHAQ8BNbtodOEtejMdZdyrfTZ8G7FbVw9179paqPgDsBqaq6lT3uGScPRAOB2qA7wLH4ex3cI2IHOkeVwD8HzDGfXwdZ2WDm/hirWoxcHIbn9X0YhYsTLj4GCdQ+ILFfL/X89xjprq/5Ffg/CL3NVU9jxMkcP/7HDAamAC8IyKfAr/gi4un0YljfItFLsHZWwCcP46zAFR1Jc7yHe3ZoqqftpGHvzyc5ct9zgPeU9UanEUrL2zVNNRRmd4CDrjpK4Az3RrSyapa3k45m9xr+fJ6SVWrVPWgez3fH/4tqrpCVZuBVcActyaxotVnK8apvZgIEtPxIcaEhK/fYiJOM9QO4KdABc4+EAnAX3HWvdkhIr/CWSMInODwLxF5EaclZIOITARWqepk2icdHFPn/reJ7v1bqfN73oRTG2qths8/Bzg1iZPEWb0XnPWATgPe6WqZVHW9OFtsngPcJSJzVLWtzudaVW06VF6trg3OniV1fs/9y5KA87lMBLGahQkXH+P8qt6vqk2quh9n69HJ7nu+P6gl4uwp0jKaSVU34fzx/CVO4ABnIbkcEZkMTpOTiPh3mnf2mNbmAZe6x4/DCW4+DW7TVleswWneQUTScH7F56vqMHVW8P0BX26KOlSZzsJZpBARGQhUq+pTwO9xlkUHqMTZGrgtH+LUZpLE2aDpq25aV4zCCfgmgliwMOFiBc4oqAWt0spVtURVy4CHcf4IvQ0sanX+c8CVOE1SuH0clwD3uh23n+KOuPLpzDFt+CtOgFkN3IXTHONr3pkJfObXwd0Zb+LstwzOH+b/qqr/L/hXcEYXxR8ij1/jLEW9Evgazi5qlTiB7BO3ie0Ot7y+cr7l6+D2p6pLgcdwVnNdCPxDVZd14fMATMUZFWUiiK06a0wXuP0HsapaKyIjcZaKHu0Gnu7m+RLwM1Xd0M3z44EmVW10a0kPqeoR3S1PT4hIf+AZVT3di+ub4LE+C2O6Jgln2GksTp/HdT0JFK5bcDq6uxUscEY/PS/OHJV64Joelqcn8nH6mkyEsZqFMcaYDlmfhTHGmA5ZsDDGGNMhCxbGGGM6ZMHCGGNMhyxYGGOM6dD/B96XcRuoTNREAAAAAElFTkSuQmCC\n" }, "metadata": { @@ -108,11 +108,18 @@ "execution_count": 4, "metadata": {}, "outputs": [ + { + "output_type": "stream", + "name": "stderr", + "text": [ + "/Users/lpsinger/Library/Caches/pypoetry/virtualenvs/dorado-sensitivity-RYVm8gWH-py3.8/lib/python3.8/site-packages/astropy/units/quantity.py:477: RuntimeWarning: divide by zero encountered in true_divide\n result = super().__array_ufunc__(function, method, *arrays, **kwargs)\n/Users/lpsinger/Library/Caches/pypoetry/virtualenvs/dorado-sensitivity-RYVm8gWH-py3.8/lib/python3.8/site-packages/astropy/units/quantity.py:477: RuntimeWarning: divide by zero encountered in true_divide\n result = super().__array_ufunc__(function, method, *arrays, **kwargs)\n" + ] + }, { "output_type": "execute_result", "data": { "text/plain": [ - "" + "" ] }, "metadata": {}, @@ -122,7 +129,7 @@ "output_type": "display_data", "data": { "text/plain": "
", - "image/svg+xml": "\n\n\n\n \n \n \n \n 2021-03-01T18:00:09.521334\n image/svg+xml\n \n \n Matplotlib v3.3.3, https://matplotlib.org/\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n\n", + "image/svg+xml": "\n\n\n\n \n \n \n \n 2021-03-18T19:44:40.431275\n image/svg+xml\n \n \n Matplotlib v3.3.3, https://matplotlib.org/\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n\n", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA7UElEQVR4nO3dd3xW9b3A8c83e5IdCAkZrIBMAVEB98LWVqtVi3tb63W1tbXXtmq9t9d729pha92jLY62KCoOtFWGikjYI4SZQEKAkIQMkkDG9/5xTiDEJztPnozv+/U6r+c56znfkwfyze/8lqgqxhhjTHN+vg7AGGNM72QJwhhjjEeWIIwxxnhkCcIYY4xHliCMMcZ4FODrALpTfHy8pqen+zoMY4zpM1auXHlAVRM87etXCSI9PZ2srCxfh2GMMX2GiOS1tM8eMRljjPHIEoQxxhiPLEEYY4zxqF/VQRhjTHvV1taSn59PTU2Nr0PpESEhIaSkpBAYGNjucyxBGGMGpPz8fCIjI0lPT0dEfB2OV6kqxcXF5Ofnk5GR0e7z7BGTMWZAqqmpIS4urt8nBwARIS4ursOlJUsQxpgBayAkh0adudcBnyDq6hv40yfbWLKlyNehGGNMrzLgE4Q/9QQs/h9yl73h61CMMQPYww8/zK9//Wtfh3GcAZ8gxD+Qq2QhQ/d94utQjDGmVxnwCQKgODSNmKpdvg7DGDPA/Pd//zejR49m1qxZ5OTkAPDss89y0kknMWnSJC677DKqqqqoqKggIyOD2tpaAMrLy49b9xZr5gocjhpO6qHFlFXVEhXW/jbCxpj+4ZF3NrJpT3m3fuYJQwfx0DfGtbh/5cqVvPbaa6xZs4a6ujqmTJnC1KlTufTSS7n11lsB+OlPf8rzzz/PXXfdxZlnnsm7777LJZdcwmuvvcall17aoT4NnWElCMAvIZMEKSO3oMDXoRhjBoilS5fyrW99i7CwMAYNGsQ3v/lNADZs2MBpp53GhAkTmDt3Lhs3bgTglltu4cUXXwTgxRdf5MYbb/R6jFaCAAaljIV1UJS7EUal+zocY0wPa+0v/Z52ww03MH/+fCZNmsRLL73EokWLAJg5cya5ubksWrSI+vp6xo8f7/VYrAQBxKU5P+jqws0+jsQYM1CcfvrpzJ8/n+rqaioqKnjnnXcAqKioICkpidraWubOnXvcOddddx1XXXVVj5QewBIEAAHxw6nDH7/irb4OxRgzQEyZMoUrr7ySSZMmceGFF3LSSScB8Oijj3LyySczc+ZMxowZc9w5V199NaWlpcyZM6dHYrRHTAD+gRwIHEpkZa6vIzHGDCAPPvggDz744Fe233HHHR6P//TTT/n2t79NdHS0lyNzWIJwVUZkMKR4O0fqGggKsIKVMaZ3ueuuu3j//fd57733euyaXv1NKCIviMh+EdnQZNtkEflCRNaISJaITG/h3OtFZKu7XO/NOAEa4kaRJvvIKyrz9qWMMabDnnjiCbZt28bo0aN77Jre/lP5JWB2s23/BzyiqpOBn7vrxxGRWOAh4GRgOvCQiMR4M9CwoWMIljoK83K8eRljjOkzvJogVHUJUNJ8MzDIfR8F7PFw6gXAR6paoqqlwEd8NdF0q8aWTOX52d68jDHG9Bm+qIO4F1goIr/GSVAzPByTDOxusp7vbvsKEbkNuA0gNTW100GFJjmtBer3b+n0ZxhjTH/ii9rYO4D7VHUYcB/wfFc+TFWfUdVpqjotISGh8x8UFku5XxQhZdu7Eo4xxvQbvkgQ1wONY2v/A6eOobkCYFiT9RR3m1eVhqYTV5OHqnr7UsYY49HPf/5z/vWvf7V6TEtDgx88eJAnn3yy22LxRYLYA5zhvj8b8NQ7bSFwvojEuJXT57vbvOpIzHDS2MO+8sPevpQxxnj0i1/8gnPPPbdT5/apBCEirwLLgEwRyReRm4Fbgd+IyFrgl7j1ByIyTUSeA1DVEuBRYIW7/MLd5lWBiZkkSDk7d9ugfcYY78rNzWXs2LHceuutjBs3jvPPP5/q6mpuuOEG/vnPfwLw3nvvMWbMGKZOncrdd9/NRRdddPT8TZs2ceaZZzJ8+HD+8Ic/APDAAw+wfft2Jk+ezP3339/lGL1aSa2qLfUHn+rh2CzglibrLwAveCk0j6KHnQCroGTXRhg/oicvbYzxpfcfgL3ru/czh0yACx9r9ZCtW7fy6quv8uyzz3LFFVcwb968o/tqamq4/fbbWbJkCRkZGV8ZXmPz5s188sknVFRUkJmZyR133MFjjz3Ghg0bWLNmTbfcgnUZbiJqmDOi4+G91tTVGON9GRkZTJ48GYCpU6eSm5t7dN/mzZsZPnw4GRkZAF9JEF//+tcJDg4mPj6exMRE9u3b1+3x2VAbTUhMOrUEEFBqLZmMGVDa+EvfW4KDg4++9/f3p7q6utPn1tXVdWtsYCWI4/kHUByUzKDKnb6OxBgzwGVmZrJjx46jpYrXX3+9zXMiIyOpqKjothgsQTRTNSiD5PoCKmq8O9erMca0JjQ0lCeffJLZs2czdepUIiMjiYqKavWcuLg4Zs6cyfjx47ulklr6U5v/adOmaVZWVpc+Y8dr95OS/TzZN21hUlp8N0VmjOltsrOzGTt2rK/DaFVlZSURERGoKnfeeSejRo3ivvvu6/TnebpnEVmpqtM8HW8liGYiho4hSOrZm2uD9hljfOvZZ59l8uTJjBs3jrKyMm6//fYevb5VUjcT6w7aV7knG5jp22CMMQPafffd16USQ1dZCaKZgERnrHU9YIP2GdPf9adH7G3pzL1agmguNIZy/2jCyq0lkzH9WUhICMXFxQMiSagqxcXFhISEdOg8e8TkQVlYOglledTWNxDobznUmP4oJSWF/Px8ioqKfB1KjwgJCSElJaVD51iC8KAuZiQZ5R+yq6SKEQkRvg7HGOMFgYGBR3spG8/sz2MPgodkEi/l5OXn+zoUY4zxGUsQHkSnOmMyHdy1yceRGGOM71iC8CDMnX70yD7rC2GMGbgsQXgSnUYtgQTaoH3GmAHMEoQn/gGUhiQTU71zQDSBM8YYTyxBtKB60AhSG/ZQVGnTjxpjBiZLEC3wSxhNmuxj+95SX4dijDE+YQmiBZEpYwmUevbtsiE3jDEDkyWIFkSlnABAdYFNP2qMGZi8liBE5AUR2S8iG5psmywiX4jIGhHJEpHpLZxb7x6zRkTe9laMrZH4Uc6b4q2+uLwxxvhcm0NtiEgizrjXQ4FqYAOQpaoNbZz6EvBH4C9Ntv0f8Iiqvi8iX3PXz/RwbrWqTm4rNq8KjabcP5aw8h2oKiLi03CMMaantViCEJGzRGQh8C5wIZAEnAD8FFgvIo+IyKCWzlfVJUBJ881A4zlRwJ4uxO51NVHDGVqfz5Z9lb4OxRhjelxrJYivAbeq6q7mO0QkALgIOA+Y14Hr3QssFJFf4ySnGS0cFyIiWUAd8Jiqzm/pA0XkNuA2gNTU1A6E0rbwYeMZW/war24uIHPImG79bGOM6e1aLEGo6v2ekoO7r05V56tqR5IDwB3Afao6DLgPeL6F49LcOVKvAn4nIiNaifMZVZ2mqtMSEhI6GE7rwjPPIUJq2LPx0279XGOM6QtaraQWkTNEZKL7/goR+aOI3CciwZ283vXAG+77fwAeK6lVtcB93QEsAk7s5PW6JuN0GvAjbu9nHDpc55MQjDHGV1qrg/gT8F/AcyLyN5y/5jcAU4AXOnm9PcAZ7vuzga80ERKRmMYEJCLxOBXkvhlWNTSayvhJzJR1fLGj2CchGGOMr7RWB3GWqp4gIiFAAZCoqvUi8jSwrq0PFpFXcVooxYtIPvAQcCvwe7cOowa37kBEpgHfVdVbgLHA0yLSgJPAHlNVn427HTb2PCYW/Yr/27Sdc8YO9lUYxhjT41pLEDUAqlojInmqWu+uq4jUtvXBqjqnhV1TPRybBdzivv8cmNDW5/eUgFHnwtL/o3rLJ7Rcp26MMf1PawkiUUS+D0iT97jr3Vsb3JslT+VIQARjDq0g98Ah0uPDfR2RMcb0iNYqqZ8FIoGIJu8b15/zfmi9hH8Adamncbr/ehbn7Pd1NMYY02NaLEGo6iMt7RORk7wTTu8UNvY8wna8z+ZNa2CmTXJujBkY2j0Wk4icICKPisg24M9ejKn3GXE2AGG7F3O4rt7HwRhjTM9odSwmEUkH5rhLLZAGTFPVXK9H1pvEZlAdPoxTyteSlVvKzJHxvo7IGGO8rrV+EMtwxmEKAC5T1alAxYBLDq6A0edwqt8mlm7u1cNHGWNMt2ntEdM+nErpwRxrtTRgJ2gOHH0uEVJDUbYNu2GMGRhaG4vpEpz+CCuBh0VkJxDT0hwO/V76aTSIP+llyyksq/Z1NMYY43WtVlKrapmqvqiq5wMnAz8Dfisiu3skut4kNJrDiSdymt96lmwp8nU0xhjjde1uxaSq+1X1j6o6E5jlxZh6rZAx5zLRbwcrsrf7OhRjjPG61iqpnxWRloa8OCAiN4nI1V6Kq1eSkefgh6I7FlNX39aEesYY07e11sz1T8DP3CSxASgCQoBROLPCvQDM9XqEvcnQKdQGRjKteg1rdh9kWnqsryMyxhivaa0n9RrgChGJAKbhTDlaDWSrak7PhNfL+AdAxhmcnvMFr+XstwRhjOnXWu0oB6CqlTiT9hggcPQ5JG9ZwLbsNXCBTUNqjOm/2l1JbVzusBuDiz6nuPKwj4MxxhjvsQTRUTHpHBmUziy/dby1xnpVG2P6r44M1hfmzUD6kqDMc5npn83rX2xDdcB2LjfG9HNtJggRmSEim4DN7vokEXnS65H1ZiPOIZQahpSsYJnNVW2M6afaU4L4LXABUAygqmuB070ZVK838hw0NJZrghYx94tdvo7GGGO8ol2PmFS1+dAaA3tShIBgZPJVnC1ZZG3czP7yGl9HZIwx3a49CWK3iMwAVEQCReSHQHZbJ4nICyKyX0Q2NNk2SUSWich6EXlHRAa1cO5sEckRkW0i8kC776YnTbkef63nW7KY11cMvKGpjDH9X3sSxHeBO4FkoACY7K635SVgdrNtzwEPqOoE4E3g/uYniYg/Ti/uC4ETgDkickI7rtezEkZD2kxuDF3Cq8tzbegNY0y/02aCUNUDqnq1qg5W1URVvUZV26yZVdUlQEmzzaOBJe77j4DLPJw6HdimqjtU9QjwGnBxW9fziSnXM7huD2mVq/l4835fR2OMMd2qxZ7UIvIErUwQpKp3d+J6G3F+2c8HLgeGeTgmGWj6zCYfZ6jxluK8DbgNIDU1tRMhdcEJ30Tf/xE3+S3ir8vP5PxxQ3r2+sYY40WtlSCycCYLCgGmAFvdZTIQ1Mnr3QR8T0RW4sxWd6STn3OUqj6jqtNUdVpCQkLbJ3SnwFBk0hzO1uVs2LKdvOJDPXt9Y4zxotZmlHtZVV8GJgJnquoTqvoEcA5OkugwVd2sque781u/CniaWKGA40sWKe623mnq9fhrHd8OWMory63JqzGm/2hPJXUMzvDejSLcbR0mIonuqx/wU+ApD4etAEaJSIaIBAHfAd7uzPV6ROJYGHYyN4Uu5u8rdlFTO7BbABtj+o/2JIjHgNUi8pKIvAysAn7Z1kki8iqwDMgUkXwRuRmnRdIWnF7Ze4AX3WOHish7AKpaB/wHsBCnOe3fVXVjx2+tB025niG1+YyuWc/7Gwp9HY0xxnQLac9YQiIyhGMVxctVda9Xo+qkadOmaVZWVs9f+EgV+ptM/lV/Ik/FPcC8O2b0fAzGGNMJIrJSVad52teesZhOx2meWuouo91tplFQGDLxSs5qWMa2vN1kF5b7OiJjjOmy9jxiur/J8jPgHeBhL8bUN029noCGI1wR+CnPLt3h62iMMabL2tNR7htNlvOA8TglCdPUkAmQPJXbIpby5up8K0UYY/q8zkwYlA+M7e5A+oUp15NQvZPTg3fw2PubfR2NMcZ0SZtzUjfrUe2H0wdilRdj6rvGXwYLH+ThqEWctWUEn249wKxR8b6OyhhjOqU9JYjGHtUrcZqt/lhVr/FqVH1VcASc+j0yiv7NOVF7+J/3s2losBnnjDF9U3sSRHRjr2pVnauqn4nIPV6PrK869U4IjeGXUfPZuKect9favNXGmL6pPQnieg/bbujmOPqPkCiY9X0G7/+UKxPy+NXCHA7XWe9qY0zf02KCEJE5IvIOkCEibzdZPuGrw3ibpqbfCpFJPBjyDwoOVvHXZXm+jsgYYzqstUrqz4FCIB74TZPtFcA6bwbV5wWGwhk/YtCC+7h72E6e+DiIy6cOIyos0NeRGWNMu7U2mmueqi5S1VNVdXGTZZU7XpJpzYnXQkwG32t4hYqawzy5aJuvIzLGmA5p7RHTp+5rhYiUN1kqRMR6gbXFPxDOepCQ4k38YsRWXvw8l/zSKl9HZYwx7dZaCWKW+xqpqoOaLJGqOqil80wT4y+DxHF8p/IvBFLHbz7c4uuIjDGm3drVk1pE/N0huVMbF28H1i/4+cE5PyPg4E4eH72RN1cXsCjH5q42xvQN7RnN9S5gH/AR8K67LPByXP3H6NmQMp3zi15ifGIgP563joNVXZ5p1RhjvK49JYh7gExVHaeqE9xlorcD6zdE4JyfIxWFPDt2HcWVR/j5W717/iNjjIH2JYjdQJm3A+nXMk6D4WeRtPYJfnJaDG+v3cM71sPaGNPLtSdB7AAWichPROT7jYu3A+t3LvxfqK3mxrInmTQsmp+9tYH95TW+jsoYY1rUngSxC6f+IQiIbLKYjkjIhDMfwG/z2zw9ZRc1tfX8eN462jPlqzHG+EKbw32r6iM9EciAMONu2PQWQ5b+lIfO/gc/WVjIayt2M2e6NQozxvQ+7WnF9E6zsZjeFpG/isg9IhLSynkviMh+EdnQZNskEVkmIuvdz/XYn0JEct1j1ohIVudurRfyD4BLnoSaMr5T/Cdmjozj0QWb2FVsHeiMMb1Pe+sgKoFn3aUcZzym0e56S14CZjfb9hzwgKpOAN7Emee6JWep6mRVndaOGPuOwePg9PuRDf/kiRML8fcTfvCPNdTbvBHGmF6mPQlihqpeparvuMs1wEmqeicwpaWTVHUJXx31dTSwxH3/EXBZZ4Lu82bdB4PHE/vJA/xydgorckv51cIcX0dljDHHaU+CiGjac9p9H+GudrTH10bgYvf95cCwFo5T4EMRWSkit7X2gSJym4hkiUhWUVFRB8PxkYAguPhPcKiIiwqf4OqTU3lq8Xbmrcz3dWTGGHNUexLED4BPReQTEVkELAV+KCLhwMsdvN5NwPdEZCVOS6iWEswsVZ0CXAjcKSKnt/SBqvqMqk5T1WkJCQkdDMeHhk6GWfcia1/lkRP2MGNEHD95Yz0r82yqDWNM7yDtaWYpIsHAGHc1R1Xb1YBfRNKBBao63sO+0cDfVHV6G5/xMFCpqr9u63rTpk3TrKw+VKdddxieOg2OVFJ2w2Iufn4DlYfreOs/ZpEcHerr6IwxA4CIrGyprrddg/UBo4BMYBJwhYhc18lAEt1XP+CnwFMejgkXkcjG98D5wIbmx/ULAcFwyZ+hch9RH9zFc9dN5XBdA7e8nMWhwzblhjHGt9rTzPUh4Al3OQv4P+Cb7TjvVWAZkCki+SJyMzBHRLYAm4E9wIvusUNF5D331ME4j7TWAl8C76rqBx2+s74iZSpc8EvY8j4js5/kiTknkrO3nPteX0ODtWwyxvhQm4+YRGQ9TslhtapOEpHBOI+GzuuJADuizz1iaqQK8++Ata/CnNd4oWgMv1iwiTvPGsH9F4xp+3xjjOmkrj5iqlbVBqDO7di2n5ZbH5nOEIGLfgtJk+CN27hxTB1zpg/jT59YyyZjjO+0J0FkiUg0Tqe4lcAqnEdHpjsFhsKVfwP/QOS1q3nkgjRmjozj/n+uZcE6G/nVGNPz2kwQqvo9VT2oqk8B5wHXq+qN3g9tAIpOhW+/CMVbCVpwJ89eO5VpabHc89oaPthQ6OvojDEDTHunHJ0oIt/E6Tk9UkQu9W5YA9jwM+C8RyH7HcK+/AMv3HgSk4dF8x+vrOajTft8HZ0xZgBpTyumF4AXcIbF+Ia7XOTluAa2U++E8d+Gfz9KRN7HvHTjSYxLjuJ7c1fyyWab09oY0zPa04ppk6qe0EPxdEmfbcXkyZEqeOECKN4G186nLGEK1zy3nJx9FTx33TROH92Heo0bY3qtrrZiWiYifSJB9CtBYXDNPIhMgrmXE1W2mb/ePJ2RCRHc+pcsPtt2wNcRGmP6ufYkiL/gJIkcEVnnztOwztuBGSAiEa6bD8ER8NdvEV29m7/dcjIZ8eHc9NIKq7g2xnhVexLE88C1OHM7NNY/fMObQZkmolPh2vmgDfCXi4mt288rt57CCUMHccfcVbz02U5fR2iM6afakyCKVPVtVd2pqnmNi9cjM8ckjIZr3oCaMvjLJcRqGa/ccgrnjR3Mw+9s4pfvZduwHMaYbteeBLFaRF4RkTkicmnj4vXIzPGGToarXoeyfPjbpYQ2VPLna6Zy3alpPLNkB/e8vobDdfW+jtIY04+0J0GEAodxRlW1Zq6+lDYDrvwr7N8Ecy/H/0g5j3xzHA9cOIZ31u7huue/pKy61tdRGmP6iXbNB9FX9Ktmrq3ZOB/m3QyJY51HTxGJvLWmgB/+Yy0Z8eE8d91JpMaF+TpKY0wf0B3zQZjeZNwlMOd1OLANXpgNpXlcPDmZl2+azt6yGi56Yikfb7Ze18aYrrEE0VeNOheuewuqDjhJYv9mZoyIZ8Fdp5ESE8ZNL2Xx+Ic51FvltTGmkyxB9GWpJ8MN74HWw4uzIX8lqXFhvPG9GVw+NYU/fLyNG178kpJDLU39bYwxLWvPUBvf97C5DFipqmu8EVRnDZg6iOZKdsBfvwWVRfCduTDiLFSV11bs5qG3NpIQGcyTV09h0rBoX0dqjOlluloHMQ34LpDsLrfjdJp7VkR+1G1Rms6LHQ43LYSYdHjlCljzKiLCnOmp/POOUwG4/KllvPx5Lv2pUYIxxrvakyBSgCmq+gNV/QEwFUgETgdu8GJspiMih8CN78Kwk2H+d2Hhg1Bfx8SUaBbcNYsZI+N46O2NXP/iCvaV1/g6WmNMH9CeBJGI0w+iUS0wWFWrm203vhYaA9e+CdNvg2V/hFcuh+pSYsKDePGGk3j0kvGs2FnC+b9dwjtrbZY6Y0zr2pMg5gLLReQhEXkI+Ax4RUTCgU0tnSQiw0TkExHZJCIbReQed3usiHwkIlvd15gWzr/ePWariFzfiXsbmPwD4Wu/gm/8AXYuhWfPhqIcRIRrT0nj3btnkR4fzl2vruae11ZTVmUd64wxnrWro5yInATMcFc/U9U2a4JFJAlIUtVVIhKJM5/1JTiPpUpU9TEReQCIUdUfNzs3FsjCqf9Q99ypqlra2jUHbCV1S3Z9Aa9fA7U1cNlzkDkbgLr6Bp5ctJ0//Hsr8RHB/OryiZw2yuaXMGYg6o6OcquAfwBvAvtFJLWtE1S1UFVXue8rgGycSu6LgZfdw17GSRrNXQB8pKolblL4CKdi3HRE6ilw2yKIGwGvfgcWPQYN9QT4+3H3OaN443szCA/259rnv+T7f19DcaU9MTTGHNOeKUfvAvbh/JJeALzrvrabiKQDJwLLceovGicy2AsM9nBKMrC7yXq+u83TZ98mIlkiklVUVNSRsAaGqBS46QOYeCUs+h94+RvOgH/AxJRo3r37NP7jrJG8s3YP5zy+mL+v2G0tnYwxQPtKEPcAmao6TlUnquoEVZ3Y3guISAQwD7hXVcub7lPnN1GXfhup6jOqOk1VpyUk2GMSjwJD4dKn4VtPQ+Fa+PNMyH4HgJBAf354QSbv3X0aoxIj+NG8dVz5zBds21/h46CNMb7WngSxG6djXIeJSCBOcpirqm+4m/e59RON9RT7PZxaAAxrsp7ibjNdMek7cPsSp7/E69fAgu9DbTUAowZH8vptp/K/l00gZ28FF/5+Kb9emEP1ERtC3JiBqj09qZ8HMnEeLR19SK2qj7dxnuDUMZSo6r1Ntv8KKG5SSR2rqj9qdm4sTsX0FHfTKpxK6pLWrmmV1O1UdwQ+/gV8/gQkjIVvvwCDj007fqDyML98N5s3VheQFBXCj2ZncvGkZPz8xIdBG2O8oauV1Ltw6h+CgMgmS1tm4kxVeraIrHGXrwGPAeeJyFbgXHcdEZkmIs8BuIngUWCFu/yireRgOiAgCM7/L7hmnjPY3zNnwtLfQH0dAPERwTx+5WT+fvupxEcEc9/ra/nWnz8nK9e+AmMGEpsPYqCr3A/v/gCy34akSXDxn2DIhKO7GxqUN1YX8KuFm9lXfpivT0zigdljGBZr800Y0x+0VoJoMUGIyO9U9V4ReQcPFcmq+s3uDbPrLEF0wcb58N4PoboUZn0fTv8hBAQf3V11pI6nFu/gmSXbaVC4YUY63z1jBLHhQb6L2RjTZZ1NEFNVdaWInOFpv6ou7sYYu4UliC6qKoEPHoB1rzt1Exf/CVKmHndIYVk1v1qYw5urCwgPCuCmWRnccloGg0ICfRS0MaYrOpUgmpx8j6r+vq1tvYEliG6y5UNYcC9UFMK0m+HsB51xnprYuq+Cxz/awvsb9hIVGsjtZwznhhnphAUF+CZmY0yndDVBrFLVKc22rVbVE7sxxm5hCaIb1ZTBx/8FK55zksO5D8Pka8Dv+HYNGwrK+M2HOXySU0R8RDDfO3MEc6anEhrk75u4jTEd0tlHTHOAq4BZwNImuyKBBlU9p7sD7SpLEF6wdz28dz/sWgbJU52BAJOnfuWwrNwSfv1hDl/sKCEuPIibZmVw7alp9ujJmF6uswkiDcgA/gd4oMmuCmCdqtZ1d6BdZQnCS1Rh3d/ho585rZ6mXAfn/BzC479y6Jc7S3hy0TYW5RQRGRzAdTPSuGlmBnERwR4+2Bjja116xNSXWILwsppyWPy/sPwpCAyDmffAKd+DoK82ed1QUMaTi7bx/oa9BAf48Z2TUrl5VoY1jzWml+lsCeJTVZ0lIhUc38xVcIZRGtT9oXaNJYgeUrQF/vUw5LwLkUlw5k9g8tXg/9UK6m37K3lq8Xbmry6gQZXZ44dw86zhTE3zOA2IMaaHWQnCeEfeMvjo55D/JSSMgXMegswLQb46JEdhWTUvf57HK8vzKK+pY/KwaG6elcGF44cQ4N/eUeeNMd2tywnCnfVtGHD0T8TGuR56E0sQPqAKmxc4JYribZB6Kpz1n5B+msdEcehwHfNW5fPCpzvJLa5iaFQI156azhXTUqyewhgf6Goz10dxZoHbATS4m1VVz+7OILuDJQgfqq+FVX+BJb9y+k+kzXQePWWc5vHwhgbl4837ef7TnSzbUUyQvx9fmzCEa05JY2paDOIhuRhjul9XE0QOMEFVj3gjuO5kCaIXqK2BVS/D0sehci+kzYIzH2gxUYDT6W7u8l3MW5lPxeE6xgyJ5JpT0rjkxGQigq3jnTHe1NUEMQ+4Q1U9zdvQq1iC6EU8JYrTvg8jzvb46Amc8Z7eWrOHv32Rx8Y95YQF+XPRxCSuPGkYU1KtVGGMN3Q1QUwD3gI2cPx8EDZYn2lbbTWsfBk+/a2TKIZMdJrHnnCJx1ZPAKrK6t0Hee3LXSxYV0jVkXqGJ4RzxbRhXDolmcTIkJ69B2P6sa4miI3A08B6jtVB2GB9pmPqDjuDAH72ByjeCtFpMOMup3msh34UjQ4druPd9YX8fcVusvJK8fcTzspM5LIpyZw9NpHgABvSw5iu6GqCWKGqJ3klsm5mCaIPaGiAnPfgs99B/goIi4OTbnEGBYwc3Oqp24sq+UdWPvNW5VNUcZhBIQF8feJQvnViMtPSYmzGO2M6oasJ4nGcR0tvc/wjJmvmajpP1Rnf6dPfwdaF4BcI4y+Fk78LyVNaPbWuvoHPtxfz5uoCPtiwl+raelJiQrlkcjKXnDiUkYntmfDQGANdTxCfeNhszVxN9yneDsufhjVz4UglpEyHU74LY78J/q0P9nfocB0fbdrHG6sL+HRrEQ0KmYMjuWhiEhdNGkpGfHgP3YQxfZP1pDZ9Q025kySWPw2lO51hPKZc5yxRKW2evr+ihvfX72XBuj2syC0FYNzQQVw0cSgXTUyycaCM8aCzYzFdo6p/E5Hve9qvqo93Y4zdwhJEP9HQAFs/dOai2PYvp1nsqPNh2k0w8lzwa7ties/Bat5bX8g76wpZu/sg4CSLC8cPYfb4IfYYyhhXZxPE7ar6tIg85GG3quov2rjoMOAvwGCcwf6eUdXfi0gs8DqQDuQCV6hqqYfz63FaTgHsak+zWksQ/VBpntOfYtVf4dB+iBrmlCgmX9WuUgXA7pIq3ltfyMKNe1m16yAAIxLCmT1+CLPHJTE+eZD1sTADVrc/YhKRe1X1d20ckwQkqeoqEYkEVgKX4AzbUaKqj4nIA0CMqv7Yw/mVqhrRkbgsQfRj9bVO66esF2DHIkBgxFlOM9kxF0Fg+/pG7C2r4cNNe/lgw16W7yyhvkFJigrhnLGJnDt2MKeOiLOms2ZA8UaC2KWqqR085y3gj+5ypqoWuklkkapmejjeEoTxrGQnrHkF1r4KZbshJArGfxtOvBqGTmmxp/ZXPubQEf6VvY9/Z+9jyZYDVNfWExbkz+mjEjj3hMGclZlgAwiafs8bCWK3qg7rwPHpwBJgPM7jomh3uwCljevNzqkD1gB1wGOqOr+Fz74NuA0gNTV1al5eXgfuxPRpDQ2QuwRWz4Xst6GuBuIzYeLlMOFyiElv90fV1NazbHsxH7kJY1/5YURgYko0Z2UmcFZmIhOSo6yvhel3fFqCEJEIYDHw36r6hogcbJoQRKRUVb8ye4yIJKtqgYgMBz4GzlHV7a1dy0oQA1hNGWx4w5kaddfnzrZhJzuJYtylEB7X7o9SVdYXlPHJ5iI+ydnP2vyDqEJceBBnZCZwZmYip42MJyY8yEs3Y0zP6WwldfOZ5I7uAkJVtc1hNkUkEFgALGxs9eSODtvmI6Zmn/MSsEBV/9nacZYgDAAHd8H6fzrJoigb/AJgxDlOR7zMC51HUh1QXHmYJVuLWJRTxOItRRysqkUEJiRHcfqoBE4bFc+JqTEEBdjER6bv8Uk/CPfx0cs4FdL3Ntn+K6C4SSV1rKr+qNm5MUCVqh4WkXhgGXCxqm5q7ZqWIMxxVGHfBidRbHgDyvPBP8hpKjvuWzB6NoR0bObc+gZlbf5Blm45wNKtRazefZD6BiU8yJ9TR8Qxc2Q8M0fGMyoxwlpGmT7BVwliFrCU4wf5+09gOfB3IBXIw2nmWuKOGvtdVb1FRGbgDBDYAPgBv1PV59u6piUI06KGBijIgo3zYeObULEH/INh1HlOj+3RF0BodIc/trymlmXbi1mypYilWw+wq6QKgITIYGaOiGOGmzCSo0O7936M6SbWk9qYphoanIECN74Jm+Y7M+D5BUDG6U6T2TEXtTlwYEt2l1Tx+fYDfLatmM+3H+BApTPPVlpcGKdkxHHKiFhOGR5HUpQlDNM7WIIwpiUNDbBnldMKKvsdKNkBCAybDmO+Dplfg/hRnfpoVSVnXwWfbSvmix3FLN9RTHlNHQDpcWGcMjyOk4fHMj0jzkoYxmcsQRjTHqqwPxs2L3ASxl63I3/sCKdyO/NCGHZKixMdtaW+Qdm8t5xl24v5YkcJX+48ljCSo0M5KT2G6RlxTM+IYUSC1WGYnmEJwpjOOLgLtiyEnPchdynUH4GQaGdcqFHnOS2jOtB8trnGhLFiZwlf5pbw5c5SDlQ6I+rHhgcxNS2GaWkxTE2LYUJKlPXwNl5hCcKYrjpcAds/hpwPnIEEqw4AAslT3YRxLiSdCH6db+qqquQWV/HlzmK+3FnKyrwScoudSu8gfz8mpEQxLS2GKWkxTEmNISHSenmbrrMEYUx3amiAwtWw9SNnKVgJKITFw4iz3eUsiBzS5UsVVRxmZZ6TLLLyStlQUEZtvfN/NiUmlCmpMUxJjebE1BjGJg2yvhimwyxBGONNhw44pYutH8L2T9zSBTB4vJMoRpwNqadCYNcromtq69m4p4xVeQdZvbuUVXkH2VteA0BQgB/jhw5i0rBoJrtLamyY1WWYVlmCMKanNDTAvvVOwtj+Mez6wqm78A+G1JMh4wwYfiYkTe50ZXdzhWXVrN51kFV5pazNP8j6gjJqap2uR9FhgUxKiWZSShQTU6KZmBJF4qD2jXxrBgZLEMb4ypFDkPsZ7FwMOxY7yQMgeBCkz3ISRvosSDyhS/UXTdXVN7BlXyVrdh9k7e6DrNl9kK37K2hw/6sPGRTChJQoJiZHMSElignJUTZq7QBmCcKY3qKyyBmBdsdiJ2mU5jrbQ2MhfSakn+YkjISx3ZYwAKqO1LFpTznr8stYl3+QdQVl7Cg6dHR/UlQI45OjGD80igkpgxg/1EoaA4UlCGN6q4O7nBJG7lJnObjL2R4aC2kznCX1VBgysdseSTUqr6llY0E5G/eUsb6gjA0FZew4cIjGXwnxEcGMGzqIE4YOYtzQQYwbGkVabJgNed7PWIIwpq8ozYO8zyD3U+e1sYQRFOH07k6dAWmnOhMjBYV1++UrD9eRXVjO+vwyNhWWs3FPOVv3VVDnPp8KD/JnTNIgxiZFMjZpEGOTBjFmSCRhQd2bvEzPsQRhTF9VvgfyPoddyyBvGezf6Gz3C4CkSU7P7tSTnddOjh/VlsN19WzdV8mmPU5pI7uwguzCcioOO73ARSA9LpyxSZGMGTKIzCGRjB0yiJSYUCtt9AGWIIzpL6pKYPeXsPsL2LXcGUeqzmnmSkw6pEyHlJNg2ElOM1v/QK+Eoarkl1aTXVjOpsJysgvL2by3gl0lVUcfUYUH+TN6SCRjhkSSOTiS0e6rVYj3LpYgjOmv6o5A4VrYvdxJGrtXQOVeZ19AKAw90UkWydMgZRoMGurVcA4drmPLvgpy9laweW8Fm/c6ieNgVe3RY+IjgskcEsHowZHuEsHIxEiiQr2TzEzrLEEYM1CoQlk+5H8J+VlOaaNwLTS4v6Ajh0LyFCdZJE91EkhwpJdDUooqDpPjJo4t+yrI2VfJ1n0VVB2pP3rc4EHBjB4cyajESEYNjmBUYgQjEyOIDrOpXb3JEoQxA1ltjTOzXn6WM2lSfhaU7nR3CsSPdpLG0BOdyu8hEyDQ+01cGxqUgoPVbNlXwdb9lc7rvkq27a+kuvZY4oiPCGZkYjijEiMZ6SaNEQkRDB4UbL3Eu4ElCGPM8apKnDGkClY59RgFq+DQfmefXwAkjnV6ew+d7AxCOHhcjyQNOJY4tu2vZOt+N2kUVbJtX+XRinGAiOAARiSEMyIhghGJEYxICGd4QgRpcWE28m0HWIIwxrRO1WkxtWcV7FntJIzCNVBd6uz3C3A67yVNcpeJTiV4cEQPhqjsKz/MjiInYWzfX8n2okNsL6qksKzm6HF+AikxYUcTxvCEcDLiwxkeb6UOTyxBGGM6ThXKdsOeNU6yKFzrvG8cjBCBuBFOJ76kic7rkAkQkdjjoVYermNn0SF2HDiWNHYUHWLngcqj41IBhAX5kxHfmDDCSXeXjLhwYsIHZl2HJQhjTPdQdebwLlwLhetg7zrntWzXsWMiBjuliyETnGXweIgb2e09wdujoUEpLK9hp5ssdhw4xE532V1SdXR8KoCo0EA3WYQ5iSMunLS4MNL7efLwSYIQkWHAX4DBgALPqOrvReRy4GFgLDBdVT3+RheR2cDvAX/gOVV9rK1rWoIwxkeqSpyK8L3rYe8GZ1DC/ZuPtZ7yD4bEMU6yGDzOXcZDeLzPQj5cV8/ukmpyDxwit9hJGrnFh8g9UMWesmq0efKICyM1Lpy02DDS4sJIcxNIYmTffmzlqwSRBCSp6ioRiQRWApfgJIsG4Gngh54ShIj4A1uA84B8YAUwR1U3tXZNSxDG9CJ1R+DAFidp7N8I+9ylct+xY8ITnQrxweOcEW0TT3ASSVC47+LGmXcjv7SK3ANV5BYfIq/Yec0tPsSegzXUNyl6hAT6kRobRmpsOKlu8kiNDWNYbBjDYkN7fYV5awnCa2U+VS0ECt33FSKSDSSr6kduUK2dPh3Ypqo73GNfAy4GWk0QxpheJCAIhox3lqYqi5yEsXcD7M923me9CHXVx46JTjuWLBLGOq/xo7tl0qX2CAn0Z2RiJCMTv9pHpLa+gYLSavJKqthVfIjc4iryiqvYXVLFZ9sOHNdEV8QZXn1YTBgpsaFO4ogJIzXOeU2MDO7Vw5H0yENBEUkHTgSWt/OUZGB3k/V84OQWPvs24DaA1NTUzgdpjOkZEQkQcaYzcVKjhgY4mAv7NsH+TU7iKNoM2/517DGV+EFMBiSMgYRM93W0kzh6sMQR6O93tHIbEo7bp6oUVR5md4mTNPKKq9hdWkV+STWfbyvmzYqC4x5dBQX4kRIdSnJMqFPiiHFKHSkxYaTEhBIXHuTTx1deTxAiEgHMA+5V1fLu/nxVfQZ4BpxHTN39+caYHuDnB7HDnWXsRce219dC8XYoynbqNIqyoSgHti6EhmN9IohOhfhMJ3HEjz72Ghbbo7chIiRGhpAYGcLUtK9e+3BdPQWl1ewqqWJ3aTX5JW4CKa1mw/pCSpsMSQLO46vGZJEc7SSO5KPvQ0mI8G4JxKsJQkQCcZLDXFV9owOnFgDDmqynuNuMMQOJf6DzeClxDIxrsr2+Fkp2OqWMohw3cWxx5tSoO9YngrB4N1mMgrhRTtKIH+k8wvLr+bqB4AB/t2+G5/4jFTW15JdWU1BaTb6bOPJLq8k/WMWa3QePG9MKIMjfj6HRIYxIiOD5G07q9ni9liDEKRc9D2Sr6uMdPH0FMEpEMnASw3eAq7o5RGNMX+Uf6DxeShh9/PaGBqfJbdEWOJDjJI8DWyH7HagqbnJ+sFNaiR/pJI64kW4SGdnjpY6mIkMCGZsUyNikQR73Vx6uY89BJ3kUlFaTf9BJJt56dOLNVkyzgKXAepxWSwD/CQQDT+A8vDsIrFHVC0RkKE5z1q+5538N+B1OM9cXVPW/27qmtWIyxrSoqsRJFge2QPFW9/1WZ1yqpo+rQmOdRBE3wlliRzjrscN7tOd4T7GOcsYY05L6OjiYB8XbnIRRvNWp9yjeDhV7jj82YoibNDKcxBE73FmPyeizycMnzVyNMaZP8A84VloYfcHx+44cgpIdbsLY5rwv2QFbPjw2uGGjiMFOwojJcCvcM5wlJsOnj626whKEMca0JCj82JAhzR2ucCrKS9zSRulOKMmFHYtg7SvHHxsS5SSKmPRjSSMm3VkGJftkGJL26J1RGWNMbxcc6QxSmDTxq/uOVDmPrUp2OEmkdKfzuncdbF5wfJ2HXwBEDTuWMGLSISbNaWkVkw6hMU6POx+wBGGMMd0tKMwZQiRx7Ff31ddBeYGTNErzoDT32LLpLaguafZZkU0ShvsaneosMWlenRHQEoQxxvQk/wDnF3tMmuf9NeVO6aM07/jXku2w4xOorTr++NAYp1f5TR90e6iWIIwxpjcJGdRyvYeq05/jaOLY5SwNtV89thtYgjDGmL5CxBkiPTwekqd6/XJ+Xr+CMcaYPskShDHGGI8sQRhjjPHIEoQxxhiPLEEYY4zxyBKEMcYYjyxBGGOM8cgShDHGGI/61XwQIlIE5LVySDxwoIfC6W0G6r3bfQ8sdt8dl6aqCZ529KsE0RYRyWppYoz+bqDeu933wGL33b3sEZMxxhiPLEEYY4zxaKAliGd8HYAPDdR7t/seWOy+u9GAqoMwxhjTfgOtBGGMMaadLEEYY4zxaMAkCBGZLSI5IrJNRB7wdTzdSUSGicgnIrJJRDaKyD3u9lgR+UhEtrqvMe52EZE/uD+LdSIyxbd30DUi4i8iq0VkgbueISLL3ft7XUSC3O3B7vo2d3+6TwPvAhGJFpF/ishmEckWkVMHwvctIve5/8Y3iMirIhLSX79vEXlBRPaLyIYm2zr8HYvI9e7xW0Xk+o7EMCAShIj4A38CLgROAOaIyAm+japb1QE/UNUTgFOAO937ewD4t6qOAv7troPzcxjlLrcBf+75kLvVPUB2k/X/BX6rqiOBUuBmd/vNQKm7/bfucX3V74EPVHUMMAnn/vv19y0iycDdwDRVHQ/4A9+h/37fLwGzm23r0HcsIrHAQ8DJwHTgocak0i6q2u8X4FRgYZP1nwA/8XVcXrzft4DzgBwgyd2WBOS4758G5jQ5/uhxfW0BUtz/KGcDCwDB6VEa0Py7BxYCp7rvA9zjxNf30Il7jgJ2No+9v3/fQDKwG4h1v78FwAX9+fsG0oENnf2OgTnA0022H3dcW8uAKEFw7B9Wo3x3W7/jFqNPBJYDg1W10N21Fxjsvu9PP4/fAT8CGtz1OOCgqta5603v7eh9u/vL3OP7mgygCHjRfbT2nIiE08+/b1UtAH4N7AIKcb6/lfT/77upjn7HXfruB0qCGBBEJAKYB9yrquVN96nz50O/atMsIhcB+1V1pa9j6WEBwBTgz6p6InCIY48agH77fccAF+MkyKFAOF99BDNg9MR3PFASRAEwrMl6irut3xCRQJzkMFdV33A37xORJHd/ErDf3d5ffh4zgW+KSC7wGs5jpt8D0SIS4B7T9N6O3re7Pwoo7smAu0k+kK+qy931f+IkjP7+fZ8L7FTVIlWtBd7A+TfQ37/vpjr6HXfpux8oCWIFMMpt7RCEU7H1to9j6jYiIsDzQLaqPt5k19tAY6uF63HqJhq3X+e2fDgFKGtSbO0zVPUnqpqiquk43+nHqno18Anwbfew5vfd+PP4tnt8n/srW1X3ArtFJNPddA6wiX7+feM8WjpFRMLcf/ON992vv+9mOvodLwTOF5EYtwR2vrutfXxdCdODlT1fA7YA24EHfR1PN9/bLJyi5jpgjbt8Ded567+BrcC/gFj3eMFp1bUdWI/TKsTn99HFn8GZwAL3/XDgS2Ab8A8g2N0e4q5vc/cP93XcXbjfyUCW+53PB2IGwvcNPAJsBjYAfwWC++v3DbyKU9dSi1NqvLkz3zFwk/sz2Abc2JEYbKgNY4wxHg2UR0zGGGM6yBKEMcYYjyxBGGOM8cgShDHGGI8sQRhjjPHIEoTps0SkXkTWNFn63Ci9IvKfzdY/99J1ThSR51vZnyAiH3jj2qbvsmaups8SkUpVjfB1HM25nbhEVRvacWyP3IOI/AP4L1Vd28oxLwLPqepn3o7H9A1WgjD9iohEiTPvR6a7/qqI3Oq+rxSR37rzCfxbRBLc7ZNF5At3HP03m4yxf7c4c2ysE5HX3G0Pi8gPm1xvg4iku0uOiPwFpxPXMBG5X0RWuOc/4iHWx4BQt/QztzFG9/VMEVksIm+JyA4ReUxErhaRL0VkvYiMcI9LEJF57nVWiMhMD9eJBCY2JgcROaNJqWu1ux+cDndXd8PXYPoLX/cWtMWWzi5APcd6jq8BrnS3nwcswxl+44Mmxytwtfv+58Af3ffrgDPc978Afue+38OxXrnR7uvDwA+bfOYGnCGZ03FGlD3F3X4+zkTygvOH2ALgdA/3UOlpHadn+EGcIZuDccbPecTdd0+TGF8BZrnvU3GGW2l+jbOAeU3W3wFmuu8jODZUdjKw3tffqy29Z2kc4MqYvqhaVSc336iqH4nI5ThDD0xqsqsBeN19/zfgDRGJwvnlv9jd/jLO8AzgJI65IjIf56/rtuSp6hfu+/PdZbW7HoEzmcuSdnxOoxXqjpkkItuBD93t63F+6YMzgN0JzlMtAAaJSISqVjb5nCSc4cEbfQY87pZa3lDVfHf7fpxRUo0BsARh+h8R8QPGAlU4YxTlt3BoWxVwXwdOB74BPCgiE3Bm72v6aDakyftDTcMA/kdVn+5A6M0dbvK+ocl6A8f+7/rhlFpqWvmc6qZxqupjIvIuznhdn4nIBaq62T2mugvxmn7G6iBMf3QfzhScV+FMqhPobvfj2KifVwGfqmoZUCoip7nbrwUWu0lmmKp+AvwYZ6joCCAXZ2htxJn3N6OFGBYCN4kzRwcikiwiiR6Oq20SX2d8CNzVuCIikz0ckw2MbHLMCFVdr6r/izPS8Rh312icR2bGAFaCMH1bqIisabL+AfAicAswXVUrRGQJ8FOceXkPAdNF5Kc4j1OudM+7HnhKRMKAHcCNOPMd/819BCXAH1T1oIjMwxlWeSPOrH1bPAWmqh+KyFhgmfv4pxK4hmPj9zd6BlgnIqvUGaq8o+4G/iQi63D+Py8Bvtssls1u5X2kqlYA94rIWTglkY3A++6hZwHvdiIG009ZM1czYPTWZrE9QUTuAypU9blWjlkCXKyqpT0XmenN7BGTMQPDnzm+TuM4bpPfxy05mKasBGGMMcYjK0EYY4zxyBKEMcYYjyxBGGOM8cgShDHGGI8sQRhjjPHo/wEfd1xZxx+91wAAAABJRU5ErkJggg==\n" }, "metadata": { @@ -142,7 +149,7 @@ "\n", "for night in [False, True]:\n", " limmags = dorado.sensitivity.get_limmag(\n", - " synphot.ConstFlux1D, snr=5, exptime=exptimes, coord=coord, time=time, night=night)\n", + " synphot.SourceSpectrum(synphot.ConstFlux1D, amplitude=0 * u.ABmag), snr=5, exptime=exptimes, coord=coord, time=time, night=night)\n", " ax.plot(exptimes, limmags, label='night' if night else 'day')\n", "\n", "ax.legend()" @@ -159,13 +166,15 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 7, "metadata": {}, "outputs": [ { "output_type": "stream", - "name": "stdout", + "name": "stderr", "text": [ + "/Users/lpsinger/Library/Caches/pypoetry/virtualenvs/dorado-sensitivity-RYVm8gWH-py3.8/lib/python3.8/site-packages/astropy/units/quantity.py:477: RuntimeWarning: invalid value encountered in multiply\n", + " result = super().__array_ufunc__(function, method, *arrays, **kwargs)\n", "nan\n", "5.000000000000003\n", "5.000000000000003\n", @@ -225,13 +234,24 @@ " synphot.SourceSpectrum(synphot.ConstFlux1D, amplitude=limmag),\n", " exptime=exptime, coord=coord, time=time, night=night))" ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { "kernelspec": { - "display_name": "Python 3", - "language": "python", - "name": "python3" + "name": "python3", + "display_name": "Python 3.8.8 64-bit", + "metadata": { + "interpreter": { + "hash": "e056a4a7562dea2d16314dfb63f3098ac154a31025c0798943f6e61d571635c6" + } + } }, "language_info": { "codemirror_mode": {