From 303d8b35013b7c2f82a5e41defd2fe2263f63237 Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 08:19:38 -0400 Subject: [PATCH 01/10] Test footprint defect flags and final-catalogue propagation first Co-Authored-By: Claude Opus 5.5 --- tests/module/test_make_cat.py | 7 +++ tests/module/test_ngmix_defect_flags.py | 73 +++++++++++++++++++++++++ 2 files changed, 80 insertions(+) create mode 100644 tests/module/test_ngmix_defect_flags.py diff --git a/tests/module/test_make_cat.py b/tests/module/test_make_cat.py index 4a725f019..528b43ac0 100644 --- a/tests/module/test_make_cat.py +++ b/tests/module/test_make_cat.py @@ -60,6 +60,7 @@ def warning(self, message, *_args, **_kwargs): "mcal_types_fail", "neighbour_flag", "n_epoch_interp", "min_dist_interp", "min_dist_noisefill", + "imaflags_iso", "n_epoch_interp_iso", "nfev_fit", "g1", "g1_err", "g2", "g2_err", "T", "T_err", @@ -90,6 +91,7 @@ def _ngmix_row(obj_id): "n_epoch_failed": 1, "mcal_types_fail": 0, "neighbour_flag": 1, + "imaflags_iso": 1028, "n_epoch_interp_iso": 1, "n_epoch_interp": 2, "min_dist_interp": 7.5, "min_dist_noisefill": 14.25, "nfev_fit": 7, @@ -124,6 +126,7 @@ def _write_ngmix_cat(path, obj_ids): name=key, format="K" if key in ("id", "n_epoch_model", "n_epoch_failed", "mcal_types_fail", + "imaflags_iso", "n_epoch_interp_iso", "n_epoch_interp", "nfev_fit", "flags", "mcal_flags") else "D", array=np.array([row[key] for row in rows]), @@ -286,6 +289,8 @@ def test_save_ngmix_data_fills_sentinels_for_absent_objects(tmp_path): for col, key, never_fit in ( ("NGMIX_N_EPOCH_FAILED", "n_epoch_failed", 0.0), ("NGMIX_N_EPOCH_INTERP", "n_epoch_interp", 0.0), + ("NGMIX_IMAFLAGS_ISO", "imaflags_iso", 0), + ("NGMIX_N_EPOCH_INTERP_ISO", "n_epoch_interp_iso", 0), ("NGMIX_MIN_DIST_INTERP", "min_dist_interp", -1.0), ("NGMIX_MIN_DIST_NOISEFILL", "min_dist_noisefill", -1.0), ): @@ -334,6 +339,8 @@ def _metacal_result(obj_id): "n_epoch_failed": row["n_epoch_failed"], "neighbour_flag": row["neighbour_flag"], "n_epoch_interp": row["n_epoch_interp"], + "imaflags_iso": row["imaflags_iso"], + "n_epoch_interp_iso": row["n_epoch_interp_iso"], "min_dist_interp": row["min_dist_interp"], "min_dist_noisefill": row["min_dist_noisefill"], } diff --git a/tests/module/test_ngmix_defect_flags.py b/tests/module/test_ngmix_defect_flags.py new file mode 100644 index 000000000..bf11b2a8b --- /dev/null +++ b/tests/module/test_ngmix_defect_flags.py @@ -0,0 +1,73 @@ +"""Adaptive defect flags on the coadd object's footprint, over kept epochs.""" + +import numpy as np +import pytest + +from shapepipe.modules.ngmix_package.ngmix import Ngmix, prepare_postage_stamps +from tests.module.test_ngmix_defect_fill import N_STAMP, _fake_inputs + + +def _run(kind="inside", seg_kind="valid", megacam_flip=True): + seg = np.zeros((N_STAMP, N_STAMP), dtype=np.int32) + seg[23:28, 23:37] = 1 # reaches 11 px right of centre, no stamp edge + if seg_kind == "empty": + seg[:] = 0 + if seg_kind == "clipped": + seg[25, 36:] = 1 + epochs = {} + for i, (ccd, bit) in enumerate(((10, 4), (20, 1024), (21, 0))): + flag = np.zeros(seg.shape, dtype=np.float32) + weight = np.ones(seg.shape) + pixel = (25, 35) if kind != "outside" else (35, 35) + if kind != "clean": + flag[pixel] = bit + if kind == "zero_weight" and bit: + flag[pixel] = 0 + weight[pixel] = 0 + if megacam_flip: + flag = Ngmix.MegaCamFlip(flag, ccd) + weight = Ngmix.MegaCamFlip(weight, ccd) + epochs[f"210000{i}-{ccd}"] = (flag, weight) + # A centrally vetoed epoch and a fraction-cut epoch carry otherwise + # unique bits on the footprint. Neither may contribute. + for i, fraction in enumerate((False, True), start=3): + flag = np.zeros(seg.shape, dtype=np.float32) + flag[25, 35] = 64 + if fraction: + flag[:10] = 64 + else: + flag[25, 25] = 64 + epochs[f"210000{i}-20"] = (flag, np.ones(seg.shape)) + vignet, tile, psf, gal = _fake_inputs(epochs) + tile.seg = None if seg_kind == "missing" else seg[None] + return prepare_postage_stamps( + vignet, 1, 0, tile, bkg_sub=False, psf_obj=psf, gal_obj=gal, + megacam_flip=megacam_flip, + ) + + +@pytest.mark.decision("shape_measurement.defect_flagging") +@pytest.mark.parametrize("flip", [True, False]) +@pytest.mark.parametrize("kind, expected", [ + ("inside", (1028, 2)), ("outside", (0, 0)), + ("clean", (0, 0)), ("zero_weight", (0, 2)), +]) +def test_footprint_flags_only_count_kept_epochs(kind, expected, flip): + """Keep bit identities across epochs, ignore outside/vetoed defects. + + The final clean epoch must not overwrite the OR of earlier epochs. + A zero-weight pixel has no instrument bit but counts if interpolated. + """ + stamp = _run(kind, megacam_flip=flip) + assert len(stamp.flags) == 3 + diagnostics = stamp.defect_diagnostics() + assert (diagnostics["imaflags_iso"], + diagnostics["n_epoch_interp_iso"]) == expected + assert diagnostics["n_epoch_interp"] == (0 if kind == "clean" else 2) + + +@pytest.mark.parametrize("seg_kind", ["missing", "empty", "clipped"]) +def test_unknown_or_incomplete_footprint_is_not_clean(seg_kind): + diagnostics = _run(seg_kind=seg_kind).defect_diagnostics() + assert diagnostics["imaflags_iso"] == -1 + assert diagnostics["n_epoch_interp_iso"] == -1 From face514d519f6181af674f38b9fde3692273d817 Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 08:35:17 -0400 Subject: [PATCH 02/10] Propagate kept-epoch instrument bits and interpolated defects on object footprints Co-Authored-By: Claude Opus 5.5 --- .../modules/make_cat_package/make_cat.py | 9 ++- src/shapepipe/modules/ngmix_package/ngmix.py | 59 ++++++++++++++++++- workflow/config/cfis/final_cat.param | 3 + .../config/cfis_image_sims/final_cat.param | 3 + 4 files changed, 72 insertions(+), 2 deletions(-) diff --git a/src/shapepipe/modules/make_cat_package/make_cat.py b/src/shapepipe/modules/make_cat_package/make_cat.py index 7be609147..cf756b2be 100644 --- a/src/shapepipe/modules/make_cat_package/make_cat.py +++ b/src/shapepipe/modules/make_cat_package/make_cat.py @@ -319,7 +319,10 @@ def _save_ngmix_data(self, ngmix_cat_path, moments=False): observation raised, see ``ngmix.do_ngmix_metacal`` — ``NGMIX_MCAL_TYPES_FAIL``, ``NGMIX_NEIGHBOUR_FLAG`` — a blend flag set when the coadd seg stamp held a non-central footprint, - shapepipe#776 — and the defect diagnostics ``NGMIX_N_EPOCH_INTERP``, + shapepipe#776 — the adaptive defect diagnostics + ``NGMIX_IMAFLAGS_ISO`` and ``NGMIX_N_EPOCH_INTERP_ISO`` (defined by + ``Postage_stamp.record_footprint_defects``; unknown footprints -1, + never-fit rows 0), and the stamp diagnostics ``NGMIX_N_EPOCH_INTERP``, ``NGMIX_MIN_DIST_INTERP``, ``NGMIX_MIN_DIST_NOISEFILL``, described by ``Postage_stamp.record_defects``; -1 is their absent distance, and a never-fit row carries 0 and -1). The galaxy is the implicit default object @@ -380,6 +383,8 @@ def _save_ngmix_data(self, ngmix_cat_path, moments=False): defect_columns = { "NGMIX_N_EPOCH_FAILED": "n_epoch_failed", "NGMIX_N_EPOCH_INTERP": "n_epoch_interp", + "NGMIX_IMAFLAGS_ISO": "imaflags_iso", + "NGMIX_N_EPOCH_INTERP_ISO": "n_epoch_interp_iso", "NGMIX_MIN_DIST_INTERP": "min_dist_interp", "NGMIX_MIN_DIST_NOISEFILL": "min_dist_noisefill", } @@ -411,6 +416,8 @@ def _save_ngmix_data(self, ngmix_cat_path, moments=False): self._add2dict("NGMIX_NEIGHBOUR_FLAG", np.zeros(n_obj)) self._add2dict("NGMIX_N_EPOCH_FAILED", np.zeros(n_obj)) self._add2dict("NGMIX_N_EPOCH_INTERP", np.zeros(n_obj)) + self._add2dict("NGMIX_IMAFLAGS_ISO", np.zeros(n_obj, dtype=np.int64)) + self._add2dict("NGMIX_N_EPOCH_INTERP_ISO", np.zeros(n_obj, dtype=np.int64)) self._add2dict("NGMIX_MIN_DIST_INTERP", np.full(n_obj, -1.0)) self._add2dict("NGMIX_MIN_DIST_NOISEFILL", np.full(n_obj, -1.0)) diff --git a/src/shapepipe/modules/ngmix_package/ngmix.py b/src/shapepipe/modules/ngmix_package/ngmix.py index 43995f926..6b6516388 100644 --- a/src/shapepipe/modules/ngmix_package/ngmix.py +++ b/src/shapepipe/modules/ngmix_package/ngmix.py @@ -272,6 +272,8 @@ def empty_metacal_output(): 'neighbour_flag', # defect diagnostics (Postage_stamp.record_defects) 'n_epoch_interp', + 'imaflags_iso', + 'n_epoch_interp_iso', 'min_dist_interp', 'min_dist_noisefill', 'nfev_fit', @@ -706,6 +708,8 @@ def __init__( # Defect diagnostics over the epochs kept (see # :meth:`record_defects`). self.n_epoch_interp = 0 + self.imaflags_iso = 0 + self.n_epoch_interp_iso = 0 self.min_dist_interp = np.inf self.min_dist_noisefill = np.inf self.bkg_sub = bkg_sub @@ -742,10 +746,61 @@ def record_defects(self, defect, interpolated): self.min_dist_noisefill, distance[noisefill].min() ) + def record_footprint_defects(self, flag, interpolated, seg, obj_id): + """Accumulate adaptive defect diagnostics for one kept epoch. + + @sc [decision:shape_measurement.defect_flagging,label:schema] footprint-defect-flags + ``imaflags_iso`` is the bitwise OR of all instrument bits in the + exposure flag stamp on ``seg == obj_id``, OR'ed across epochs kept + after the masked-fraction cut and central veto. ``seg`` is the coadd + SEG_VIGNET, already MegaCam-flipped to the epoch frame; it is not an + exposure detection or a reprojected isophote. ``n_epoch_interp_iso`` + counts kept epochs with any interpolated defect on that footprint. + Zero-weight/off-tile pixels with flag 0 contribute no instrument bit + (SExtractor semantics), but count if interpolated. Neither quantity + includes clean quarter-turn copies introduced by defect weighting. + Missing segmentation, an empty footprint, or a footprint touching + any stamp edge makes both quantities -1 (unknown/incomplete), not + clean; this sentinel persists across epochs. Never-fit rows in the + final catalogue carry 0, as do the other count diagnostics. + These diagnostics do not change the fit or impose a catalogue cut. + + Parameters + ---------- + flag : numpy.ndarray + Exposure instrument flags; integral floats are accepted because + vignet extraction can store flags in floating-point arrays. + Non-finite, fractional, negative or int64-overflow values on + the footprint raise ValueError rather than silently losing bits. + interpolated : numpy.ndarray of bool + Interpolated defects, as used by the fill and central veto. + seg : numpy.ndarray or None + Coadd segmentation stamp in the epoch's orientation. + obj_id : int + Tile object's segmentation label. + """ + footprint = None if seg is None else seg == obj_id + if (footprint is None or not footprint.any() + or footprint[[0, -1], :].any() + or footprint[:, [0, -1]].any()): + self.imaflags_iso = self.n_epoch_interp_iso = -1 + return + if self.imaflags_iso == -1: + return + values = np.asarray(flag)[footprint] + if (not np.isfinite(values).all() or (values < 0).any() + or (values >= 2**63).any() + or (values != np.floor(values)).any()): + raise ValueError("Non-integral or out-of-range footprint flags") + self.imaflags_iso |= int(np.bitwise_or.reduce(values.astype(np.int64))) + self.n_epoch_interp_iso += int((interpolated & footprint).any()) + def defect_diagnostics(self): """The defect diagnostic columns, -1 for an absent distance.""" return { "n_epoch_interp": self.n_epoch_interp, + "imaflags_iso": self.imaflags_iso, + "n_epoch_interp_iso": self.n_epoch_interp_iso, "min_dist_interp": ( float(self.min_dist_interp) if np.isfinite(self.min_dist_interp) else -1.0 @@ -1082,7 +1137,8 @@ def compile_results(self, results): # Defect diagnostics (Postage_stamp.record_defects), likewise # object-level and replicated. for key in ( - "n_epoch_interp", "min_dist_interp", "min_dist_noisefill" + "n_epoch_interp", "min_dist_interp", "min_dist_noisefill", + "imaflags_iso", "n_epoch_interp_iso" ): output_dict[name][key].append(results[idx][key]) # ngmix 2.x reports the solver's function-evaluation count @@ -1733,6 +1789,7 @@ def prepare_postage_stamps( stamp.epoch_cuts["central_veto"] += 1 continue stamp.record_defects(defect, interpolated) + stamp.record_footprint_defects(flag_vign, interpolated, tile_seg, obj_id) # One unpickle per exposure (all CCDs), reused across this object's # epochs; the cache is per-call, so bounded by the object's exposures diff --git a/workflow/config/cfis/final_cat.param b/workflow/config/cfis/final_cat.param index bd2f9ba98..261dd0792 100644 --- a/workflow/config/cfis/final_cat.param +++ b/workflow/config/cfis/final_cat.param @@ -178,6 +178,9 @@ NGMIX_NEIGHBOUR_FLAG # min distance (px) from the centre to an interpolated / noise-filled defect # pixel, -1 when none (ngmix Postage_stamp.record_defects) NGMIX_N_EPOCH_INTERP +# @sc [decision:shape_measurement.defect_flagging] +NGMIX_IMAFLAGS_ISO +NGMIX_N_EPOCH_INTERP_ISO NGMIX_MIN_DIST_INTERP NGMIX_MIN_DIST_NOISEFILL diff --git a/workflow/config/cfis_image_sims/final_cat.param b/workflow/config/cfis_image_sims/final_cat.param index 2be87f372..cddeff1a2 100644 --- a/workflow/config/cfis_image_sims/final_cat.param +++ b/workflow/config/cfis_image_sims/final_cat.param @@ -35,6 +35,9 @@ NGMIX_NEIGHBOUR_FLAG # min distance (px) from the centre to an interpolated / noise-filled defect # pixel, -1 when none (ngmix Postage_stamp.record_defects) NGMIX_N_EPOCH_INTERP +# @sc [decision:shape_measurement.defect_flagging] +NGMIX_IMAFLAGS_ISO +NGMIX_N_EPOCH_INTERP_ISO NGMIX_MIN_DIST_INTERP NGMIX_MIN_DIST_NOISEFILL From 73c43b0775ef5c64433b0e016688c2f44209cf33 Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 08:35:18 -0400 Subject: [PATCH 03/10] Record adaptive footprint defect flag semantics in ASTRA Co-Authored-By: Claude Opus 5.5 --- astra.yaml | 26 +++++++++++++++++++++++++- universes/committed.yaml | 1 + 2 files changed, 26 insertions(+), 1 deletion(-) diff --git a/astra.yaml b/astra.yaml index 70d7150dd..30c856690 100644 --- a/astra.yaml +++ b/astra.yaml @@ -1351,7 +1351,7 @@ analyses: psf_epoch_averaging, galaxy_pixel_weights, psf_likelihood_noise, megacam_ccd_flip, defect_fill, blend_handling, defect_weighting, epoch_masked_fraction_cut, - central_defect_veto] + central_defect_veto, defect_flagging] decisions: ngmix_seed_mode: label: Per-object RNG seeded from sky position @@ -1546,6 +1546,30 @@ analyses: options: megapipe_flip: label: Flip CCDs < 18 and 36/37 (MegaPipe orientation) + defect_flagging: + label: Adaptive exposure defects on tile objects + rationale: >- + [HARDCODED] We retain every instrument bit by OR'ing exposure flag + pixels on the object's coadd footprint (seg == obj_id), flipped to + each epoch, across epochs kept after the fraction cut and central + veto. NGMIX_IMAFLAGS_ISO invents no bit for unflagged zero-weight + or off-tile pixels. NGMIX_N_EPOCH_INTERP_ISO counts kept epochs + with an interpolated defect on that footprint, including those + unflagged pixels. Missing, empty, or stamp-edge-touching footprints + give -1 in both columns because completeness is unknown. Never-fit + rows retain the existing count convention of 0. The Napoli + collaboration decision requests adaptive exposure flags propagated + to tile objects. Stamp-wide diagnostics count a distant defect + equally for small and large objects. The coadd footprint avoids + exposure detection matching; it is not reprojected, consistent with + the existing segmentation overlay. These columns impose no cut + and do not alter the fits. + default: footprint_or_and_interp_epochs + options: + footprint_or_and_interp_epochs: + label: Footprint IMAFLAGS_ISO OR plus interpolated-epoch count + stamp_wide_only: + label: Stamp-wide interpolated count and minimum distances only defect_fill: label: Image content of defect pixels before metacal rationale: >- diff --git a/universes/committed.yaml b/universes/committed.yaml index 79f62f356..13b5e3503 100644 --- a/universes/committed.yaml +++ b/universes/committed.yaml @@ -56,6 +56,7 @@ analyses: psf_likelihood_noise: psf_noise_1em5 megacam_ccd_flip: megapipe_flip defect_fill: interpolate + defect_flagging: footprint_or_and_interp_epochs blend_handling: uberseg defect_weighting: fourfold_zero epoch_masked_fraction_cut: ten_percent From 6edc7bea7b93df8e85ea53862b663e72b8a12ed2 Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 08:51:46 -0400 Subject: [PATCH 04/10] Cover footprint sentinels, safe casts, and ngmix FITS round trips Co-Authored-By: Claude Opus 5.5 --- tests/module/test_ngmix.py | 3 ++ tests/module/test_ngmix_defect_fill.py | 2 + tests/module/test_ngmix_defect_flags.py | 40 ++++++++++++++++++- tests/module/test_psf_averaging_properties.py | 6 ++- tests/module/test_psf_grammar_properties.py | 9 ++++- 5 files changed, 55 insertions(+), 5 deletions(-) diff --git a/tests/module/test_ngmix.py b/tests/module/test_ngmix.py index 37b3e58b5..c83bba8c2 100644 --- a/tests/module/test_ngmix.py +++ b/tests/module/test_ngmix.py @@ -332,6 +332,7 @@ def _fake_metacal_result(T, T_err, T_psf, T_psf_err): "mcal_types_fail": 0, "neighbour_flag": 0, "n_epoch_interp": 0, + "imaflags_iso": 0, "n_epoch_interp_iso": 0, "min_dist_interp": -1.0, "min_dist_noisefill": -1.0, # original image PSF (psfex/mccd) family @@ -638,6 +639,7 @@ def test_process_counts_flagged_fits_across_batches(tmp_path, monkeypatch, flags epoch_cuts=Counter(considered=1), epoch_failures=[], jacobs=[galsim.JacobianWCS(.186, 0., 0., .186)], defect_diagnostics=lambda: dict( + imaflags_iso=0, n_epoch_interp_iso=0, n_epoch_interp=0, min_dist_interp=-1.0, min_dist_noisefill=-1.0, ), ) @@ -1334,6 +1336,7 @@ def test_process_logs_and_counts_failed_epochs(tmp_path, monkeypatch): epoch_failures=[], jacobs=[galsim.JacobianWCS(.186, 0., 0., .186)] * 2, defect_diagnostics=lambda: dict( + imaflags_iso=0, n_epoch_interp_iso=0, n_epoch_interp=0, min_dist_interp=-1.0, min_dist_noisefill=-1.0, ), ) diff --git a/tests/module/test_ngmix_defect_fill.py b/tests/module/test_ngmix_defect_fill.py index 96de7e44c..b24bf971f 100644 --- a/tests/module/test_ngmix_defect_fill.py +++ b/tests/module/test_ngmix_defect_fill.py @@ -764,6 +764,7 @@ def test_defect_diagnostics_summarise_the_kept_epochs(): ) assert stamp.defect_diagnostics() == { "n_epoch_interp": 2, + "imaflags_iso": -1, "n_epoch_interp_iso": -1, "min_dist_interp": float(EPOCH_INTERPOLATED_DEFECT_RADIUS), "min_dist_noisefill": float(EPOCH_CENTRAL_DEFECT_RADIUS), } @@ -774,6 +775,7 @@ def test_defect_diagnostics_summarise_the_kept_epochs(): psf_obj=psf_obj, gal_obj=gal_obj, ) assert stamp.defect_diagnostics() == { + "imaflags_iso": -1, "n_epoch_interp_iso": -1, "n_epoch_interp": 0, "min_dist_interp": -1.0, "min_dist_noisefill": -1.0, } diff --git a/tests/module/test_ngmix_defect_flags.py b/tests/module/test_ngmix_defect_flags.py index bf11b2a8b..53118cf5a 100644 --- a/tests/module/test_ngmix_defect_flags.py +++ b/tests/module/test_ngmix_defect_flags.py @@ -15,7 +15,8 @@ def _run(kind="inside", seg_kind="valid", megacam_flip=True): if seg_kind == "clipped": seg[25, 36:] = 1 epochs = {} - for i, (ccd, bit) in enumerate(((10, 4), (20, 1024), (21, 0))): + ccds = (10, 20, 21) if megacam_flip else (20, 21, 22) + for i, (ccd, bit) in enumerate(zip(ccds, (4, 1024, 0))): flag = np.zeros(seg.shape, dtype=np.float32) weight = np.ones(seg.shape) pixel = (25, 35) if kind != "outside" else (35, 35) @@ -42,7 +43,6 @@ def _run(kind="inside", seg_kind="valid", megacam_flip=True): tile.seg = None if seg_kind == "missing" else seg[None] return prepare_postage_stamps( vignet, 1, 0, tile, bkg_sub=False, psf_obj=psf, gal_obj=gal, - megacam_flip=megacam_flip, ) @@ -71,3 +71,39 @@ def test_unknown_or_incomplete_footprint_is_not_clean(seg_kind): diagnostics = _run(seg_kind=seg_kind).defect_diagnostics() assert diagnostics["imaflags_iso"] == -1 assert diagnostics["n_epoch_interp_iso"] == -1 + + +@pytest.mark.parametrize("value", [np.nan, np.inf, -1, 1.5, float(2**63)]) +def test_non_integral_or_unsafe_flags_raise(value): + """Do not silently truncate float flags or overflow the integer OR.""" + from shapepipe.modules.ngmix_package.ngmix import Postage_stamp + + stamp = Postage_stamp() + seg = np.zeros((5, 5), dtype=int) + seg[2, 2] = 7 + flag = np.zeros(seg.shape) + flag[2, 2] = value + with pytest.raises(ValueError, match="footprint flags"): + stamp.record_footprint_defects(flag, flag != 0, seg, 7) + + +def test_footprint_columns_survive_ngmix_fits(tmp_path): + """Every shear extension keeps integer diagnostics, including -1.""" + from astropy.io import fits + from tests.module.test_ngmix import _fake_metacal_result + + inst = object.__new__(Ngmix) + inst._zero_point = 30.0 + inst._output_dir = str(tmp_path) + inst._file_number_string = "-001-001" + rows = [] + for i, (bits, count) in enumerate(((1028, 2), (-1, -1)), start=1): + row = _fake_metacal_result(.18, .02, .09, .001) + row.update(obj_id=i, imaflags_iso=bits, n_epoch_interp_iso=count) + rows.append(row) + inst.save_results(inst.compile_results(rows)) + with fits.open(inst.get_output_path(str(tmp_path))) as hdus: + for hdu in hdus[1:]: + assert hdu.data["imaflags_iso"].tolist() == [1028, -1] + assert hdu.data["n_epoch_interp_iso"].tolist() == [2, -1] + assert hdu.data["imaflags_iso"].dtype.kind == "i" diff --git a/tests/module/test_psf_averaging_properties.py b/tests/module/test_psf_averaging_properties.py index d6799bc29..49895b3e4 100644 --- a/tests/module/test_psf_averaging_properties.py +++ b/tests/module/test_psf_averaging_properties.py @@ -260,6 +260,8 @@ def test_all_epochs_failed_raises_zero_division(flagged_specs): "NGMIX_NEIGHBOUR_FLAG": 0.0, "NGMIX_N_EPOCH_FAILED": 0.0, "NGMIX_N_EPOCH_INTERP": 0.0, + "NGMIX_IMAFLAGS_ISO": 0, + "NGMIX_N_EPOCH_INTERP_ISO": 0, "NGMIX_MIN_DIST_INTERP": -1.0, "NGMIX_MIN_DIST_NOISEFILL": -1.0, } @@ -270,6 +272,7 @@ def test_all_epochs_failed_raises_zero_division(flagged_specs): "id", "n_epoch_model", "n_epoch_failed", "mcal_types_fail", "neighbour_flag", "n_epoch_interp", "min_dist_interp", "min_dist_noisefill", "nfev_fit", + "imaflags_iso", "n_epoch_interp_iso", "g1", "g1_err", "g2", "g2_err", "T", "T_err", "flux", "flux_err", "s2n", "mag", "mag_err", "flags", "mcal_flags", "g1_psf_orig", "g2_psf_orig", "g1_err_psf_orig", "g2_err_psf_orig", @@ -280,7 +283,8 @@ def test_all_epochs_failed_raises_zero_division(flagged_specs): _INT_KEYS = { "id", "n_epoch_model", "n_epoch_failed", "mcal_types_fail", "neighbour_flag", - "n_epoch_interp", "nfev_fit", "flags", "mcal_flags" + "n_epoch_interp", "nfev_fit", "flags", "mcal_flags", + "imaflags_iso", "n_epoch_interp_iso" } _SHEAR_EXTS = ["1M", "1P", "2M", "2P", "NOSHEAR"] diff --git a/tests/module/test_psf_grammar_properties.py b/tests/module/test_psf_grammar_properties.py index 4a6442f13..d35f6a462 100644 --- a/tests/module/test_psf_grammar_properties.py +++ b/tests/module/test_psf_grammar_properties.py @@ -67,11 +67,13 @@ def info(self, *_args, **_kwargs): "id", "n_epoch_model", "n_epoch_failed", "mcal_types_fail", "neighbour_flag", "n_epoch_interp", "nfev_fit", "flags", "mcal_flags", + "imaflags_iso", "n_epoch_interp_iso", } NGMIX_KEYS = [ "id", "n_epoch_model", "n_epoch_failed", "mcal_types_fail", "neighbour_flag", "n_epoch_interp", "min_dist_interp", "min_dist_noisefill", "nfev_fit", + "imaflags_iso", "n_epoch_interp_iso", "g1", "g1_err", "g2", "g2_err", "T", "T_err", "flux", "flux_err", "s2n", "mag", "mag_err", @@ -113,6 +115,7 @@ def _base_row(obj_id): return { "id": obj_id, "n_epoch_model": 3, "n_epoch_failed": 1, "mcal_types_fail": 0, "neighbour_flag": 1, + "imaflags_iso": 1028, "n_epoch_interp_iso": 1, "n_epoch_interp": 2, "min_dist_interp": 7.5, "min_dist_noisefill": 14.25, "nfev_fit": 7, "g1": 0.10, "g1_err": 0.011, "g2": -0.20, "g2_err": 0.022, @@ -184,7 +187,8 @@ def _run_save_ngmix(ngmix_path, obj_ids, cat_size_target=None): GRAMMAR_RE = re.compile( rf"^NGMIXm?_(?:{_COMPONENT})(?:_ERR)?(?:_(?:{_OBJECT}))?_(?:{_SHEAR})$" rf"|^NGMIXm?_(?:MCAL_FLAGS|MCAL_TYPES_FAIL|N_EPOCH|NEIGHBOUR_FLAG" - rf"|N_EPOCH_FAILED|N_EPOCH_INTERP|MIN_DIST_INTERP|MIN_DIST_NOISEFILL)$" + rf"|N_EPOCH_FAILED|N_EPOCH_INTERP|MIN_DIST_INTERP|MIN_DIST_NOISEFILL" + rf"|IMAFLAGS_ISO|N_EPOCH_INTERP_ISO)$" ) # The shipped final-catalogue param files, two levels up from tests/module/. @@ -424,7 +428,8 @@ def test_param_file_ngmix_tokens_are_producible(param_path, obj_ids): r"NGMIXm?_(?:G1|G2|T|SNR|FLUX|MAG|FLAGS)(?:_ERR)?" r"(?:_(?:PSF_ORIG|PSF_RECONV))?_(?:NOSHEAR|1P|1M|2P|2M)" r"|NGMIXm?_(?:MCAL_FLAGS|MCAL_TYPES_FAIL|N_EPOCH|NEIGHBOUR_FLAG" - r"|N_EPOCH_FAILED|N_EPOCH_INTERP|MIN_DIST_INTERP|MIN_DIST_NOISEFILL)" + r"|N_EPOCH_FAILED|N_EPOCH_INTERP|MIN_DIST_INTERP|MIN_DIST_NOISEFILL" + r"|IMAFLAGS_ISO|N_EPOCH_INTERP_ISO)" # HSM: g-type, explicit PSF/STAR object, singular FLAG; the multi-epoch # sink in make_cat._save_psf_data appends a bare epoch index. r"|HSM_(?:G1|G2|T|M4_1|M4_2|RHO4)_(?:PSF|STAR)(?:_\d+)?" From c5fcdc3ab008230a1063d458859355d9efc3376e Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 08:54:01 -0400 Subject: [PATCH 05/10] Add read-only no-fit tile defect validation and document pending survey check Co-Authored-By: Claude Opus 5.5 --- .../masked_pixels/tile_defect_flags.md | 64 ++++++ .../masked_pixels/tile_defect_flags.py | 201 ++++++++++++++++++ 2 files changed, 265 insertions(+) create mode 100644 scripts/validation/masked_pixels/tile_defect_flags.md create mode 100644 scripts/validation/masked_pixels/tile_defect_flags.py diff --git a/scripts/validation/masked_pixels/tile_defect_flags.md b/scripts/validation/masked_pixels/tile_defect_flags.md new file mode 100644 index 000000000..02046278d --- /dev/null +++ b/scripts/validation/masked_pixels/tile_defect_flags.md @@ -0,0 +1,64 @@ +# Tile-object defect flags + +We propagate exposure instrument bits onto each tile object's coadd segmentation footprint, using only epochs that pass the masked-fraction cut and central veto. +The implementation does not change fits or impose a catalogue cut. + +- `NGMIX_IMAFLAGS_ISO` is the bitwise OR of exposure flag pixels on `tile_seg == obj_id`, OR'ed across kept epochs, with every instrument bit preserved. +- `NGMIX_N_EPOCH_INTERP_ISO` counts kept epochs containing at least one interpolated defect pixel on that footprint. +- Both are `-1` when segmentation is missing, the footprint is empty, or it touches a stamp edge. + We cannot establish completeness in those cases. + Never-fit final-catalogue rows retain the existing count convention of `0`. +- A zero-weight/off-tile pixel with flag `0` contributes no instrument bit, but contributes to the epoch count if interpolated. + Clean quarter-turn copies made by defect weighting do not count as defects. + +The coadd segmentation is MegaCam-flipped to each exposure's stamp orientation, without reprojection, matching the existing segmentation overlay. +The columns supplement the stamp-wide diagnostics rather than replacing the tile SExtractor's `IMAFLAGS_ISO`. + +## Fixture validation + +The synthetic path runs `prepare_postage_stamps` with multiple epochs and both flipped and unflipped CCDs. +It checks inside/outside-footprint defects, clean objects, unflagged zero-weight defects, rejected epochs, missing/empty/clipped footprints, safe float-to-integer conversion, integer FITS round trips in every shear extension, and final-catalogue values and never-fit sentinels. + +| Mutation | Observed failure | Result | +| --- | --- | --- | +| Reduce/count over the entire stamp | Outside-footprint case becomes `1028 / 2`, expected `0 / 0` | 2 assertion failures | +| Overwrite the OR on each epoch | Final clean epoch erases `1028`, leaving `0` | 2 assertion failures | +| Accumulate before either epoch cut | Rejected epochs contribute bit `64` and interpolation counts | 8 assertion failures | +| Drop the make_cat write mapping | Final-catalogue value stays `0`, expected `1028` | 1 assertion failure | + +All four mutations were restored before running the suite. +The runtime SIF lacks pytest; validation uses `shapepipe-develop-240b37e4.sif` on Slurm allocation `24917684` under `def-mjhudson`. +The params-pin test runs inside that image with the existing `/project/6001537/cdaley/snakemake-env/lib/python3.12/site-packages` prepended and pytest plugin autoload disabled, because the image lacks Snakemake. +Regeneration leaves the pin unchanged. +`completeness.py` counts files, not columns; this change adds no product file. + +## Flag dtype + +`split_exp.py` casts flag images to `int16`, and `get_stamps` preserves the input dtype. +On real exposure `2254826p.flag.fits.fz`, CCD HDU 1 is `uint8`, with values `0, 1, 3, 8, 10, 11`; applying the pipeline's split cast and stamp extraction produces `int16` stamps. +This checks the extraction path, **not a retained production SQLite flag vignet**, which was not available. +The fixture also exercises integral `float32` flag stamps. +The new reduction validates finite, nonnegative, integral, int64-representable values before casting to `int64`, so floats cannot silently truncate or overflow. + +## Real-tile validation remains pending + +**No real-tile fractions or spatial plot have been measured for this change.** +The retained final catalogues and `nflagged` maps are available, but a complete tile with its image, flag, weight, PSF and segmentation stores has not been located. +For example, the `smk-ab925-us` tile `181.308` manifest points to `/local/scratch/sp-181.308`, and the shared tile directory contains `cleaned.json` and only the exposure-finder output. +`fir-b0001.yaml` likewise points to node-local inputs/stores and requests `clean_tiles: true`. +The `ngmix-fixtest` inputs for `196.305` are broken symlinks to removed `p3-first10` stores. +A bounded search of the shared products, scratch run directories and test runs did not locate a complete alternative; this is not proof that no retained tile exists elsewhere. + +`tile_defect_flags.py` opens stores read-only and calls the real `prepare_postage_stamps` path without fitting. +Its module docstring gives the full command. +It writes per-object arrays, counts, flag-bit breakdowns, explicit denominators and two Poisson reference estimates under the requested scratch output prefix. +With `--nflagged` and `--plot`, it overlays flagged objects on the exposure defect map. +A synthetic one-object end-to-end smoke test produced bit `4` and one interpolated epoch; those are fixture results, not survey measurements. + +For the real check, compare the measured instrument-bit fraction with the mean of `1 - exp(-0.0028 * A * N)` over known footprints, using each object's area and kept-epoch count. +That expression ignores spatial correlations and the central veto: every kept interpolated defect is at least 7 px from the stamp centre, so compact footprints can be protected completely. +The script also reports a reference using only footprint area outside that radius; even this is not a prediction calibrated for the epoch-selection correlations or invalid-RMS/noise-filled defects. +Explain the discrepancy rather than treating agreement with an independent-pixel approximation as a requirement. + +Back-filling would require retained or regenerated stores, but would not require shear refits. +No back-fill has been attempted. diff --git a/scripts/validation/masked_pixels/tile_defect_flags.py b/scripts/validation/masked_pixels/tile_defect_flags.py new file mode 100644 index 000000000..e82bce3fa --- /dev/null +++ b/scripts/validation/masked_pixels/tile_defect_flags.py @@ -0,0 +1,201 @@ +"""Measure footprint defect flags without fitting, using retained tile stores. + +Run inside the ShapePipe container on a compute allocation, with PYTHONPATH=src:: + + python scripts/validation/masked_pixels/tile_defect_flags.py \ + --sexcat SEXCAT --segcat SEG_VIGNET_CAT --image IMAGE_SQLITE \ + --flag FLAG_SQLITE --weight WEIGHT_SQLITE --psf GALAXY_PSF_SQLITE \ + --headers HEADERS_SQLITE --background BACKGROUND_SQLITE \ + --rms BACKGROUND_RMS_SQLITE --blend-handling uberseg \ + --nflagged NFLAGGED.hsp --out /scratch/cdaley/defectflag-check/TILE \ + --plot scripts/validation/masked_pixels/tile_defect_flags-TILE.png + +All stores open read-only; nothing is written beneath an input directory. +Omit --segcat only if SEXCAT itself carries SEG_VIGNET. Omit --background +only for already background-subtracted images; omit --rms only when the run +used stamp-estimated RMS. Pass the production blend/fill choices explicitly. +The output contains every tile row, including unavailable PSFs, zero-kept- +epoch objects and unknown footprints, distinguished by status and -1. +The Poisson expectation uses 0.0028 as the exposure DATASEC pixel flag rate +from #887. It is an unconditioned reference, not a calibration: spatially +correlated defects and the central veto violate its independence assumption. +""" + +import argparse +from collections import Counter +from contextlib import ExitStack +import json +from pathlib import Path +from types import SimpleNamespace + +import fitsio +import matplotlib +matplotlib.use("Agg") +import matplotlib.pyplot as plt +import numpy as np +from sqlitedict import SqliteDict + +from shapepipe.modules.ngmix_package.ngmix import ( + EPOCH_INTERPOLATED_DEFECT_RADIUS, Tile_cat, centre_distance, + prepare_postage_stamps, +) + + +def measure(args): + """Read production stamps and invoke the pipeline's epoch preparation.""" + tile = Tile_cat(args.sexcat, seg_cat_path=args.segcat) + if tile.seg is None: + with fitsio.FITS(args.sexcat) as cat: + objects = next(h for h in cat if h.get_extname() == "LDAC_OBJECTS") + if "SEG_VIGNET" in objects.get_colnames(): + tile.seg = objects.read_column("SEG_VIGNET") + if tile.seg is None: + raise ValueError("A real-tile validation requires segmentation stamps") + n = len(tile.obj_id) + data = dict(number=tile.obj_id, ra=tile.ra, dec=tile.dec) + for key in ("imaflags_iso", "n_epoch_interp_iso", "n_epoch_interp", + "n_kept", "area", "area_outside_veto"): + data[key] = np.full(n, -1, dtype=np.int64) + status, dtypes = Counter(), Counter() + with ExitStack() as stack: + def store(path): + if path is None: + return None + if not Path(path).is_file(): + raise FileNotFoundError(path) + db = SqliteDict(path, flag="r") + stack.callback(db.close) + return db + + vignet = SimpleNamespace( + gal_vign_cat=store(args.image), flag_vign_cat=store(args.flag), + weight_vign_cat=store(args.weight), psf_vign_cat=store(args.psf), + f_wcs_file=store(args.headers), bkg_vign_cat=store(args.background), + bkg_rms_vign_cat=store(args.rms), + ) + radius = centre_distance(tile.seg[0].shape) + for i, obj_id in enumerate(tile.obj_id): + footprint = tile.seg[i] == obj_id + data["area"][i] = footprint.sum() + data["area_outside_veto"][i] = ( + footprint & (radius >= EPOCH_INTERPOLATED_DEFECT_RADIUS) + ).sum() + key = str(obj_id) + if key not in vignet.psf_vign_cat or key not in vignet.gal_vign_cat: + status["missing_store_row"] += 1 + continue + psf, gal = vignet.psf_vign_cat[key], vignet.gal_vign_cat[key] + if isinstance(psf, str) or isinstance(gal, str): + status["empty_store_row"] += 1 + continue + stamp = prepare_postage_stamps( + vignet, obj_id, i, tile, bkg_sub=args.background is not None, + psf_obj=psf, gal_obj=gal, blend_handling=args.blend_handling, + defect_fill=args.defect_fill, + ) + data["n_kept"][i] = len(stamp.flags) + for flag in stamp.flags: + dtypes[str(flag.dtype)] += 1 + if not stamp.flags: + status["no_kept_epoch"] += 1 + continue + for name, value in stamp.defect_diagnostics().items(): + if name in data: + data[name][i] = value + status["unknown_footprint" if stamp.imaflags_iso < 0 else "known"] += 1 + if i % 1000 == 0: + print(f"{i}/{n}", flush=True) + valid = data["imaflags_iso"] >= 0 + kept = data["n_kept"] > 0 + interp = data["n_epoch_interp_iso"] > 0 + flagged = data["imaflags_iso"] > 0 + wide = data["n_epoch_interp"] > 0 + summary = dict( + inputs=vars(args), n_objects=n, status=dict(status), + kept_flag_vignet_dtypes=dict(dtypes), + n_interp_iso=int(interp.sum()), n_imaflags_iso=int(flagged.sum()), + n_stamp_interp=int(wide.sum()), + fractions_all_rows=dict(interp_iso=float(interp.mean()), + imaflags_iso=float(flagged.mean()), + stamp_interp=float(wide.mean())), + bits={str(1 << bit): int(((data["imaflags_iso"][valid] & (1 << bit)) != 0).sum()) + for bit in range(63) + if ((data["imaflags_iso"][valid] & (1 << bit)) != 0).any()}, + ) + if valid.any(): + summary["fractions_known_footprint"] = dict( + interp_iso=float(interp[valid].mean()), + imaflags_iso=float(flagged[valid].mean()), + ) + for area in ("area", "area_outside_veto"): + summary[f"poisson_{area}"] = float(np.mean( + -np.expm1(-0.0028 * data[area][valid] * data["n_kept"][valid]) + )) + summary["mean_area_known"] = float(data["area"][valid].mean()) + summary["mean_kept_epochs_known"] = float(data["n_kept"][valid].mean()) + if kept.any(): + summary["fraction_stamp_interp_kept"] = float(wide[kept].mean()) + out = Path(args.out) + out.parent.mkdir(parents=True, exist_ok=True) + np.savez_compressed(str(out) + ".npz", **data) + Path(str(out) + ".json").write_text(json.dumps(summary, indent=2) + "\n") + if args.plot: + plot(args, data, summary) + print(json.dumps(summary, indent=2)) + + +def plot(args, data, summary): + """Overlay interpolated-footprint objects on the retained nflagged map.""" + import healsparse + from figstyle import INK2, ORANGE, PAPER + + if not args.nflagged: + raise ValueError("--plot requires the exposure --nflagged map") + sky = healsparse.HealSparseMap.read(args.nflagged) + ra, dec = data["ra"], data["dec"] + # Unwrap a tile crossing RA=0 before laying out its sky grid. + centre = float(np.rad2deg(np.angle(np.mean(np.exp(1j * np.deg2rad(ra)))))) + x = (ra - centre + 180) % 360 - 180 + extent = [x.min(), x.max(), dec.min(), dec.max()] + xx, yy = np.meshgrid(np.linspace(*extent[:2], 700), + np.linspace(*extent[2:], 700)) + image = sky.get_values_pos((xx + centre) % 360, yy, lonlat=True).astype(float) + image[image < 0] = np.nan + at_objects = sky.get_values_pos(ra, dec, lonlat=True) + selected = data["n_epoch_interp_iso"] > 0 + summary["map_nflagged_positive_at_interp_iso"] = int( + ((at_objects > 0) & selected).sum() + ) + fig, ax = plt.subplots(figsize=(9, 8), facecolor=PAPER) + ax.set_facecolor(PAPER) + im = ax.imshow(image, origin="lower", extent=extent, cmap="Blues", + interpolation="nearest", aspect="auto") + ax.scatter(x[~selected], dec[~selected], s=0.6, c=INK2, alpha=.12, + linewidths=0, label="Other tile detections") + ax.scatter(x[selected], dec[selected], s=13, c=ORANGE, linewidths=.3, + edgecolors="white", label=f"Interpolated defect on footprint ({selected.sum():,})") + ax.set(xlabel=f"RA offset from {centre % 360:.3f}° [deg]", + ylabel="Declination [deg]", title="Kept-epoch defects on tile-object footprints") + ax.invert_xaxis() + ax.legend(loc="upper right", markerscale=2) + fig.colorbar(im, ax=ax, label="Flagged exposures per sky pixel (nflagged)") + fig.tight_layout() + Path(args.plot).parent.mkdir(parents=True, exist_ok=True) + fig.savefig(args.plot, dpi=180) + plt.close(fig) + Path(args.out + ".json").write_text(json.dumps(summary, indent=2) + "\n") + + +def main(): + parser = argparse.ArgumentParser(description=__doc__) + for name in ("sexcat", "image", "flag", "weight", "psf", "headers", "out"): + parser.add_argument("--" + name, required=True) + for name in ("segcat", "background", "rms", "nflagged", "plot"): + parser.add_argument("--" + name) + parser.add_argument("--blend-handling", choices=("noisefill", "uberseg"), required=True) + parser.add_argument("--defect-fill", choices=("interpolate", "noise"), default="interpolate") + measure(parser.parse_args()) + + +if __name__ == "__main__": + main() From 33a26b11b6e9e45c6c4e572f8af886819d228b25 Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 09:32:35 -0400 Subject: [PATCH 06/10] Guard committed SEG_VIGNET delivery to adaptive defect flags under every handling Co-Authored-By: Claude Opus 5.5 --- tests/workflow/test_defect_flags.py | 50 +++++++++++++++++++++++++++++ 1 file changed, 50 insertions(+) create mode 100644 tests/workflow/test_defect_flags.py diff --git a/tests/workflow/test_defect_flags.py b/tests/workflow/test_defect_flags.py new file mode 100644 index 000000000..72334e49a --- /dev/null +++ b/tests/workflow/test_defect_flags.py @@ -0,0 +1,50 @@ +"""The committed detection config supplies the adaptive-flag footprint.""" + +from pathlib import Path +import re + +import numpy as np +import pytest + +from tests.workflow.harness import Campaign +from tests.workflow.test_dag import _committed_value + + +@pytest.mark.parametrize("input_type", ["data", "image_sims"]) +@pytest.mark.parametrize("blend", [None, "noisefill", "uberseg"]) +def test_committed_detection_supplies_defect_flag_footprint( + tmp_path, resolve_dag, monkeypatch, input_type, blend): + """Every handling enables SEG_VIGNET and Tile_cat consumes that column. + + This guards the shipped path, rather than assigning tile.seg by hand: + removing the workflow export, disabling the committed ini option, or + dropping the loader's SEG_VIGNET read must not silently turn an ordinary + clean object into an unknown-footprint (-1) object. The SExtractor runner + column-writing behavior is covered by test_sextractor_seg_vignet. + """ + from shapepipe.modules.ngmix_package.ngmix import Postage_stamp, Tile_cat + from tests.module.test_ngmix_uberseg import _write_tile_cat + + campaign = Campaign(tmp_path / "campaign", input_type, "psfex") + if blend is not None: + campaign.config["blend_handling"] = blend + campaign.write_config() + with resolve_dag(campaign) as dag: + job = next(job for job in dag.jobs if job.rule.name == "tile_detect") + exports = dict(re.findall(r"^export (SP_SEG_VIGNET)='([^']*)'$", + job.params.pre, re.MULTILINE)) + config_dir = (Path(__file__).parents[2] / "workflow" / "config" + / ("cfis" if input_type == "data" else "cfis_image_sims")) + enabled = _committed_value(job.shellcmd, config_dir, "SEG_VIGNET", + exports, monkeypatch) + assert enabled == "True", "committed detection disables SEG_VIGNET" + seg = np.zeros((3, 5, 5), dtype=np.int32) + seg[:, 2, 2] = [4, 5, 9] + path = _write_tile_cat(tmp_path / "tile.fits", seg) + tile = Tile_cat(str(path)) + stamp = Postage_stamp() + stamp.record_footprint_defects( + np.zeros((5, 5), dtype=np.int16), np.zeros((5, 5), dtype=bool), + None if tile.seg is None else tile.seg[0], tile.obj_id[0], + ) + assert stamp.imaflags_iso == stamp.n_epoch_interp_iso == 0 From 869c072c84282b6749bd29b0b41c37d73d6bb18f Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 09:32:36 -0400 Subject: [PATCH 07/10] Read the stacked tile SEG_VIGNET and report the blocked private regeneration plan Co-Authored-By: Claude Opus 5.5 --- .../masked_pixels/tile_defect_flags.md | 23 ++++++++++++-- .../masked_pixels/tile_defect_flags.py | 30 ++++++++++++------- 2 files changed, 39 insertions(+), 14 deletions(-) diff --git a/scripts/validation/masked_pixels/tile_defect_flags.md b/scripts/validation/masked_pixels/tile_defect_flags.md index 02046278d..4c5c05e43 100644 --- a/scripts/validation/masked_pixels/tile_defect_flags.md +++ b/scripts/validation/masked_pixels/tile_defect_flags.md @@ -2,6 +2,8 @@ We propagate exposure instrument bits onto each tile object's coadd segmentation footprint, using only epochs that pass the masked-fraction cut and central veto. The implementation does not change fits or impose a catalogue cut. +This change stacks on #925 (`feat/uberseg-seg-vignet`), which puts the int32 `SEG_VIGNET` column in the committed tile catalogue under every blend handling and makes `Tile_cat` read it. +Stacking only on #922 would leave the committed workflow without the footprint input and produce `-1` for fitted objects. - `NGMIX_IMAFLAGS_ISO` is the bitwise OR of exposure flag pixels on `tile_seg == obj_id`, OR'ed across kept epochs, with every instrument bit preserved. - `NGMIX_N_EPOCH_INTERP_ISO` counts kept epochs containing at least one interpolated defect pixel on that footprint. @@ -25,9 +27,15 @@ It checks inside/outside-footprint defects, clean objects, unflagged zero-weight | Overwrite the OR on each epoch | Final clean epoch erases `1028`, leaving `0` | 2 assertion failures | | Accumulate before either epoch cut | Rejected epochs contribute bit `64` and interpolation counts | 8 assertion failures | | Drop the make_cat write mapping | Final-catalogue value stays `0`, expected `1028` | 1 assertion failure | +| Remove the workflow's SEG_VIGNET export | Committed detection disables the footprint | 6 assertion failures | +| Disable SEG_VIGNET in the committed ini | Committed detection disables the footprint | 6 assertion failures | +| Drop Tile_cat's SEG_VIGNET read | Ordinary clean objects receive `-1` instead of `0` | 6 assertion failures | -All four mutations were restored before running the suite. -The runtime SIF lacks pytest; validation uses `shapepipe-develop-240b37e4.sif` on Slurm allocation `24917684` under `def-mjhudson`. +All seven mutations were restored before running the suite. +`test_committed_detection_supplies_defect_flag_footprint` resolves the actual workflow and ini for data and image simulations under default, noisefill and uberseg handling, then checks the `SEG_VIGNET` load and a clean object's zero-valued diagnostics. +The runtime SIF lacks pytest; validation uses `shapepipe-develop-240b37e4.sif` on Slurm allocations under `def-mjhudson`. +The rebased branch's ngmix, make_cat, unit and PSF property suite passes 830 tests with one skipped; the workflow/DAG/params-pin suite passes 31 tests with one skipped. +The seven mutations were rerun on the rebased branch on allocation `24920131`, with source edits and pytest on the same compute node and bytecode caches cleared between mutations. The params-pin test runs inside that image with the existing `/project/6001537/cdaley/snakemake-env/lib/python3.12/site-packages` prepended and pytest plugin autoload disabled, because the image lacks Snakemake. Regeneration leaves the pin unchanged. `completeness.py` counts files, not columns; this change adds no product file. @@ -49,9 +57,18 @@ For example, the `smk-ab925-us` tile `181.308` manifest points to `/local/scratc The `ngmix-fixtest` inputs for `196.305` are broken symlinks to removed `p3-first10` stores. A bounded search of the shared products, scratch run directories and test runs did not locate a complete alternative; this is not proof that no retained tile exists elsewhere. +A private regeneration config is saved at `/scratch/cdaley/defectflag-check/run/config.yaml`, with separate work, products, index and tile-store directories under that prefix and both cleaners disabled. +It selects tile `181.308`, imports its retained PSF tar directories through read-only symlinks, and reads the shared raw mirrors. +The dry run targets only `manifests/tile_vignets.json`; its command and log are `run/dryrun.sh` and `run/dryrun.log` beneath the same prefix. +The plan contains **35 jobs: 9 exposure retrievals, 9 exposure splits, 9 exposure PSF chains, and 8 tile-level jobs**. +All nine source exposure output directories are absent, consistent with the recorded cleanup; the PSF tarballs do not include the split images, headers or SExtractor background maps needed for the vignets. +The run stopped at the dry run, as required when any exposure-level job is planned. +No split, detection, PSF fitting, vignette generation or ngmix job was executed. + `tile_defect_flags.py` opens stores read-only and calls the real `prepare_postage_stamps` path without fitting. Its module docstring gives the full command. -It writes per-object arrays, counts, flag-bit breakdowns, explicit denominators and two Poisson reference estimates under the requested scratch output prefix. +It reads `SEG_VIGNET` through the rebased branch's `Tile_cat`, without a separate segmentation path or a script-only fallback. +It writes per-object arrays, counts, flag-bit breakdowns, explicit denominators, the edge-touching and empty-footprint counts, the unknown-footprint fraction among objects with kept epochs, and two Poisson reference estimates under the requested scratch output prefix. With `--nflagged` and `--plot`, it overlays flagged objects on the exposure defect map. A synthetic one-object end-to-end smoke test produced bit `4` and one interpolated epoch; those are fixture results, not survey measurements. diff --git a/scripts/validation/masked_pixels/tile_defect_flags.py b/scripts/validation/masked_pixels/tile_defect_flags.py index e82bce3fa..375f130ee 100644 --- a/scripts/validation/masked_pixels/tile_defect_flags.py +++ b/scripts/validation/masked_pixels/tile_defect_flags.py @@ -3,7 +3,7 @@ Run inside the ShapePipe container on a compute allocation, with PYTHONPATH=src:: python scripts/validation/masked_pixels/tile_defect_flags.py \ - --sexcat SEXCAT --segcat SEG_VIGNET_CAT --image IMAGE_SQLITE \ + --sexcat SEXCAT --image IMAGE_SQLITE \ --flag FLAG_SQLITE --weight WEIGHT_SQLITE --psf GALAXY_PSF_SQLITE \ --headers HEADERS_SQLITE --background BACKGROUND_SQLITE \ --rms BACKGROUND_RMS_SQLITE --blend-handling uberseg \ @@ -11,8 +11,8 @@ --plot scripts/validation/masked_pixels/tile_defect_flags-TILE.png All stores open read-only; nothing is written beneath an input directory. -Omit --segcat only if SEXCAT itself carries SEG_VIGNET. Omit --background -only for already background-subtracted images; omit --rms only when the run +SEXCAT must carry the committed workflow's SEG_VIGNET column. Omit +--background only for already background-subtracted images; omit --rms when the run used stamp-estimated RMS. Pass the production blend/fill choices explicitly. The output contains every tile row, including unavailable PSFs, zero-kept- epoch objects and unknown footprints, distinguished by status and -1. @@ -28,7 +28,6 @@ from pathlib import Path from types import SimpleNamespace -import fitsio import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt @@ -43,12 +42,7 @@ def measure(args): """Read production stamps and invoke the pipeline's epoch preparation.""" - tile = Tile_cat(args.sexcat, seg_cat_path=args.segcat) - if tile.seg is None: - with fitsio.FITS(args.sexcat) as cat: - objects = next(h for h in cat if h.get_extname() == "LDAC_OBJECTS") - if "SEG_VIGNET" in objects.get_colnames(): - tile.seg = objects.read_column("SEG_VIGNET") + tile = Tile_cat(args.sexcat) if tile.seg is None: raise ValueError("A real-tile validation requires segmentation stamps") n = len(tile.obj_id) @@ -56,6 +50,7 @@ def measure(args): for key in ("imaflags_iso", "n_epoch_interp_iso", "n_epoch_interp", "n_kept", "area", "area_outside_veto"): data[key] = np.full(n, -1, dtype=np.int64) + data["footprint_clipped"] = np.zeros(n, dtype=bool) status, dtypes = Counter(), Counter() with ExitStack() as stack: def store(path): @@ -77,6 +72,9 @@ def store(path): for i, obj_id in enumerate(tile.obj_id): footprint = tile.seg[i] == obj_id data["area"][i] = footprint.sum() + data["footprint_clipped"][i] = ( + footprint[[0, -1], :].any() or footprint[:, [0, -1]].any() + ) data["area_outside_veto"][i] = ( footprint & (radius >= EPOCH_INTERPOLATED_DEFECT_RADIUS) ).sum() @@ -113,6 +111,15 @@ def store(path): summary = dict( inputs=vars(args), n_objects=n, status=dict(status), kept_flag_vignet_dtypes=dict(dtypes), + footprints=dict( + n_edge_touching=int(data["footprint_clipped"].sum()), + fraction_edge_touching=float(data["footprint_clipped"].mean()), + n_empty=int((data["area"] == 0).sum()), + n_edge_touching_with_kept_epoch=int( + (data["footprint_clipped"] & kept).sum() + ), + n_unknown_with_kept_epoch=int((~valid & kept).sum()), + ), n_interp_iso=int(interp.sum()), n_imaflags_iso=int(flagged.sum()), n_stamp_interp=int(wide.sum()), fractions_all_rows=dict(interp_iso=float(interp.mean()), @@ -135,6 +142,7 @@ def store(path): summary["mean_kept_epochs_known"] = float(data["n_kept"][valid].mean()) if kept.any(): summary["fraction_stamp_interp_kept"] = float(wide[kept].mean()) + summary["fraction_unknown_footprint_kept"] = float((~valid)[kept].mean()) out = Path(args.out) out.parent.mkdir(parents=True, exist_ok=True) np.savez_compressed(str(out) + ".npz", **data) @@ -190,7 +198,7 @@ def main(): parser = argparse.ArgumentParser(description=__doc__) for name in ("sexcat", "image", "flag", "weight", "psf", "headers", "out"): parser.add_argument("--" + name, required=True) - for name in ("segcat", "background", "rms", "nflagged", "plot"): + for name in ("background", "rms", "nflagged", "plot"): parser.add_argument("--" + name) parser.add_argument("--blend-handling", choices=("noisefill", "uberseg"), required=True) parser.add_argument("--defect-fill", choices=("interpolate", "noise"), default="interpolate") From 27009cc38c0239edf21d1ebab27f287068648e8b Mon Sep 17 00:00:00 2001 From: Cail Daley Date: Sat, 24 Oct 2026 11:53:09 -0400 Subject: [PATCH 08/10] Validate footprint defect flags on real tile 188.307 Measure real epoch preparation without fitting, record denominators and bit counts, and add a spatial overlay and conditional Poisson interpretation. Cache read-only stores and exposure WCS decoding for bounded validation runs. Co-Authored-By: Claude Opus 5.5 --- .../masked_pixels/tile_defect_flags.md | 132 ++++++++++++++---- .../masked_pixels/tile_defect_flags.py | 118 ++++++++++++++-- .../tile_defect_flags_188307.json | 79 +++++++++++ .../tile_defect_flags_188307.png | Bin 0 -> 687384 bytes 4 files changed, 287 insertions(+), 42 deletions(-) create mode 100644 scripts/validation/masked_pixels/tile_defect_flags_188307.json create mode 100644 scripts/validation/masked_pixels/tile_defect_flags_188307.png diff --git a/scripts/validation/masked_pixels/tile_defect_flags.md b/scripts/validation/masked_pixels/tile_defect_flags.md index 4c5c05e43..456cee404 100644 --- a/scripts/validation/masked_pixels/tile_defect_flags.md +++ b/scripts/validation/masked_pixels/tile_defect_flags.md @@ -44,38 +44,112 @@ Regeneration leaves the pin unchanged. `split_exp.py` casts flag images to `int16`, and `get_stamps` preserves the input dtype. On real exposure `2254826p.flag.fits.fz`, CCD HDU 1 is `uint8`, with values `0, 1, 3, 8, 10, 11`; applying the pipeline's split cast and stamp extraction produces `int16` stamps. -This checks the extraction path, **not a retained production SQLite flag vignet**, which was not available. +The real-tile check below additionally inspects all 112,734 kept flag stamps in newly regenerated SQLite stores: every stamp is `int16`. The fixture also exercises integral `float32` flag stamps. The new reduction validates finite, nonnegative, integral, int64-representable values before casting to `int64`, so floats cannot silently truncate or overflow. -## Real-tile validation remains pending - -**No real-tile fractions or spatial plot have been measured for this change.** -The retained final catalogues and `nflagged` maps are available, but a complete tile with its image, flag, weight, PSF and segmentation stores has not been located. -For example, the `smk-ab925-us` tile `181.308` manifest points to `/local/scratch/sp-181.308`, and the shared tile directory contains `cleaned.json` and only the exposure-finder output. -`fir-b0001.yaml` likewise points to node-local inputs/stores and requests `clean_tiles: true`. -The `ngmix-fixtest` inputs for `196.305` are broken symlinks to removed `p3-first10` stores. -A bounded search of the shared products, scratch run directories and test runs did not locate a complete alternative; this is not proof that no retained tile exists elsewhere. - -A private regeneration config is saved at `/scratch/cdaley/defectflag-check/run/config.yaml`, with separate work, products, index and tile-store directories under that prefix and both cleaners disabled. -It selects tile `181.308`, imports its retained PSF tar directories through read-only symlinks, and reads the shared raw mirrors. -The dry run targets only `manifests/tile_vignets.json`; its command and log are `run/dryrun.sh` and `run/dryrun.log` beneath the same prefix. -The plan contains **35 jobs: 9 exposure retrievals, 9 exposure splits, 9 exposure PSF chains, and 8 tile-level jobs**. -All nine source exposure output directories are absent, consistent with the recorded cleanup; the PSF tarballs do not include the split images, headers or SExtractor background maps needed for the vignets. -The run stopped at the dry run, as required when any exposure-level job is planned. -No split, detection, PSF fitting, vignette generation or ngmix job was executed. - -`tile_defect_flags.py` opens stores read-only and calls the real `prepare_postage_stamps` path without fitting. -Its module docstring gives the full command. -It reads `SEG_VIGNET` through the rebased branch's `Tile_cat`, without a separate segmentation path or a script-only fallback. -It writes per-object arrays, counts, flag-bit breakdowns, explicit denominators, the edge-touching and empty-footprint counts, the unknown-footprint fraction among objects with kept epochs, and two Poisson reference estimates under the requested scratch output prefix. -With `--nflagged` and `--plot`, it overlays flagged objects on the exposure defect map. -A synthetic one-object end-to-end smoke test produced bit `4` and one interpolated epoch; those are fixture results, not survey measurements. - -For the real check, compare the measured instrument-bit fraction with the mean of `1 - exp(-0.0028 * A * N)` over known footprints, using each object's area and kept-epoch count. -That expression ignores spatial correlations and the central veto: every kept interpolated defect is at least 7 px from the stamp centre, so compact footprints can be protected completely. -The script also reports a reference using only footprint area outside that radius; even this is not a prediction calibrated for the epoch-selection correlations or invalid-RMS/noise-filled defects. -Explain the discrepancy rather than treating agreement with an independent-pixel approximation as a requirement. +## Real-tile validation: 188.307 + +We regenerated inputs for tile **188.307** from the same campaign as 181.308. +It needs seven exposures rather than nine: `2608415`, `2254827`, `2608414`, `2254828`, `2254835`, `2114047`, and `2114046`. +All 21 raw exposure image/weight/flag files, the tile image and weight, and the external catalogue were already local. +Retrieval used symlinks; nothing was downloaded. +The committed detection exports `SEG_VIGNET`, and the check reads it through `Tile_cat`, without a script-only segmentation fallback. +We ran the real `prepare_postage_stamps` path with `blend_handling=uberseg` and `defect_fill=interpolate`, **without running galaxy fits, tile_ngmix, or make_cat**. +This diagnostic path stops before applying the fit-time defect-weighting scheme. +These are detections presented to epoch preparation, not a shear-selected sample. + +### Measured fractions + +Of **36,051 objects**, 374 have no kept epoch, 35,056 have a known footprint and kept epochs, and 621 have kept epochs but an edge-clipped footprint. +There are **112,734 kept epochs** out of 120,421 considered: the masked-fraction cut rejects 5,122 and the central veto rejects another 2,565. + +| Diagnostic | Count | Denominator | Fraction | +| --- | ---: | ---: | ---: | +| Interpolated defect on the footprint, `NGMIX_N_EPOCH_INTERP_ISO > 0` | 173 | 35,056 known footprints | **0.4935%** | +| Instrument bits on the footprint, `NGMIX_IMAFLAGS_ISO > 0` | 174 | 35,056 known footprints | **0.4963%** | +| Interpolation anywhere in a kept stamp | 12,943 | 35,677 objects with kept epochs | **36.2783%** | +| Unknown footprint, both new columns `-1` | 621 | 35,677 objects with kept epochs | **1.7406%** | +| Footprint touches a stamp edge, before epoch cuts | 633 | 36,051 detections | **1.7558%** | +| Empty footprint | 0 | 36,051 detections | **0%** | + +Relative to all detections, the first two rates are 0.4799% and 0.4826%, and the stamp-wide rate is 35.9019%. +All 173 interpolation-positive objects also have instrument bits; one additional object has instrument bits without footprint interpolation. +The diagnostics distinguish flagged pixels from the subset actually interpolated, so those counts need not be identical. +The 1.74% edge-sentinel share is below the proposed few-percent threshold for reconsidering the truncation policy; this tile does not justify changing it. +The script leaves no-kept-epoch array entries at `-1` to distinguish them from measured rows, whereas final-catalogue never-fit rows use `0`. +After that final-catalogue convention, a literal `NGMIX_IMAFLAGS_ISO != 0` would select **795/36,051 (2.2052%)**, including the 621 unknowns; **174** is the count with positive instrument bits. +Consumers must handle the unknown sentinel and fit status explicitly. + +| Instrument bit | Objects carrying the bit | +| --- | ---: | +| `1` | 82 | +| `2` | 24 | +| `8` | 96 | + +These counts overlap. +The exact positive OR values are `1: 55`, `3: 23`, `8: 92`, `9: 3`, and `11: 1` objects. +No other instrument bits occur on known kept footprints. + +### Poisson comparison + +The requested independent-pixel reference is the object-averaged `1 - exp(-0.0028 * A * N)`, using footprint area `A` and epoch count `N` for each known footprint. +Here the mean area is **66.77 pixels**, the median is **31 pixels**, and the mean number of kept epochs is **3.160**. +Using each object's values rather than plugging their means into the exponential gives **31.75%**, compared with the measured **0.496%** bit-bearing rate. + +The central veto protects compact footprints: **78.22%** of known footprints lie entirely inside its radius, and only **19.12 pixels** lie outside it on average. +Replacing `A` by that outside-veto area lowers the reference to **6.81%**, still well above the measured rate. +The veto is defined on interpolable defects, so this remains an explanatory reference rather than an exact conditional model for instrument bits. +For the same known objects with kept epochs, **6.05%** touch an instrument flag in any candidate epoch before either cut; the corresponding unconditioned reference is **33.00%**. +Thus the disagreement already exists before the cuts, and epoch rejection reduces the observed fraction further. + +The 0.28% rate from #887 is an exposure-level DATASEC average, not a local independent probability calibrated for this tile's detected objects. +Bad columns and cosmic rays are spatially correlated, detection is not random with respect to defects, and conditioning on the epoch cuts changes the distribution again. +In the kept stamps here, **0.1448%** of all pixels carry instrument flags, while only **0.01491%** of known-footprint pixels do. +Those measurements and the before/after-cut comparison explain why treating every footprint pixel as an independent 0.28% trial overpredicts the object rate. +We have not measured a separate full-DATASEC rate for these seven exposures, so the relative contributions of local defect density, correlations and detection selection are not isolated. + +### Spatial check + +![Kept-epoch footprint defects on tile 188.307 over the exposure defect map](tile_defect_flags_188307.png) + +Orange points are the 173 interpolation-positive objects; grey points are other detections. +The blue map is `/project/6001537/cdaley/sp-products/smk-g17-int/nflagged_smk-g17-int.hsp`, sampled on the tile at NSIDE 131072. +Its value is **flagged CCD pixels per sky pixel**, summed over exposures, not a count of flagged exposures. +The map shows the expected thin bad-column and other structured defects. +Only **23/173** affected-object centroids land directly in a positive map cell, because the flagged pixel can be displaced from the centroid within the footprint. +Every affected centroid (**173/173**) is within 8 arcsec of a positive map cell, compared with **1,644/2,000 (82.2%)** deterministic sampled controls. +This is a coarse spatial consistency check, not an independent false-positive estimate: the map combines campaign exposures and its any-touch rasterization broadens thin defects. +The grid fraction with any positive map value is 6.97%; the corresponding fraction at other object centroids is 7.36%. + +### Reproduction and resource scope + +The private run root is `/scratch/cdaley/defectflag-check/run/188307/`. +`config.yaml`, `dryrun.log` and `path-audit.txt` record a **29-job** local Snakemake plan: seven retrievals, seven splits, seven exposure PSF chains and eight tile jobs. +All outputs were audited below `/scratch/cdaley/defectflag-check/` before execution. +We mounted `/project`, `/scratch` and `/home` read-only, then overlaid only that private prefix writable. +Neither production products/run directories nor `feeder/` or `latest/` were modified. + +The coadd-weight CLASSIC vignet FITS copy was too slow on shared scratch and was stopped; it is not an input to the diagnostic path. +`EpochSnakefile`, `epoch_vignets.ini`, `epoch-dryrun.log` and `epoch-run.log` record the private one-rule restart of **only the committed MULTI-EPOCH vignetmaker section**, using the completed detection, merged headers, exposure forest and PSF store. +It completed with zero runner errors and **36,051 rows in each of the five epoch stores**. +We did not mark the original full tile-vignet stage complete. +The private completion record verifies the actual five store row counts after the writer exits successfully. + +`measure.sh` and `measure-cached.log` record the no-fitting diagnostic command. +`--cache-stores` reads those SQLite stores into RAM without modifying them; exposure WCS objects are decoded once rather than repeatedly per object. +A synthetic end-to-end check verifies identical per-object arrays with and without this cache. +The final targeted pytest rerun was interrupted by the allocation time limit, so it adds no completed-suite result beyond the earlier 830-test and 31-test runs above. +The real measurement, spatial plot and cache-equivalence smoke had already completed. + +The single `def-mjhudson` allocation was **24921492**, with 16 CPUs, 40 GiB and a two-hour limit; it has ended, and we also issued `scancel`. +`diskusage_report` showed 398/1024 GiB scratch before and 369/1024 GiB afterward; these are shared, non-transactional quota snapshots rather than the private run's incremental size. +The requested runtime lacks pytest, so we used `/project/6001537/cdaley/containers/shapepipe-develop-240b37e4.sif`. + +The measured per-object arrays, original JSON and PNG are `/scratch/cdaley/defectflag-check/188.307.{npz,json,png}`. +The committed [JSON summary](tile_defect_flags_188307.json) preserves the full input paths, counts, denominators and spatial statistics. +Run `bash /scratch/cdaley/defectflag-check/run/188307/measure-exec.sh` in an equivalent allocation to repeat the measurement after changing its job ID, or use the command in `tile_defect_flags.py` with those paths. +`--plot-only --out /scratch/cdaley/defectflag-check/188.307 --nflagged --plot ` redraws the figure from existing arrays and refreshes the spatial statistics without rereading the stores. Back-filling would require retained or regenerated stores, but would not require shear refits. No back-fill has been attempted. diff --git a/scripts/validation/masked_pixels/tile_defect_flags.py b/scripts/validation/masked_pixels/tile_defect_flags.py index 375f130ee..bc76583dc 100644 --- a/scripts/validation/masked_pixels/tile_defect_flags.py +++ b/scripts/validation/masked_pixels/tile_defect_flags.py @@ -19,6 +19,8 @@ The Poisson expectation uses 0.0028 as the exposure DATASEC pixel flag rate from #887. It is an unconditioned reference, not a calibration: spatially correlated defects and the central veto violate its independence assumption. +After measuring, --plot-only --out PREFIX --nflagged MAP --plot PNG regenerates +just the figure and spatial comparison from PREFIX.npz and PREFIX.json. """ import argparse @@ -35,7 +37,7 @@ from sqlitedict import SqliteDict from shapepipe.modules.ngmix_package.ngmix import ( - EPOCH_INTERPOLATED_DEFECT_RADIUS, Tile_cat, centre_distance, + EPOCH_INTERPOLATED_DEFECT_RADIUS, Ngmix, Tile_cat, centre_distance, prepare_postage_stamps, ) @@ -48,10 +50,12 @@ def measure(args): n = len(tile.obj_id) data = dict(number=tile.obj_id, ra=tile.ra, dec=tile.dec) for key in ("imaflags_iso", "n_epoch_interp_iso", "n_epoch_interp", - "n_kept", "area", "area_outside_veto"): + "n_kept", "area", "area_outside_veto", "n_candidate", + "precut_imaflags_iso"): data[key] = np.full(n, -1, dtype=np.int64) data["footprint_clipped"] = np.zeros(n, dtype=bool) - status, dtypes = Counter(), Counter() + status, dtypes, cuts = Counter(), Counter(), Counter() + stamp_pixels = stamp_flagged = footprint_pixels = footprint_flagged = 0 with ExitStack() as stack: def store(path): if path is None: @@ -59,6 +63,12 @@ def store(path): if not Path(path).is_file(): raise FileNotFoundError(path) db = SqliteDict(path, flag="r") + if getattr(args, "cache_stores", False): + print(f"Caching {Path(path).name}", flush=True) + cached = dict(db.items()) + db.close() + print(f"Cached {len(cached)} rows", flush=True) + return cached stack.callback(db.close) return db @@ -68,8 +78,12 @@ def store(path): f_wcs_file=store(args.headers), bkg_vign_cat=store(args.background), bkg_rms_vign_cat=store(args.rms), ) + # Decode exposure-level WCS objects once rather than for every object. + vignet.f_wcs_file = dict(vignet.f_wcs_file.items()) radius = centre_distance(tile.seg[0].shape) for i, obj_id in enumerate(tile.obj_id): + if i % 1000 == 0: + print(f"{i}/{n}", flush=True) footprint = tile.seg[i] == obj_id data["area"][i] = footprint.sum() data["footprint_clipped"][i] = ( @@ -86,14 +100,36 @@ def store(path): if isinstance(psf, str) or isinstance(gal, str): status["empty_store_row"] += 1 continue + flags = vignet.flag_vign_cat[key] + data["n_candidate"][i] = len(psf) + if footprint.any() and not data["footprint_clipped"][i]: + bits = 0 + for name in psf: + own = Ngmix.MegaCamFlip(footprint, int(name.split("-")[1])) + values = np.asarray(flags[name]["VIGNET"])[own] + bits |= int(np.bitwise_or.reduce(values.astype(np.int64))) + data["precut_imaflags_iso"][i] = bits + # Reuse the fetched flag row when the production preparation reads + # it; no second SQLite fetch is needed for the pre-cut comparison. + per_object = SimpleNamespace(**(vars(vignet) | { + "flag_vign_cat": {key: flags}, + })) stamp = prepare_postage_stamps( - vignet, obj_id, i, tile, bkg_sub=args.background is not None, + per_object, obj_id, i, tile, bkg_sub=args.background is not None, psf_obj=psf, gal_obj=gal, blend_handling=args.blend_handling, defect_fill=args.defect_fill, ) data["n_kept"][i] = len(stamp.flags) - for flag in stamp.flags: + cuts.update(stamp.epoch_cuts) + for flag, epoch_seg in zip(stamp.flags, stamp.segs): dtypes[str(flag.dtype)] += 1 + flagged_pixels = flag != 0 + stamp_pixels += flag.size + stamp_flagged += int(flagged_pixels.sum()) + if stamp.imaflags_iso >= 0: + own = epoch_seg == obj_id + footprint_pixels += int(own.sum()) + footprint_flagged += int((flagged_pixels & own).sum()) if not stamp.flags: status["no_kept_epoch"] += 1 continue @@ -101,8 +137,6 @@ def store(path): if name in data: data[name][i] = value status["unknown_footprint" if stamp.imaflags_iso < 0 else "known"] += 1 - if i % 1000 == 0: - print(f"{i}/{n}", flush=True) valid = data["imaflags_iso"] >= 0 kept = data["n_kept"] > 0 interp = data["n_epoch_interp_iso"] > 0 @@ -110,7 +144,13 @@ def store(path): wide = data["n_epoch_interp"] > 0 summary = dict( inputs=vars(args), n_objects=n, status=dict(status), - kept_flag_vignet_dtypes=dict(dtypes), + kept_flag_vignet_dtypes=dict(dtypes), epoch_cuts=dict(cuts), + kept_stamp_flagged_pixel_fraction=( + stamp_flagged / stamp_pixels if stamp_pixels else None + ), + kept_known_footprint_flagged_pixel_fraction=( + footprint_flagged / footprint_pixels if footprint_pixels else None + ), footprints=dict( n_edge_touching=int(data["footprint_clipped"].sum()), fraction_edge_touching=float(data["footprint_clipped"].mean()), @@ -134,11 +174,24 @@ def store(path): interp_iso=float(interp[valid].mean()), imaflags_iso=float(flagged[valid].mean()), ) + summary["fraction_precut_instrument_flags_known_kept"] = float( + (data["precut_imaflags_iso"][valid] > 0).mean() + ) + summary["poisson_before_cuts_known_kept"] = float(np.mean( + -np.expm1(-0.0028 * data["area"][valid] * data["n_candidate"][valid]) + )) for area in ("area", "area_outside_veto"): summary[f"poisson_{area}"] = float(np.mean( -np.expm1(-0.0028 * data[area][valid] * data["n_kept"][valid]) )) summary["mean_area_known"] = float(data["area"][valid].mean()) + summary["median_area_known"] = float(np.median(data["area"][valid])) + summary["mean_area_outside_veto_known"] = float( + data["area_outside_veto"][valid].mean() + ) + summary["fraction_known_footprints_inside_veto"] = float( + (data["area_outside_veto"][valid] == 0).mean() + ) summary["mean_kept_epochs_known"] = float(data["n_kept"][valid].mean()) if kept.any(): summary["fraction_stamp_interp_kept"] = float(wide[kept].mean()) @@ -155,6 +208,8 @@ def store(path): def plot(args, data, summary): """Overlay interpolated-footprint objects on the retained nflagged map.""" import healsparse + import hpgeom + from matplotlib.colors import LogNorm from figstyle import INK2, ORANGE, PAPER if not args.nflagged: @@ -174,19 +229,41 @@ def plot(args, data, summary): summary["map_nflagged_positive_at_interp_iso"] = int( ((at_objects > 0) & selected).sum() ) + summary["map_grid_fraction_positive"] = float((image > 0).mean()) + summary["map_fraction_positive_at_other_centroids"] = float( + (at_objects[~selected] > 0).mean() + ) + # A flagged footprint pixel need not share the centroid's 1.6-arcsec + # HealSparse cell. Eight arcsec covers a 51-px stamp's half-diagonal + # (6.6 arcsec) plus one map-cell half-diagonal. + def nearby(rows): + return sum(bool((sky.get_values_pix(hpgeom.query_circle( + sky.nside_sparse, ra[i], dec[i], 8 / 3600, inclusive=True, + )) > 0).any()) for i in rows) + + selected_rows = np.flatnonzero(selected) + controls = np.flatnonzero(~selected & (data["imaflags_iso"] >= 0)) + rng = np.random.default_rng(0) + controls = rng.choice(controls, size=min(len(controls), 2000), replace=False) + summary["map_positive_within_8_arcsec"] = dict( + flagged_n=len(selected_rows), flagged_matches=nearby(selected_rows), + control_n=len(controls), control_matches=nearby(controls), + ) fig, ax = plt.subplots(figsize=(9, 8), facecolor=PAPER) ax.set_facecolor(PAPER) im = ax.imshow(image, origin="lower", extent=extent, cmap="Blues", + norm=LogNorm(vmin=1, vmax=max(2, np.nanmax(image))), interpolation="nearest", aspect="auto") ax.scatter(x[~selected], dec[~selected], s=0.6, c=INK2, alpha=.12, linewidths=0, label="Other tile detections") ax.scatter(x[selected], dec[selected], s=13, c=ORANGE, linewidths=.3, edgecolors="white", label=f"Interpolated defect on footprint ({selected.sum():,})") ax.set(xlabel=f"RA offset from {centre % 360:.3f}° [deg]", - ylabel="Declination [deg]", title="Kept-epoch defects on tile-object footprints") + ylabel="Declination [deg]", + title=f"{Path(args.out).name}: kept-epoch defects on object footprints") ax.invert_xaxis() ax.legend(loc="upper right", markerscale=2) - fig.colorbar(im, ax=ax, label="Flagged exposures per sky pixel (nflagged)") + fig.colorbar(im, ax=ax, label="Flagged CCD pixels per sky pixel (nflagged)") fig.tight_layout() Path(args.plot).parent.mkdir(parents=True, exist_ok=True) fig.savefig(args.plot, dpi=180) @@ -197,12 +274,27 @@ def plot(args, data, summary): def main(): parser = argparse.ArgumentParser(description=__doc__) for name in ("sexcat", "image", "flag", "weight", "psf", "headers", "out"): - parser.add_argument("--" + name, required=True) + parser.add_argument("--" + name, required=(name == "out")) for name in ("background", "rms", "nflagged", "plot"): parser.add_argument("--" + name) - parser.add_argument("--blend-handling", choices=("noisefill", "uberseg"), required=True) + parser.add_argument("--blend-handling", choices=("noisefill", "uberseg")) parser.add_argument("--defect-fill", choices=("interpolate", "noise"), default="interpolate") - measure(parser.parse_args()) + parser.add_argument("--plot-only", action="store_true") + parser.add_argument("--cache-stores", action="store_true", + help="Cache read-only stores in RAM before measuring") + args = parser.parse_args() + if args.plot_only: + if not args.plot or not args.nflagged: + parser.error("--plot-only requires --plot and --nflagged") + summary = json.loads(Path(args.out + ".json").read_text()) + with np.load(args.out + ".npz") as data: + plot(args, data, summary) + else: + for name in ("sexcat", "image", "flag", "weight", "psf", "headers", + "blend_handling"): + if getattr(args, name) is None: + parser.error(f"--{name.replace('_', '-')} is required for measurement") + measure(args) if __name__ == "__main__": diff --git a/scripts/validation/masked_pixels/tile_defect_flags_188307.json b/scripts/validation/masked_pixels/tile_defect_flags_188307.json new file mode 100644 index 000000000..3d7288aef --- /dev/null +++ b/scripts/validation/masked_pixels/tile_defect_flags_188307.json @@ -0,0 +1,79 @@ +{ + "inputs": { + "sexcat": "/scratch/cdaley/defectflag-check/run/188307/work/tiles/18/188.307/output/run_sp_tile_Sx/sextractor_runner/output/sexcat-188-307.fits", + "image": "/scratch/cdaley/defectflag-check/run/188307/epoch-output/run_sp_epoch_vignets/vignetmaker_runner/output/image_vignet-188-307.sqlite", + "flag": "/scratch/cdaley/defectflag-check/run/188307/epoch-output/run_sp_epoch_vignets/vignetmaker_runner/output/flag_vignet-188-307.sqlite", + "weight": "/scratch/cdaley/defectflag-check/run/188307/epoch-output/run_sp_epoch_vignets/vignetmaker_runner/output/weight_vignet-188-307.sqlite", + "psf": "/scratch/cdaley/defectflag-check/run/188307/tile-store/sp-c9b3bf5e-188.307/output/run_sp_tile_PiViVi/psfex_interp_runner/output/galaxy_psf-188-307.sqlite", + "headers": "/scratch/cdaley/defectflag-check/run/188307/tile-store/sp-c9b3bf5e-188.307/wcs/log_exp_headers-188-307.sqlite", + "out": "/scratch/cdaley/defectflag-check/188.307", + "background": "/scratch/cdaley/defectflag-check/run/188307/epoch-output/run_sp_epoch_vignets/vignetmaker_runner/output/background_vignet-188-307.sqlite", + "rms": "/scratch/cdaley/defectflag-check/run/188307/epoch-output/run_sp_epoch_vignets/vignetmaker_runner/output/background_rms_vignet-188-307.sqlite", + "nflagged": "/project/6001537/cdaley/sp-products/smk-g17-int/nflagged_smk-g17-int.hsp", + "plot": "/scratch/cdaley/defectflag-check/188.307.png", + "blend_handling": "uberseg", + "defect_fill": "interpolate", + "plot_only": false, + "cache_stores": true + }, + "n_objects": 36051, + "status": { + "no_kept_epoch": 374, + "known": 35056, + "unknown_footprint": 621 + }, + "kept_flag_vignet_dtypes": { + "int16": 112734 + }, + "epoch_cuts": { + "considered": 120421, + "masked_fraction": 5122, + "central_veto": 2565 + }, + "kept_stamp_flagged_pixel_fraction": 0.0014477094273839075, + "kept_known_footprint_flagged_pixel_fraction": 0.00014910042069069652, + "footprints": { + "n_edge_touching": 633, + "fraction_edge_touching": 0.017558458849962554, + "n_empty": 0, + "n_edge_touching_with_kept_epoch": 621, + "n_unknown_with_kept_epoch": 621 + }, + "n_interp_iso": 173, + "n_imaflags_iso": 174, + "n_stamp_interp": 12943, + "fractions_all_rows": { + "interp_iso": 0.004798757316024521, + "imaflags_iso": 0.004826495797620039, + "stamp_interp": 0.3590191672907825 + }, + "bits": { + "1": 82, + "2": 24, + "8": 96 + }, + "fractions_known_footprint": { + "interp_iso": 0.004934961204929256, + "imaflags_iso": 0.004963486992240986 + }, + "fraction_precut_instrument_flags_known_kept": 0.06053172067549065, + "poisson_before_cuts_known_kept": 0.33003992762146317, + "poisson_area": 0.3174972118531098, + "poisson_area_outside_veto": 0.06808699656730423, + "mean_area_known": 66.76725810132359, + "median_area_known": 31.0, + "mean_area_outside_veto_known": 19.11758329529895, + "fraction_known_footprints_inside_veto": 0.7822056138749429, + "mean_kept_epochs_known": 3.1601152441807394, + "fraction_stamp_interp_kept": 0.36278274518597414, + "fraction_unknown_footprint_kept": 0.017406172043613534, + "map_nflagged_positive_at_interp_iso": 23, + "map_grid_fraction_positive": 0.06971020408163266, + "map_fraction_positive_at_other_centroids": 0.07358269691733095, + "map_positive_within_8_arcsec": { + "flagged_n": 173, + "flagged_matches": 173, + "control_n": 2000, + "control_matches": 1644 + } +} diff --git a/scripts/validation/masked_pixels/tile_defect_flags_188307.png b/scripts/validation/masked_pixels/tile_defect_flags_188307.png new file mode 100644 index 0000000000000000000000000000000000000000..f70042ce1493a9c7b89d959f5cc7d38f85ce509d GIT binary patch literal 687384 zcmeFZcT`hvw>26Nq!*=&QUs+-lioo@6zNhWfK=(d6ObYu6cI!^C?JI1YotpJNC~|P zp@$Yi3*3#r_dVyH_k8!?@2@)s$RK00^W=HfnrqIv);3C0U4@MJ4lxJ>B2#_#L>mOU zsS5((XA|B8Ug2CVRffnJ(>xHvd_IJ|ws;$`LL{?^$^_>rLC zBOzWETMrKxcS(MJ$N#$Gk+Yi(e{kF`a1#)b%QHiF5Qwqt`Wpv)=i^-v4hW?BL_yal zeLKTDfZFP^M^1RAa)Xs9C-fPWGPVDW&mHgIINL~}qy6T6MICIK4x=+O{ghGCXd4N} zWSurJB@Hz`UYnxU#)E41wC!Do)KjdW@58Dp_*ljq9=6q#|3UANG}hPp*f{)n(hDC? zf$~4T!nodeiv_a~hWyX3+yW(Hn*VqW_?7Wlj^cmb0nhvfm+pW5eF%#KOE~fWemAmj znymlp+VXeQ-2bsG;8)QPisAorQNRM)Z)z9*uLUSf%K!h1{lB`iA#-z}ztPdWKcd3C z@}nRN2lU#f?UY{13*+Q;tAb6~ea-)5A%-6Z1o7?dqL*-cB38$lLTBp7fNGVZ0=cX> zJIh;GD2;46*xPs5m$k(6!!C$R934McG300xde8Z)jBGt@l{;yuqvbszO^`mL2Q6|S zGh$p z<|DirvOd>z%K}3m6G|aTnWOpZRqDRa-7?tz!`4ebE+ZqO;bF~$iY8Qi2Zh`~EwqO@ z2sdje4cW+CknSeb&Fogbbv}E2Kk17-Hc)T}xl|H0_Y;MH%hNZGG8sPY`!eb|2fecz zHOD5$Tc$K7IR~4?WlKm{Jv|>^i9}Gsll$z@2d|fq^Ji-_&NY?f5`Pw%t3RY2+q&k| zJpG#mZ4rMnXKsK@Y&H@B4x-}Gg9h{V`(~k<;IgIE+*LO4wpaS5p9aSP53j zg0LiB3?fCAR#O$`9XHd1+Ta4t^8xIiR7bb&+_~d2{yAl3IVR|xJd2rOelF(;Jd|Q7 zKg~(m-`{^8+2W7smhoN)=hBH6;ah1jtlAq8hD)NJarnpHGxiAY!`wP(Wc9=4(JyV+ z7`ozG-uLkq>k_Q?j?8>2i1-s9jVNiZY+7v}#Qf2|VKw zN6u!keCjm-W7-Iq-z;~s31C<0DiGk6XmB!Bt>Z#Gg5&lvYHK^LP(hihxr7bj467gS47S2$fd>c^b)6Rpl$Di^=%4j|7yAjn3j zeExh?5mDI5R=DBX%}|i1k@YY#H_|NIrd4)y{?xo2QSVxh%xtR8@)zoQ{Yh$>3EvAE zXw+9z%N9rS7uKG@183L?Lr&glfx+=}UUM`SUJGO|QXaW~dp$h*`**bI<$=Bg@dqF> zNd0$2$jQ*?$$D2d*$bC=+SZ~HzV}&I%eM#4kJt!9G^;CTwM|UIWiYc;N)fwZ`S}mW z5D2zg2`gJ=Da$no)1KoUK2R{7X#f$$1ybSlQB0P5zwgHP4T-*_wornlKcCc}=w!(n z)Y#Ev1|Iu5%r`TDT)JrVEy`S1%6phB4*MjK=~b>jan5II{Z9{QPLOgag6#P5V4NGX zUKuQrOzr6jvNu6^p$g6WXFT-cE@onNku19J$>)fG<&4b@u6x-g&@vm@|AWy1e`~|CX!@eZ5#^5OULbiZs|S~u_Nv#)KS{3qZ^@nZh6A|` z>`xa*R=`&}?E!dY`fr$NQ1kgv+N6D%Ykp4i^z^jL<*8g3qjM`h$Y-zZYK2?3a*Eff z^;B;M)7x`5PKC*raKu*kd`3j4tXq{MEBf@1jwNkyIX}(0tDnE+gZ@MgRrks;(Td|t z(D*pQRxZnZ$gmbcg*W0pRc?L*WcD-olD+UjBopC|)<4jZ+a!?i>AIpNg#p<~&88Ax zPHM6mX4Ri!FK6m4lOR|>KuX> zS|szc(yG6nwRiq%Q3-mo#s~t?jqUOEET_cmK zww!8{j|1`>A77g6?eXD`V2R9Bkj?~wbYH#lbbmn@vZ&yn^iKg?fAmWUsD;b^A+lo= zA2MzlgQYwQSRe58JKwppsby-K@G?zQA=&CzWxDrZEJK;C+|?frWNu^QSL2Bftf8(K z`#nnyNL_f+hA_D60o_%*XBs$Q{{bKM)P_y^`)BNL{C29}!k0Uy0CJ+j?8zQ4d&Y@C zH+;{_yHU4rNmpH6iwla+a=ZIbU0NE|C}KZhHUXvBo9p{*>*%Qb>J`mXQs(#9+ZT9t zB*+)EHh|fK(V60cE)Ei71l|4+az=+>=F`G`57(|iwXnFj19H9_y1zm(LjA^XZ7?mN zhv|wxnQu#b4ZCn96P+r&(BDw1Uk0vUU1OZ1z>V15=1*`CjnR-jrDhJ2rW+iz6K|Yv zVXRnyoWB7M!u2(IwvlttUR#90{$Sl?YoJt7QPFUnGxK)sR_UZOEYO6@1v!MpV}^(` zGKx7W0t-H-m!egbz z6x~Yq&RauxeEq{gPLG5$kqB_#zT`uQtn-zhB9OK7l?Ek8tp*u#078<+H zkaYhrtod6c5fzauaBopQ7frgNHK?UNCL@wv5`=hlbd@VbqJR; zg`zh)sp3^0rnLcBlwj4LylNs9aN;r>w6IWdb$Pye5!e+?W5f{ay^xXWeFcCriC*DR z*ySDH)zA8_STeC?;1E_JvL5T|CmoKk9@n0iYmidZOq|$U$i}FR%lSWbD0vm12D#YcUfK{LC2?BfN z@;VuDigj4ZZ5F9_oT+VwWw`bA`nFFiB9dC~)*T`}QPs7z%YDW9Z>Fn}_0ezTUzMV& zF!Mf_r*DCP7ICbjA~sp+4%7!T@lBl_*#`-=o*cY@+jNU9K#Jo!SR?K_6lqS~LVUXy zzO4mY8@87K_jJOMaTUBE@vD%m$= zz4`quufClLWY8KQ5ufIE-`cY@O?#B;{2U%MpCCLWJ)VmM;j1)9cY7#-EAvs<-sNR6~u3PHkDfk+X z2CnXj;>Zk96P0Vux&=cbh#=z|$crh~ej#-UG^sT9=&_DDd$ylez`;gAWRDba5xh)f zh;#MqsA=cBIGcB;h@slWS0xxC#Q5U1)HPG-WWM*YKn=%C0%4uxrjdMhasRJC5m02h z79|T*KsaXxLGGx_?5OdYG*ZJ_x7}tTJ*MU7Uj*$Up=&R{$Xy9L%+!*BT!s}0mamu~ z7^&)N{`Cym9Fk$*h+yJ!bC-h5-(uCE9IK~=KvKX zI(M;7Y(BOE2(z&Lh(AeHZ@R?9x1jFB6B7^U=Yn8?u6u0+Io>){lG8I4PFTNyia#U+ z*grZm6<_@ZYi!@XRnC$P1k$u((sy066S&4@>Hgrd3{%2Nt=vULc8lC$i(&oY5X=sS zKp_BD z46+V33!9eILP6l;p=}A~rMv{$35eW!_Kj=Cl;@3_*R*6{6m#l41K_=lpsE6b_3L-7 z>nIZ`EzsSnqIcv!8Ys6=G5MsvSE+*i+M(FLZ=nW+7OZ9eA%Gb<6&1I1^6$m~yPd;Z z^!|^XgTvz%IBA)GWkHtdKZHv9FQM9(d9yy8A+Xg_w?8%tU!({+5e7~(!E?(ErnXx> zN_lboe#_;zDu*4>RE$Zo@t0NJ1Pc!4&>uQwLS*j0=jW}W3=FOKp(8x6B_*d#>JFng zvTlNa0=HDHqfyg>JQ`K++HEQZ^R66M#|NT56EOSH>w_g|rW|N*5Rm_Lhs%VTs3Z(E z?RvwmZU8)+T%K&-wSx@drUVPX%so$sbSLburkI)o@aM84q9drrQT}^$dis476-CcV z0@f6~`xSus*iY)RWY?Q@T?-dODc**Zn|FA1mzztm(G%J!eg@dE8p(@g0+}DMSmPk? zG;4Dyr`gOe#{-{yeaQf#S0|3H7QQ~=UzM@OPM7zrhXQzzU(=jgUr_~+(+@KzTr4Q% zo0|fK=nEVB>Iguy-G^jJfH2hut7jv?{`U$py&i(h&W`jys29kfuH47MvG?$quj1U27tbQC<=94C=$*$c@X@8z_-Y~IeG>DLE% z)NvhRW!k=U4P3xUd`(pv8zo08xA*Mhwz$U6c%ROY$OvywY9YQ+HDBX;`LNpsJ(V0n z-QJ!tUTW7s`XnpJ@O^-LG(u)-6#g$fXZjbOcYo55i3_@1fjCT6hM%8$FO+n)1f54D zPuA#Kopx!*8$ zx0yYL-<@oP(=8Uzfz8s?CiyCPR~okGAV;So>=F?6^^;3TXduEdr$P{dc zKJSscoYyU^9936xxeB}*-n=G3FK33Eyq!zS2|!*8fuu|)n80b2h9|Mca|MlU)aH1x zpv{m9NVkhg1HkjZ^W&~+o6NE4^NnU-1E8Gj91f<{0lUut(nVF$fqf=~Y+hVv@0(7S zQ8>;kms5DqC7*`CO@L(vQpd|~J_E@MrFLp3R~$ zU|`?I2Pw@^2)((y6v4EooQj)dts(T>6LRZsx_wdj7p_|xcEcf%TbaIVhM=O{#b0a> z%>nev1M-zo)Zym^e%4Ybg={A>i@#d}d4Z=>01i)cYR)ON^38!lNw1L+G~O;SLECY# z%Cpy<=A6^ez;NvhSi?^VHDvec$)&c3)BW>CO8JLWGu*mj9Yl*)Be?hz8tn>XXAup^ z^qoyGOtq8hWZw|^raM+21nOauVxW*a3jV=je!M;P%eo2_Tx&q8|N4iK&!$B|RTho1 zg6IrQ!9{UdtEgiS;B18QY@6SOoS4%KzQF?o@%z>O6pQoY<@hRc8c}6q?<(BTb9_{IF#R7Y}eqIeX6Zl28#HxDx@~B|O6=~uPu)HM>5PEOX_u@OuDYCi#OlR^#PWgGc1`K{yiU}!^&OEZhYUpJgAn=T&9U5#khSGcWm z?;sOo&%`;H8t(C?lP1XA^O@VOL1Y@#IxtWxw*j17nn#xj@Jb9S|B#Hr%ToPl&M+v7 zY@1{6jx(}G`+*V%hV$*Jx`*FJd%{hdtV@m$-w`mG8OKw-->v?%iVQq^{8o;Gib%r~ ziU_ClcpSgdoyx5nEdtX!d1WKG29K|-6jBR1RepU|+%Et)YxNsz#;&7-b@!=>%}@5s zC%qQ}T|2)qFa2WUaOE4ySVr6t(~d?PLTl*`V`o1iECu)M$zy70XP1%F3{^8utu-W^ zgPE2|oByUd|JDtchw}P~aDvM~W_=vo!uF3d&R1-Y>KAQ5!34|!;u=KDHinSfQBPd@ z`NYgsZ4~E9jO&(HQ*0KYRYPN?L=hdf4i23FO5Td;D)-zH<%{A=cvakpr+I^b)MYkV zh55XURs##If{-duN+11_yZ=!LRHSGH5dRPCy}!?@YX-qnmhWfZ;tY1ECYrb4LX}bu zyZM}3_i4d(govY268>RoZkY<$y?q~q$oQ*rrJ$E8D*6%?$ECA!dK3LlX6Ae8*>w=M}Os^o?uh$4`$H&h$-JjiCv5L~f^vIo|~F z8kGOk?iF^v4A?<_d>GBf24FYrDiAmLPTqNXRnT_*PS-L#sF94n_3^#Y0{Bm3=ax{e zo~$rrmVdV{y$yF*(ki$amL#JYNEtsNtXE51oi?c*QKSx{Ps*j zCm$qJ1sj)aOS!xweAEzovj`$s}i3_78{fs%8AV4H`7HYB8 zEEpRZ?H@~_7I&LaUAf9xD+O{*ggGws99bT@Ta{R&l{tmK$$fBM`4HP)Cb7-6Ri6Y$ z+zWq``JcO^XJh!a*OH~b(yG2dw}|Id^01QgsbSas35IFnMPudNx`oygfHvS1My-@k za|>HHPkS!3o{JYZ6`FqsGJK>oh z^#V*0-}C_yB8O09wa^}@8J#=oJOnfzp!eEI-H8|gF<;T4sipO9dioAX9*}r|LGI*> zb%(cPV`s~pVY#1Qk+Q3{>u6|iqj$d+g-o&||^|Rn}y6&}*jz)mpr!K!k@Y zzV;>WUd_@+6=noP7O(;D&gl^2z5`vN%Yv2xpBD#YQ15apUi^UGi4=2f$-Ie?J6O*< z*=_n9zxf4yN$|{|Q!Yz%M*<3NnLawacDT1B-EFojdM2BjMXan$_J^qU+n1qR8&TM~ zpeyWSt5v{&GK&h2KrDB~2nL>!f{GS5SWLj7Hz=e*R#sM^;JSG~0_PTVZD^*a+}Y(! z+=`BdOz3#=={Y^;{NK<6b~3F+7~mOzy;Rej-`IEu(DwKsKzt48qycI9H{8}4-gd<| zKa6>2>b{ON;wu3Nbq_~>^iVlEW4xske z(qErT^U3f2*?oO#{A;n5Dkn?ui9_QZ`**c1XJ_XgB96>x{?-d2RFzF;>igo-)~Dy# zS@)r#ru>Qut@Bi0hEssHIV&m2#BaTpOI5qck>%ECZ349!N>6b4-CMx+sTTx_mpJ1* z+89+iZ_drd1*B>G@wuT*R>YO`38nUr3$KFFx6Ow`O54L~VSJIz|L_kWz6}!KEu`g)ai29e66~^$FCFiD{uID9O3!*<4JEI!Nw{mj( ze{@d%jB2>c!p@3#1Y8z>ro&8JPpy1TS2qBId8XF6=bEko5Yf4A;-#RZbVJl(7C4ub z?1l!zrh=j(V%SL+{9yWzjiDWL=0xc@2Qk8Z0ejW1ZyFfa*l0w0Pm&t&>TPXpz<2)dAompn7D6k^~7H{1(Y5ScOrt zFJFR!*9Hf#9S9ZI6&w+*ZNPgr(m&OVR@V^R45RQ;x92ld_5Iq!vg z@8m&%>`S=95Cx#oh)Jt4&aie?kbD`q#if5x% zcsym}{21+|s3IMY-kMn2`u^ftwHZ4%5E0z-$ev+6%1Co=xaABxe8jH$yW2IZ1=GH_ z)ba6ofk4nDZouJS9cZiqT+=rz*+hS6A7ScwCK%*|83IU#>lVtX1Lb_q)*U-V*71HP z_P|!6)@d%<_mJ*SI;sb-wUB}b_Ho_>UAF!?Ht+9A9em}z2^q{WU8NR!*aeoFQ`M09 z7v1#vTfB4A3>!UCADSWefzH#K;THB!v-=@6_w`0Mp<1a@Od+51ISpF;;@*BqKR9r^ z_F{&|n!M+2AjzxAJ2S@D3+@McBtNf>;nJO>_Qh*wb7>bgA7yIlB|)vM!2P*>#u@g> z+l2jgTM%>5m=+VYvew?@%GEdowLzN=I~(WM0$K{D*fX)~UJ)RW44U`at0w(XEg+Dm zXN}KB7jkcg>DHDSk9e#L1uV6#|)YgmnyxEkrxRu;k2e!0v`)*->9x z1J9TDdzcgev-cdDlfwj95fN8dlrrdNv+toG(D>Uq1cZ_w?N4;)B{{P!KyA+ih9ocn zKra;p0$QfhDgr&6{QNy^R~R#OT=3H}u&2+`sg3LOTbnnh(KMj3mZP;a8L0~dkZf&O z4ECPf`Dqs*U4T?{mpi3)m+O6XKMQj9!*+8n42H!r=UdTOKz8t1T^-4pWWtHUfTbUv(x+q=;}z+EjqA^?7Z(encCOK zBac4@h~MqzDo77YSPKXSEDC;A!}w{)(kZ>%CFPU6wf4bh>>(@hQnb!3r*|$-dKP9L z3d&qv9PVuvk5-)TS{%qy;o+?GFo#gceu=v7?P@VD0F9kNBgnv)xpHKtzj$Mrlmx2T zeeV|=1uH4UdAt_8D3Oe?6Om%v!6w9xmw5zL)Il-r;|En=udBt*+@Hn{cfW@&^R&O z6+xrX^~(qSd%IH3x%y?2{Rh{Uhc^MFZZ1JM)LaxUwfvR$;ok%Cxp(b%XnZU9H`jX_ z1`L{fZ(;;U%)XNK#&4^WMI~?7>VzyI{J1!>WJGm!^_vl7mc?jyTx1VbVIlPAPn;X& zz@}@0!MjDRiv4~Ey#=>4I&nZ8nJyt#rHQ)S;Xv)kZSWQTnwIvl*HY`1NqAf6P4lC5 zlj^W0C&O9?iND==V`DAT1(GI|?9CfDNnKub7i0}}{DZEHQ>U1_{`T|)IU zB+R;Z8aTS+W=1_3NU_tYJ||0)ZwpxS5s-`h*Sx38Bp~@}-e}<=@i-nzqZU7r)|2HI zfWdn&ge4p35(BzVd~9s&OL8F4T<#%*wGe@T7WLH?gx;CMPe*lBNZpRriSKlZLlZ*rCojE zb#DtZ@1IvY(JUhkxD?A>FZOo47jjB!$UHqg*GCJ80A~I;jsxbi|7KA5&zHOeML^)I z^}sO@xhC7k102gYVYi#M&b^yd4@a{REXJd zkd0TFeKQAC12(g}I^eVUqhP-(%aj!4a}d=dXcZmERCdqtrDWR{6M8Eaop+bRVi|eZ zJ;;O`e2q=0p>!gU_(ge#0MN|I?02u?@IN?J{yP$Y^RJ1YNL)dxlK|7x9CQlM6YY?} zTmajwMP`1~$HQ;YexlUap!`HI2TAs?K76R+zepc&lOq>^829tThat?e4oD}5f3Cq( zxDa5KGx(>U0-yy`ANI-*#2j=XD2w%A0_^k~AYer0o8zyC|GwnnpTUZM3|g5av3VaV zVDaih3#IiJl}>j5b=z{*|DJN{B$t1Dr^^$e-eqW$F*Wz{UxQAO%vt{O<8|CO?!S=%;_m&c{Q{ zll$+Szv^(76?P844|crTdB%QweSLknJ9k*55XKCugZ6_>NI>>-Vo?u;TDOE{f&1Zphv^a&0WrCd;_#@nkv+C z33s=scRx_~&ElWJd)AA0Q@?Pq_o}I|4yxOE7XkeGu>RcZdMMVM`=`M*r*Lt+93gLk1N(Z^q_CQ@fZxAi05MC zcJEVVzg>T{I!>u*E*;P}WF&=l{IgU;rZRga0a_f4@kTDwl&lOTO zA8)LsGU9;hm?nDM=eiqMFLJ)1wB^{Fe?6b5DP9Qmg^LLPy3xU4<$qDxls3%APIp<= zlqOH949o(CSxnNQwOb#T+Ee7o8L7$H@Dr$;HbP4+Ow~awcj9O-4)kK_+8(5Vl#>Oa zwa6B;#*MP1k1xi!kwExlQgd}fUL=Oa$CHEhUX}Oaadw}9A&!Ee#gubzuDhGuZl+;_ zYt|Fb(eC_MSsOLJ;8OlPuup(%j2<~?v|@!TW;Y9l17kETBU;%WPT`B@n_fVh3J57O z+)gD0sTRxKj$%nwyPu)1l=v+q9Qm%&C@G@2LRggb_ARxUb?cpS!C{hzB!;19=Qu9? zUnL>7NDULy+Ihoya;}kgQEc{nTUclVKB(MSr(~;^hNR9^JkFxDz`!Q>Rg2zEMaM~3 zY60k$@`$Txerpxcvu;-m5gX`;*NXUqR8QBJV#Oq`VC-w|AEIISppZf$0h}7_HdZa! zs$VF*c4eIA9^v~@afsz(B&xY`XTB^pOG3hicqy@Es{$G43<9mqm75OPboatK1vu3ClTTqHQE3w~x$J<|z)_gVMcsHdaikjL`URfFc{zt-Nrb#qpFl>3N zn*v|T#@GLCV6poZB9CnL~kmHE2=3XDfy)E4L(b z=z}Dfp?kyG!Pm2uI?8S3Ms1hKZobjzb#eKDdko>7A%V3!K$L6sZq<-gLJc58XG6z5 zQBBF}f@?Jzd_NwTlnK5+ljmNYJl2?JUGyW!* zZe+;*=$+}U7onpmkBYZ|LU6oOBo@!Oq?P07_JZ5W)XBSf^6+eiYT{`{zu0jTWm|sA zS+pkwMdoiUElbYS$NC0;Xt9%1gK;e4GkH%HzKl0Nxo-u^HiX>RG}P45`DCG5(z7vH z=>;74OE}$x^KsVb{^Ml*^*ev&P~V`Moax? z@)!c+h&TBKuwf; z;&&`&lB;6AmbzRZQoRo;puG)bs1Kq9$k-q&;*Zbp0H%_>TJQRfubVD;?kmVJUOo}p z|Cybr!hp2bGLbR=nRDC*;)O&R%-X%b>gdLE={H<8@m01_vdD;;Tme zOB8q4QOM8eU2nTD4A{HGxsfj4cx_63{lCmvRF`6+rzDov91evYzJ1z9skT$D96~h{ z68+<6xNJ=D`|R$;Icw;|(vZlV8*U#nOF{1@vP)?18MPc55M?qzFe?BQ-Xe^Q$E)*# zE4Hk*R6w(28$h7FX}EV)vwrewv-@T^zu_0HamEJ^t=y;i3m&A)R??=PJ$=&O=idBy z^4yl)pz%-%wXx`RuXu&!ev?@aVUIs(oaUl#vj#- z`?CCzg{X)#?_0&#g1RHI38YoY0=iw_z~>N0!joc)50`LN`h_pkicC@kpYq8Fx~iGOw*WMs>Fku05`L$1$^^*J^|iyyFjozr?{@J4=;&w(fSJ#nczdwz zo2zQn$^iKM!OjgS#;FM6E9nIeCvdqpAc-2$K!`dlwB7>(X!g?Ik&a4qt7O+v>D}ka znT~{m7knQMwV9ornxPs%))Kx4C7oYBsQl0{f?x!~oEE=*DU#nD8>^TK)a~eK;-2>I zCz1>PO3bnk0i2^G1kg^xDW<|yoy+Ujz3g)Sr0|Khv!0>M%C&zgd3jyD`BsRI;zl3h&g0)bT#Pn2gUp+R!3g&`1ypImuBZrpYxs zqQDyR=){-I`?69)%PYTecaBdF&9LjQuW!F7`NY3vO?{9%xLHz)iU7ZCp>&?>3>@KN z-)Ky$1x|xc;yWZ`cK@;tbbzcJPBG^3LEP!*XXpBk@H%~^ch#p`u=xAZs0^6#daja8 zSSNiKuU%m@Wn12m4CwH1{t)d!mp2;BbRY4-&jZlO;)wpx@%~F1Wv%l!T&hvh6lwgm zh2|0uNq=}-Dm!x|;WVyoMx(3IrV!RBPve+0R*HVd@rO=_DTc8i zXYj~fk6Q_sUmQ*A;U?QA06JMPal5%J+$+E zFe)96uh7Yb35v!S+oeDgjd*rfdLwLq&H0HK-%vkpmZgYD!E)OON5*CZpO9bcLDfe^J!5Wtv;!My%li;Y2Q6>(23o+C`kUXBq!3Zz(VVnZ%G!3YN95sCM6*u zehXvz$AKGG=EMY3kWk?t}~Io1l{=AZ_N!XAlxKH(nAc zWSjb=@p^yBcq7ku&rF>%W{&gDmrWVqt%#cR++QC?xj3IwhD6oGaYjcb|G|2>rtW0|8cqY5Ra#>A)`433RDQobH z0(OsYvt?P;Nl$W1@S~novrw3eeJR1IV`>Lu>>pAqj6d*ZvsP1MCe7JQZ2Cz2{!n`Q z;-u9K0nSvQ=+B}h4%fAcR-vHz8IYCq{Nr|KQwzQS>cjFbT7aHzt{c~H(O&DX0!dW#(lq-IOs`%crf;`IQ=hbEVw ze(2>*Li1BDR4X^QzMg|HR8_WCR*H@q;{w_s*2x2Nzny1h*Q$Fk#k0=z4bI`0BPjU6 zcs^12AaC>F84~k?Z`jH;>-l$p=Y(eMtnz}ZT3_oT3$f zG5VL8(qB+M1xy7bDA-fSE1Z{LL~{x3!lgPP;qCLfR-a5#QgW|ps)SN;4Q=Qf`1*Kr zx;d~UVPt4p^lD~M+`{9fdKxODN{GKK`j%*Ghuzf4S_FC&VfE_f@b`dfsk+G%#Q3Wl zA-)5yfre*zrD~BDNW-O9oa+QwJ2$FfmHW4(TkG+Izp{1NX*{k|EmBl3T_)DRjm9>3W>RN+=82wdh94~*hMCi>uu+nQe?>9}!-;wk z$E9xiMZYPHvPgYMkmrfk6XS;<%M*+K4yf$Ba99m0ZH>B*r=5O!_*I#SF&d@gRUddI zOu{FH9z@ro#=U$*ZXX%%FY5btyEQXPo&gW^jEnK3o_KA2l-2|U{s3Tkn8(cfy}Syp z^|ReFWocBhHFLOOzu@#z{=Sxvk3g|dtWQOB%+oxhcCC3jTaIKgFWPzQ`amMj&A3kJ z?gQ2K3B2>%4Qc7S{IfTxat#bAD)^JtsV!cJ{s^iws)9-v;}uwg<-kccXJ~=K9O(iw zQ27-e`UM1Gvm3?(gf`Xrh*xAnNKf{Mv>HkT;L2YveB@&4_^o`L%fQuq{*JN}@qTWbs>(@-fm zN$*rg{mJ)Nh0r6W`NtXRY-xjviR`7pX@dK^7#xr@*rBXhu+l~6=C81t<1P5``I;j> z=MuJN>fbTK0Qbh)J9@-vs^r|;->q1|(TBoq&%^R8;&&*DV1DzJUmB1aC*Dvx&G zUxuRY#$mtn^Y%)G?vSda?=ER_nw>W2bdx*p&EEWJYH_>H9DanZ@RYIG&{Lxq3f|fW zs8_kWGs|tktNE^u-l=|<^Z2%*S`Wz#ivgWNA^g2el2PZb&P+Cq)4yvNGxbP9rpkJR7(Fs&akz;hiA_c8u(Ec{W}9Ah^LwcyGH2*^p4=y|`cu7G&wfi+Hz&#X zo$l_%Uf&ur%a!d{MZ3eJCt8E=a#11ds*&E`SYCzk^}u%|I_*Pegf4U~*SFN!Zv(Gq zZ<6Cm>{v~7LFb&HxY)4pcm0+^wK-CXIRf0=2kyd=?^v{BtoK*8Zpc530}`1MlLhsi zMTa3@R&Z4(SYxEfG=AlM5LvUzT1QCfJpO7peHJw(JTh>fVZBwrz-g_sEKtW*qFMBC ztubkZ@XwG{#gMiF<+~4_3));nJ3FWEZ%@}-{SFRK;WR4|r%}EI^y{K?>xGJ3qIvtT zpSqAN$=r~LzH+$vbRk$f1F0v$55l0@y%5GNbd@gpOD^=lLjoLlg>>m_y-PJElQt>5 z{hGhE!6%Q)4#As9*nSK2dNg3v&=q9A3kuEBOV2g||2?)){5G{agHW zPO;(dFb@h@sg-8*KGF9^hd8F`0mejL(Q~rAKYEiJU5bn#cJ;!ry9(hQdw2B?ZPQbL z9?n5@Ia`dJ)KYLO77}o4gao(D-=B?=bjkS>!Kw*T)PbnVjWnEkQ0gWB-VE?yNN}mI zWmdQ%tLgGH;Ypu(+gN$Z&UDD@a^nZUq9|j!J`7sFRPt`#vl0d^4(FWOG;k6EEW&B( z{dVJO(JSE-rhwqDq;vul62DETxuKOLF_xlCVAsv4A{7+c>iy`Us`K}LDwi9JM%Sl3 zWqToY?%NAA6gPJ(Zu5l;e5c-@!D?DNpq?$Q7@Gcoy%wMLY?^Fx21tx9In4dk2Q;s! zwId@SCVltt!h$n8HcM8Vjv3rU5%T;tmS=pQJZanx1ty&4z*J2;pV@?iKfo+?E~^kR z-o%G{MBxWC>k#w0hda$P!2uM)d5Ns4&kWWy1v3fG_8h%-ZYgUWPLkdXovqpqBa(9M zf0L;;K=R=BW3j?thfIWS^HIDL&n0(c8JXaOc2nMr*?9dR(q*s?>A)r@1sxOj!K4Ya zL;5T#t)fy$_h75x&r=+&Sl?Qn!{e#V^ZZo0##^wCx1YZ%u6e1fGb>O&(Evuu!uSRh zhdGDipXLFrDb$D_B^aD};n!=q`ZP~V#kD?;KAbd>o%N$u33Z3wl9jZ951nEK)H~fs zk477o;*Z`r@W-ZcErA#n3`IY4l<(a)IPCvA7_-`#)}Yy<*G}1$)YAT~lXk{WAv--< zY5Rza3$JO<&2P2V)s&v*LwGJ%>~AhS(&vW4ZRw1b$_i#bUHa{Du@Opovt4fdG~F#} zD3kR7ndIeuKKs0>*M>E9AF#`-0dLS7Ewr(#XzEz(ivr>G$R1WzF3w6^{4Ox8BKLFp zpykJ*Z?}lwQxcw0;XVr=wGQVRcw_t46H%~olZYPiHn_^%qN;#)<9wQSYZ1(L(+CIT zrsrmIxhP*$6_`Q1*j+s_@$nW>bb1&3xTNiP$58?u^G>g%RmtB5$Fk(yRx~~%hyo8c z*g7j4^5|hO8Ka4JZ}qMN?X%~`0h|KNO(MP6*J2NY&)m7k-(g6%VfwWPGv=I*X;uD6&zAyuR{n-0u|~L^YIDxpIQgW719|84dx@( zCERt{ulyK&0%aW4_&GxgKkGCLcjey8h%M3WSof7-17%C{SadukfZOgu+yR@Dd-1|! zJk#S2;9K1#H$E<=KXfq6$K7Ufj9|H04t(;eGrn%<`j#$EFlc&OQY}kNSqbO&1~s07 zyrUyoQPJ}PcHoi8$&Y^-#0WOa_eXiw3#+c5(mFnR?$+YBz%>22I;N$N*!0>tPN;VO zCraH!!YBTzMoNt}d&kdNqZAA)yQ;iOL}}JetqflkcGU_BvfI-kbe=?fwuxoXdoErD zxjLO30Veiq;)*bdiFQ@Z^Ut4cA-+h>=`?>klJ;?(nVNi_7zR8}VJQ^#_3ND;y(3Z{ z9&P^h@Yaa9%^0&UDAuJh9kO2;bl^xfE$Rj)FmG(+XU`gA1Bw2+9#%p?*o`%5iq z!nwP*LuztLgbP$WJ2|Cqf9&o1W%**|hXb^ym;T;kfoQ&#u^W1c14;*z7IXRvzr1RD zu?Az|WBvx`UZv#g5L%o_dfo9~KVZdjO2s!MtPI*1Rg^!nT2}h57(T+it)kM&c4oAr5^)089zMkC_apl@lzSIfhzlMqe8jIeV0G?H+id33GXc@yo22Pb3a$94_BM^ zV|K7fnt?6MV@Ra=bRQ!oUC@Uj$-NJ)W1cJ|FD8?0o=Z*2755`SbPi+Vw-l z?vyd%<5pJ2eoHzSJ8l;=m15fpp~TY?gAHwZ`^TTC{>Uu$qgLKX+9^cwYIDb(fomx7 z-bw2NRFxNSoy=pjKD58~laiEtzMHW~%h3*)IXO+sZO6JK&KAvfW<9YtjEg2aoF-aW zq9dP5efmA~dW-!iqQcn;W5H#S-Cr3i&+>ule2&}D{|BY=E!Ix~VooFQ7{*G_$4R>P3o_1XJ`*WxdiKsmGIY@1R4Q!z<7DY@fIqo35!yG1{>K@e zhCiK!_Poh@*Bt0b+ZyUSPY(~yJW+1qsCm0pOIzI|$wGJ^PvKSMdqvh!%ITt?FUczA zOrC#E7xmcO?o*mWBug)4^wPWHV{y0=(p85c0T5fBe zBd;$(qvgL8SOM=)Jf>hY726sW=ZMag*Be8IBCsn6aW2+|Lnz7!qrSC4^d2OAVXLpP zWh^+I$0-kp_=lfxva?%8mGX01OO4%xZUcma^O|skb9~@T&W@lvSpW2LlSHhl+w&8$ zJ#Me>2W5U!m7~U@KVH)Rdso7vW+GHBs;Gpf8(_-&OllsKVUA z>@xpoO<_q?qIpKv^z7)L*51;(UorcTHTn>opBe#x@v{FH@r}9+@4iMaPxykLL00pD zK_)Kk0vNSniLA-BdF$SORHf)^qf7!a>x&FPoH`rN;DDM$rkcjLkX^gmR63RSSvl(X zbal;7*qN0-pjdPe+&5}bzvfyXE$wTceGAqRpFX)xROO4ry%RY_2vG!8z2>*Jzl`zRp`V2P z_A-S!ejbM9^Bi%aBm)_7TxNRN8YFP!6=03!Rt7EtP6~V_mq25E{Z(-vL>0@vy72mCeO{M~*!Vmb(I&^AW@tiF8 zH}$%D#{|ZIW~+9r5r=QkGh+!5a`r+)L{Ddm%dYc9eO{a|}7LYn*ao--W+ zY%12&T^i?w*q_#n%zm24{){O9LwRjLUG-k;N!hcA##RR!ZfQ43rUrh?U5hwWVNzdZ z&&KQg(YHOS>4vGme+1#l!~#1XS~pS%0oM%hl(%CEgfgpY-@D6`Ym zu0A7LJVPq*5p$;g%Ziy^fP+Zo`4g=p%F0#_^En9fu)6sML?1G$-+5>4JL@4NyVU9n zuQJA$R06-Z+mSa9Bv5rw!nmG1+gJm8Wl$-!W+WB~vntMes&~V8@BM_{)uLu4%q>2f zWus;v_-w_XX$i{RI6J)m*h>-~p04WWhHR6aO$lT$392$SZa!{$)xHOGqb_}|IHPmd z!Ygcf&i@~_&N42_?~C>Wf~3-ol+q2-je??dcXx-hqzovj(w)-Xjgr!mLn915bcet_ z{Qmd8yD!EUd}QVsW}b7-UVH8DT9nhf{aNl{348R+Q?0cD`HC}&*>%dOUREjbb6+RE+l8&V)kK;^ z44Ozw2imi6tQDbL-5#T z&83;t;jp^{%T3`=o(!b(TNae*BhTyz)Y`eFzE^Q9)l3YZESBQ96PZAm->PB;1}Zb< zF(7WKR@A+=H80_sm-f&crq6Xn7zMwMlvljd4(IHrsm|SKGnj{bI+O~|#7`8F0*$Tc z!iQ7oU!I@bxetFa9oQCZsay4oh=&UP^%t9294t}dRG-d``Hpxkr$NTAMXFslwJR(G zgWU3^!gt7R*|cDvjp@>isx%g{x^`Nyt4q_?lfvSbSi|{#4U3)@4H>lLFcq(Kp!9sW z^%fEs>ea!Wx+W6pVi#fZ!Ab8Pux7ew=UgcFT$&L)!gq0P-wyD^R}C+#u-6>MfSXg@ zdO5is|HG51;(eXGdg)Vh-lXxH_-3(kW<9LS1bpOqs?0Bb5Zc<8TOCOM*Hr{cW@CO5 zGJURdFX#JXWu_3&hjJFQ*G)*K>ECwufUABUEM*(~+++Hc*XB2w6_tpN6p&;dyWIEN zB3~=%*N*J?mEt6(xk*U=vRF`q^bEyXFh?Y=HQlpf!%Bv=qG zTW;j=zBP4^e?y0ErTUgV1v)vm4wmZj&_8M@zyL`D`bEjc)1U9kgd_}hT3&5aC=;Sg zs4g{HNaLA{bR^YqYsUVYW3rdOHfB5oLj`4JWp4{6=5DX~)-3J4U3E;o_vQ)Fsy zn|k?un7ix>`E?SW!PYm@st}I`<8iScD?rtLJPle!b3|S6N0-?HzR8rIvk~-F>CR*C zQ+@C?9S)dpj|Zjo1JY^G=K@;oYLeM-*`hDvveZZ`yGIl^evMvn`JkVYemF6iHWGmU z)$CXPO_~~yfa&a;KO?=d`LW>U)7~Uw#P2&jTQ?PRy&%en92CfiQzy+-8<}yaNP~_H z=Dbyl;?@^<7D7)M0M0}lP+2w%2@m0}jQjltW z4o&y}xb;*L6>wTP)q-(ZDJH@J5s8Jk!IpZ$s=O#~=hlo_#rfik z_ICvn(XVWi;|EU+kLUQ!DE8)h)jD^?)G=j@+I6UBW^+G;{^RC7<--(N)TZ=g?5XZN z)8)WInAU!!S|aV)QKmW2coL{{gm5bidF_7ta(W}PPt4tcFqoe)2sdjj zRN3}DQr+}7$EILS%k2F(_*S_j!EHZVO4t>lKQ%Ko;!fHYlB_Ua9mwPvZTgH7hMLl{ zuCF7GSaIqfk(JpE#bK-J^1Xncp+UxVHltRCbAcLUE=l|k>rB0n=L<}e(IRF2mixQF zQ>wV{j+MuiRb7K8c`YL_)evkq%z4DFSTE+eK^lM%grNgx87v6{FHJiw1!Fut#IwnJ zLl$DIO(p;uh25syFd>=+p1edz5gqJNul(h{rt=6EwD-XiP|f!I+iy-0Wz7h6y6EU= z2gtt@nlUp%$rKDufu2$^g40S|@m&eORDlm6cn7#uZpnV)xV+lQR(rbX7qfdcsV>7( zZA~ze;<(kQsKt+@V)ek^adhz!8yy+AS~4yC{C%%WXL1-cG6hSCZJ2FJ%olMvCs3kF z4+2s}e^)=99u*T>(cMpN;5(pHiIv)9)r~Y8j{8>4mdniY3cKj0ULHmF$p9e)$ayoj zT3P4QBs8d>qG5CIBLrT_y%#KON0M|ACEuOSO>J_AZ*%24HJy9^PaZCWdn~UPww&$TgRc>m-PnyC*M*(`1 z+s;PRw{+E=QHh-iGTWg!>6u}FrWbw5nFH_4u3BBwL!mD*A+yW>@QX=bEgl}Sh*q~0 zXs;jdn1mC;+uo-V6ox-`D|dfx{DH?G(wH5m?dRc+y|0l$Fd)ZitvF~fj!4Lx+uw19 z=^_c=ti*nu$52EqSF}-LF2ts2C5rp^k+9BUL@e9u9gmt4K1uN3QhueTM~?iXdqt(w z`_ge=?VoiM*g(#gO(NYU7e@unv_&9W8afMqT&|Cz=RQafI=2T!p8Ot(EBnHBtrd3M zJGV`y%9M|oU*qK+BDA!B!b4wCwX0gn&8kGvvE>Mr&AvT3ULfku+B|D27mEafG`j%f z3U0^C3DC$w2g|=;r`%S$G@IwJ=#zi=A%MtKU0t2iPSa#l67T$(ardNT_iY~An<{lo zJG+9Ubwgv$F*e4K_TlUgr0WnmJWo!>Yi=;5WJFb=Jh)w~QhK zPjc|Ei8PY!5;wbI#L79lNbBEYqkPkSIXKpcMD{vYW!~Ay(8QsnxOm{G*14pyHiwY5 zh4O!aB-T@s6(}Gg-59*%a$TIKAj`u0CTh zt!9iG0-=g|pB5cNzAS;IGcqdi^HV{{^X5u?Hj%caCyvuTo~-<*UU;qh)`&7~`0$x6 zfloe6yQ#>MPDlE?G(1yCh!l97!oXBipN1$G&!Sag+U)UJaNx0+8e}pe z{jDcZWDogjL3#5h5fVs761|b&P`BZQkfbh;Mth%hi%6W7RP5pQzN~C@?h+3`57 zP(Lt1tR_~J_H`BP4zVwnKY^odghjxiziR04zZlmLqErX|Pm=Za^H_A{qyE-)978>U z9w@Qd`p4Yx?{EnLtB=4I3~lfY??H)Pfub&y3+?E#>I@iVgLyxcw zO>~?5g7hi&T6=THbc?^f_kkMJ)XH6#ZRH!aX1$`T-po|)Z_Ijot*&?HCB^+s)BP!$ zYuEs5g&Wg{7bRQgTq5h09hVmU%TrGnBJ-=TJs#rP#RvVg)1qMJ?f9E0S5ko0Q|GbB zHc3qink($BudfIgzrh2>7A>trH=LVN9(<&)H_2p*b2%4^N8uVRx8}am&q@ZkB=Qoz z4OF&xs?{8s8X|?uedI7Z_G&w-@wcxZ*lFUc@5A9sEfDf;H&!eaSBV<-?L9;zxIC)v3}0hqB|4~ z9wK~V37;V#GE+RfUZ^w3#gvayd(0R2JXM}(@iK}AJ%_QzOq}oY_f8Zl;IoCOmw3s2 z|Hzye-!rg@QIw=CuhOdXrV6cXJ`vmRX2{sV&92E}ilL;+d3@=E@Lq9;rJwV7K+Eya zz@P8m3npy$p{YvJALr7)ZnTD6){2a|p2$7eS0xk8rw1uhvZ>8Vs=fDMMh73t7e1)5jEB852(nd_Jsn(iRPO2=LWhOcdIiqXhWR(Gh|-P0J@I4LXdf}( z#d^*xrq6lbmBsy0i>wYFi+ewr3Pz?+pB_niP4{i$DWjbZ?rU^f1nvhmSQN~EuIxx! zf|K9JLc>!PlhyIL|G5sK4A?4b8yZB~-y~Q7NQmj}^%)?YRJEQLmy~w`E}IJgs0gya zX83(ke&+UNRDG=BBLYL!E^yiT_9Z~kA>}frQIxpF*2+8rXd6oh5kzT?A3JQj@v16p zAzF*pZD7PhV9yf9TIs+1VG2iPJK?0FAg!&G-Lc}LW>0h7yR2)pnpBpa=Rl&wHwt!| z$M=`fKV8PBdbZTbl{y#_m>Vt5zbO$sGxncUYq&)jB|8=>Y%iyErCLsr?eU&1DB-mW z{7*{VD7Hl0Q`yp9Qv+W?MzY5F^W(2CFuussRNz2aJcbA@fK1bHP};ZWrl-Qe^jT-k z5IrBuH;9+QYKiouG<>)>FdpmYx5M!=Tuq|`i;D6d@Zf7ZUphJOY}?#}f>on;a5H9= z=|wBXaZNw*`}c^%#5y=IOAB>QNa~brXqc9!de?5#gENOt;m*+M;~{*`HRFRB0U8t=ChM z@|y6(HBX*@%4kM_5uX2L`^BgQ;0$MMqA@u8hSNo$J3B< zJxc`!H89sj5=JeKr=c-bSb~58#avwrW;~wv-1+|Y0_;h<;}zt76fA!l+UW|!RS1t- zmT#Ye!A1ZC%FwFzKPqIy$yQC={it9EYxz^Hz|A&K*yWR&=9Ud=Nl7jPHG$=iypOR- z3(qvZ@y-Q&9Ex97JXq*$*Kk3o|k!J z6~pGbQalDeIM|$4ve`E*XtZC|7$EkMKSUdbciT;wI3qEMW%UPWIMy?XdXJ9XDLdXA zZm{5@C$$#{o%PZ`)E=*BJJ9taIXS-Z+FE13q8mD8CwgYuoz!5*c|yvJ zv&mot(3t3P@wA064C1Xq#i^fGfRTVeb^PaOd_p5_e_5=c+2y8dB3xKOoAj+XI3-mL z#a%PYs%ySnUB!;m5@sAv;Qp9lZULNm68UgqY^kl%=g-l>`G5Z-y%ASP@%w#QhPVpt zljj1_8T)Cx73R`-CoV0Q6`xs=Lq^)@26^IY&C3Ii@CvMwj)Wyn%)Nx(v}Lq8PAhV; zXZO5LUi7~wDr{Q!X8LArxK%6NC%Y`bDwUGC!hg1oG}0AGG{=uu{0DcNFm5k3P?@az zn1v}Ag8I=;wZxgD^$)!+6P1h<&0q0|mIg<9JW&fHrU+jkRar6RFINkG96<6ucLU@aH>^IL z2YE!wxJ-fTBSwjnpF79un7woPUNekgkA0RQQl$QFkN{pVJ|Jj>B!e`kN9$Sps;KDG zLOD6=Z}wDNYf5xl{azwPJq*0z5DodjZP}lzr@DsnJ^~IKCO_W38{8oqjF`RerEps9 zCQy3K5-)9F@XWDs#U~8_Q33Z>0F0giTHxWNni3pvdw4a=^}ig!`^@GBc`Gl(Z-i_U zpIMrn@9|IR*}W}m`#vLRDC&glF+6g%qg#EVd=vr5RK1f$TAR#jr=1tI3w)eC_&z7o zO*Q%y`>oiKC(MIQO}>$sNHQvJpMO{XE$o;lAd5qLl>biocA2#IeCs~pa{vziDtRs z-@cZ1QoB*xk5RYm@y1DE%kG-Hu8jR`lPYgEYIoFTaZ6qksM?e=Oj=s0IpQ zdU`xiL%cNyN4sPE;#2^vF|8u~-v7eZ*4Af;B(0R@Z!=DVsrs#9Iy=xy;2RY*NTJ zh>PMD;qdsNMv~^Q?#zZd>GQ!^Q6*_ZUreDeSLk_AT1Snmfc;D$fQFG2sbd0uShvVQ zYgG2Du#GrguR~+`bRWKG;hndX?vrCP)fQO(E=#jso+F|!PfNm>!hf6dJp6FgWeTqI z{e}HPBX>a>SeII_G(KEWf!$G*>Wy&v{yX&wK?*e$!f)~n-+(JoPtVThNfrZPOds5O z5M`{QSn5Xq`wK8jn$_Gq=W5ir$I98%68|PZ;pTDrO*1wJCw^$o+0)ZnC`EwEVo~nqesOwe?f<|PtH9(Hcj63Go`onZl;bi zw%>Y>*AQW@NMD5?c>QsGxJc)r)nI4Gnlr{Eo&A-6Glxl9*|ow6&_2Ocza$4R)Obz9 zrFLgY{!|(8-&j~#@t%4z)4#K*s`!1r^s=Y7d~bim-WO`XOr+^?xxQkpj&c#wE|iVy zl%;z5ZCx4Z^1!6`!16&1ed!Td@MNs4QhU4CYtzx}fWoAcdI>8?qa3$E@XLf_!syS> z28I}+iKPF0Sao?V%j9hm_8=Wj{`9wJbA~0}7NXAA+WPz^RyPY+zLbw?)%))$4g8`S z99=C|Mq17z*sEGQlFy*)8C?F3{}xBP`^?Z_y7k+UjEV3_#xiK0Kaa(0`uxbk?(q$j z1b#-*_-v+WySsL8J(uOyGRQIPEOd{3uK8hpYJ^P9u9;!^hX0FaLW3Jo8v7S)0fE?@ z-?%cwtDjh& zyo_$*on%Lcc?DKClAD6l+&mn*lgRO&?QyfqQ#-}fS|D31ctqJ1wW-bAO=bRo#^v_) z3yo9?Q@zt{{&)_;%ss0Y>Zw?-2m7z_7z7aaieKV8w+ zcCk?6*Y3O`jZMmhiz9xC3;CPwx%3O1^4;!Ixt~i=0#qIkeK|*`c+t74#*bIZS4@wyBGz3qh?R)5ZCI^ zI)8olDUFBPz+d2u2aGH~(K7TV@|d$>0%LcEG*)u`v8mONdz(bfSnP4sym_oI~78}F7E>52OM=x&1|g-A8=`Z0S}~W zAfB}jf%5AkJNf){izfjVpFZS}{6At|6Bb@P#+?*VI zP>1z_VffXef(&Vvsrb>Zr*Q9IHw(K5;oqs3hS~i}+B)F#bsD4}Et><~i32&VCP4LHt~^n|KL zm66jGJ%qtWaTEx>llyuwI9j1%Y<-H$-5bPNuMtKv=*~#oUt-1Z-=CLzepqZi=eo{?-XvFb zK)q~8yONUk8g2s0LwU!-j9tk)%W<`n%A}bHO4&0{-AvOEU@8&pM6SIH_E85gEvYwe z!tJ=1D|$Q9l&`NFEbJ=V6`lpaMPJk5fra&Vb*au^M9ycq=O<;3Q9z(eNl6in;Hs(g z?2$~c&~FJ)SZbelGe-Gt-a&mqWnKP)oncOQyVme!=%t&Urs+*+Yt8Vmke$1@wx1<& zyGo4s5SS@SAEx2)KnAk?AkGcZgM9`$F6$T}1-OAia z&#q4I37hXWf?wXSg%X#2+FKYSOSdMEJ-&{)OU?jzlknAv27frl7)pKbHwTD@7-R^8 zV%leLh@Po0g~v!xwZuG&7$srvj2IP-4p{I`^;qKP1&h~wAcxk*H%&g<~ybn;@(17Mr(vFmE4A4ln zBE!tvvMH7A0j^A(O*>8&B)3IO{OV}=+01Glv| zhI;u5BY?7W?tC?Za}^-0djP;R9NfpKfme;($t^T%F%%Wi&M$&BWA6~KiwPx(sLZ=l ztGTvc#<=nE@j!O?wb10o?YXaiy1z(a_EZQ{2hDA+qC8mSX_s%w<~Q+=%RjUXBLWtE z6VH7W6rKP|EPypa0TP8`B442k{BeQ}JQs<;0*Q<37TN?Th+&jGR7L@}b?fb1BS+MZ zjrdAp$64pT-rhpM2j=F+yW|Cb(t6ZOK`Q=Hus4dLpp2-5Qbm?EKnmOTw)SMkx)go= z)*vX$XJaPeR{+}aV6oHr7% z>&!DV$NJvIazA3*&+9JE;S)Ydsv)<)BeSj3RnpM{Ynv7J9OuwU_$EisE*20D&Agc? zWs_(dA9B`1JGch(E1)PHxxn33x3rQLnNBq85dAdd%)4ABcjKIZ9*9VOiy&8aW z67z61)H!Iocwiw!-R$2m0+LQJyB@#UV8Za_R0RhrWGICXYgFt)f632nbEIe1D)SRC z&|wfT7>}gQ0@`^`0()weX^$(Vn9MK+AO({cU%9r2Z8N%wcY_+j`C>7kr_*WL?5=uB zyTOqiK#c4M?_!>bOjlQci|FC*qbxHxn0mKd8Z$hcS_z&f?c;I3x1YHkaPK{Z#E2e$ z&X*Sg%*{%lGZw(mAY*5rKLQ*kKmw!$!6DU%`l=&|?Y!8p?|wM1KuK&JcqtkzSYgJq z1AO4+97f2wOU$OOEhp}3xV^xg^-^YNza{1uYej3YcIgwndV5(tJyy%nJ6sZj(}e^)e$WCDEf0@Ur`j@WV5hkqKg&9XPLl#VUyh0SauXUc*z6^1%WGT9FH%eez^J zJL2h59k{I#fFF~j@T`Zkmj9fqvjdpp@abt{00yDrwH*4HY?V>OT!AX}8>7i8H2zH} zMa&bO1X?SLz5}M+Yw3~B44=Cj=oDqr*ZgmEEhf74li1z9eI&iHJ(E4|!muiU1N%8U*Hg@jcH8Zr{#zwY~4 zVq}$ZQBcE|%f376)lqOYob#2;%W)8b1!!z6o7p}RCQ_wuBE$-M9wiuCz@TE)m}MHm z!eh3Rig4+}yTxSbQhE>Yx{nh@FEJp0Bj%e+J#Mt}g^fI>yutLu@uC|oLlL5be+Tpt zhqU#?b-d>){qm{ZYWT|UTU^>sM3?`Y$DTjguX`1+_wQ`p8T9j>84T{MrSt82io7G^ zFZt@uA%`!kj-sNOWfe%0RttKC2y3?C0hY;733K&r3AWv5F~@734#j z;fgXUq>*H+2*y63ZP%VZC7uoRKK813A-ref)$hO9ea4B>hU;^&=;`?(hHY?5%qyfk zk_aqTjb+uaVu=`uJJNjlXJr)VT-fj!-)lSrYwE(yupySF4r}sX5 z&Ay^+G@GXIS$3W#!oe*+&qV5SzlH`6aKI_{tROs?LqpTC1Zl6WjI zsr(xU-_w&vJ?s_IOoW}B{x08a$UNZprxQvp^SyTOWi^H0NdBJnwRE>v%45I5fq>&7 zDu7DJ#!y9FEC*1>#~TCN|!@_4VvgNK;o7ngRU(ftj0# zAe9G5XdqM~V=7J$1onVC0rTMq2r^Y3KxJNl5BP9aEt1=Q5)%qwzbIfm3O*N^CjpE` zya$f5m7y0I8kKv_biR_f%Ce3pl{h6}x%w*VL`6jju6~vi+$nrs*U>bLXR*s_+%{d$ zFzGpaZ7tcN2_WjBV|;d>n8LN^Q#v8kZRt^!&E6TznlJ>nPZZ<|u5P@>yCBwM0x+|L7b|XTl$}8S2o~W=0lgK$Y8&s<1?cZ zOkTzxVbL_oDuDkb3qU+ak12jPc?Phrm5e{79TrNCN{E*(BIWRbkP zY_7^ntRXY!v0=ZER_ZJ|tD(+#WZ2dcdEp+`5Sf`_Q_bedka6KTlpcxo{;zGK*R;~m zNQyLqGqk~&79k2=h)p-hd4^g^NZ@jNHkY$RM2x2m=tK2=rq=x|xaXOXhs=quW&(Kf zoRh6AF#i}c`Z_R%C@ubEq`9**oN zV3~D}A}B%en<)ZODRyM`R#3J%R1~a{W?Pe{dAN3NxryNu!U6R#2osM)K~yI3w&Q2|p;zh$IVXT~cpF@WbaW62Ywj^Z=ktO|g^az+M^`#Ycsxg$;&_kzoT&u84hFcUL zX73N-UmqDpdx&+5AMJMBTxDw(CmcZIkP=kYSJ^$)Enq+gFyx ziuVvgiRZ(8rGmb~)jjZzM*T6ZZV@vW18A|) zaW;pRHdLeTmRYD3#qnFNIK!H3N@HLQe5X-2D1DBV)yC6UzQg4c7oYqx`aQ~Q&B%P9 z;o`*YRrU;wi|#N3h54nKNfOAPrKJqD$2HwxoPw0k_C>w@Yz~+^H<)Lk)nCL4Oo@@U zS$;aXCeC?2=Qb)D)ZkB<7pbdhR(ElMn_E}EU6(?HKNt#o4v}Q}GPJ=b#3AH5$bdqQ zhP+ZuMW7k+Yx;fm0HH*i-oNW;kfT(U3t|I1`4myHWoculKEdR+%H(;XX%<+YwxYi- z4H5Oo{DJR+LgdM5W_b(zV!%`e{sq3W1b1Je9inU1HI2VGY=lJ%?6KIZ+HLI&y%cA%#EVom+M&>c4C0m4#YaWRj^ zOmk1a;7%9$j|rRisMIkWtan_U>HR*7(NVK!j)ZHAGR)Jq*`$36xUW&VEJj+aGV`aUZd6%Jfp|mVM8PN3Vtyy^RSA zb@ukks|FM>+gx&s-A62IX1q5|63)0m?n0>4ty}z)-aNPmd|j4wOP$8XqoP3Eaem`% zcC~xSPVyw+^9^3@%nL}Um#__@)(sba;O*=!X=<8vdj`(VNQEZlLnEXP9BOjtsW9u# z0Q#bg)!QRUa~n{O?YTAP{iG0l7_maL^L?1-|A7H#PZxcwxs>GLLb8H;wH0 zG~A3RCY$bRDXBa@oP(P`ZA~*#9Kv z-Ko9a)E=eQ$b&6evM(XnG${?@;mm2ZMJ2k+5*nqbrS?!gzAcO#sIWKXqepV;+7JFb*S`ZJaQ)4L+yxtf|E*MfbhLFcitZG4Y9aywFEjU%Gk!e77J&Li&A-@)Pp)U2*3}gOZyR@%?LS`rPHz{Pn>!5KLUml$anL32#c@Rpyz4tXEu{b(;XC>jVcS z+%F(Ja68PsluzY<0?7P0mP1L=*FIJU9*SFUu?o@G^)3YxNEk4|8|PbDXccyBWyS0% zLyMT&)LXTg@YA%u(H5JP-}HZy!WvRGQW>?k$Oe-wq`-}J7(yTw$Fvc8TU|5#MK%{8AcOFP1C(y7c4*PAWz{sw=C0@Wf>nW$JE%%%~Ho@0+HRpt1cC#^hKWC9%-cF@}Q zplJj6JP>67jRX^P>_b^9CRWhYc2XiUc}WDx4}3>p?M7FQXGRX8%*D$6*?&3G%{CnI z47sOx_dzLoSXv(9VTO2>;2<*WcgyI^zC6PxlFR^N2$$_7hmbe_ZXCTLG{?hiyrF`8 zqX+=NGyUtnvli;l=Qls*3gXAIvxxG(q?9Q8nm*Lj9s?8v_C`f9zw1KX8^G?SjqOir zhn4Hrv21!yZ}{I{3)xa|>Ih_Mp?`>O_{teFbrCRC953Yio8u!j)|Wb<=RLHMK;E2C zW(l$}-n#{92sZfaUag=9o+0!#h~3XBisa64=jtu6Ktxf{j&PacZO5`-13u)>pFe@h zf+*eTrAk&W&mB88v<$wKlsDs7vzZ;O_IRDoZXUuq?uY>YmKki*!lH*iAGJ?0qW1{7 z`3A3d=;098j79+bN%2~Lr{mAjqqT&1Mcw5%MQK7iea52u_UtpGpdxN!K)Ey_CQPQW zcTQoNagbhl+b!BZyYz^)>f+}+1&GdyU6eQ4YY368eS4B{WC`*kxsR4@fvIwg#+hC>b^FhWH)E4iS-KjcVDk7?r z=~}a;6pzxj0G{b3iCK-s-ZNuZ41H=~55=T!xsM-TxzD-?Q+O#9YJ&Nu#+{g{1d3Ze zEv~nmQt=slm2o?#(Bq~|9a2N;9Pk^;0(X*ECq zN`FTJ8_9vfCZW0SI<*CadO-%zrLYQTaDB6AQiCi^mDf4gq$etBQIWRw{RXJqlEwR_c_;Fn-JLGhilIN{ zBr{7}toh>oRIJ*3!AHHGf`A10&9%OUj2w!>Hr*_bibfK@N{wu>X^4 zmMpn-`D-+Kbg=xh?6|ue_P+(7g<4LFUe{S&g}w&brDnOCA=2wFizH5qA@+L4|pKFOgUiwihNa7pmd*9&o@3Bt%F9uB}BJ{e2ZiW(dPw+5`^ zd3)ylPU>wR0xqt(~sz-B@w z*T+bBbA~TX(_*aAKLZ4jav(O%^HjYZ;XN%~;;#v9eacrr0qy}*b|eUi!3XXoU#Eoa zx`_g%zB-s88iceCK}z{yvS`lAzGG^-nD{D{-(g_;(}d7C7Hy5r-PtEJoL->sJ_6PL zHS3@8@$a6)30~xS!ew_gvuwq52Encj*cFl5o8@wO7B9yCJqS#o0}7FJt8p8|;ia|L z93jHpx$&B{NOZX=WK37wjI5eUI;bAxm@$i~f&@}H{bO%mr0lxZE`Q*q-L!2qS#~?j zM`m|v@00U8ax+IH#jy%fzpU)g>~r!zNlg<@VL5TD*+!7Qw8jM{1t!Q(Kuu5Ew8_O| z8KlTev7L%=pGYY4MmAR2N|9*8x7h9_O>m4vg2+aZL)*?(yL#^S&X;q`%GfLBz}-a& zsLg1QG{18lNZ@)Zfp&Xxh2+cUFEE8negk1`y6#0)IW!=nn*XGegarOCj;Ke))Q1(E zh`s@kz_9sykKPl)U%OWD4A*0SbjXxRRo14uyS9chi)-(=*H6WMH8w(-AZ@v%vQYI;?)8tO1?fIl?O>?F@?SXrh#P=rZo~1tkJ1hlKfB5zCT7K z9nlv}nOwAAbnfsCQ==}Bd*nXYYqqCe9`#n8jBl83Fs?*#-lQ7N_DC0?^K-wx)#vm- zowvQcl&QM-jOSO#SirH^d#L=S_S#pbkj^PBit26asFy{mUxGYCT`L9~QX#Aj;jSQ< z3})r=*jcq;UR*}M<_QQ&5+B-N&b7$DoZBn=*)^NWXO>Uz$(einYv!;#ylqS8HZ2aW zquLj*>j5Wv27#E&?6;Jg_I4LsyzyA)jdUB>nervRB

BV4@;FDf|y3MuQ0)SU`dX z8Kw7{inXI-FWwu8A4C1z9Jx)_5T1i57%rt$->qm~9Y{DLiYV(@(FQ0$NBw+y@3Rhkl6I1!b7g{k0aU zguf_gZ-Nk)@KLSOp5q@zX6GUm*-p#ijhjX|b`~PTe5SBo5Zn9+QVBBa9`<*^-O8@o z%I^SPnD4iK7DK|{UmgSE#JWMc#1}GNL#gF!bFi)Y3mRlTZHOr8V!d2%xlEgB3@a_u^5na1K(`n6)v~q4u!NXP*C1flkpm z{^LEj|39w<*Fkn*2w(%9kE0>p31^srC(tLpDvDp3pE0~r_2K*nG6efhVM{+v-cWv$ zLpzALyS)Km<;NiGs2B8H!aM%0b3y_EexGAg9I$&QM}rURY4> zy3qA~A#4*Hd(&d~lt@Z4%LTVORckYJoB#!g57HZ@P8-k)3A=^Us)8~tM(a`Pm9Psj zV5I1L0(@NdE-aTd1o*;4bmUN)-3z;=MQo|Rj4vN!5~Pr*jTvi*E4{rxuJ(ZZmh(r* zqFX`regx`eC(fQ-ry*~DAn$vE4X#%Bz%8pJptoZ?_NjB~JZ-5XBSyhuCq$CFnW$jKvXz2hX3R&ix%Z{_vwAxUO>F!@$n%Pwwg~MK317A|V!B` zVVGX0V||CfgW*z<{&M@xWtY*xD;j;MIvqw9vSgz2(L->-WTDFCHKO8Z(8w^w?SOYw zbdg?iKSXHzFy#9!V(3*jHW?4fgXF6yviAXmbHP9;XHh;3WHMOMF<7YJX-1)H^IUrr z**)XNiSpN=9fmWKS5w_Tj}cWxAm|hS@oe%M8kPMOVr|;X#Va}EwHc3b-b@epHKnE0 z6RafhtEnL4>&eycZ_izR)1tvaD#GfL(cAe`yKX^8C9yL?Meu;RKNuc`?5rPB=HO8< z7r?VO*_*3tyt_uYfGgTu6gUfCE)0)UG}+I4^|bqUdYI|9!16{&?W0Md#V@=CB7AY4 zcXv3u)`?H{za=2SFQF~yDgeMu~^WfvY3iiJI_)B!b z&UF|_+b$qYZF8~U9>%FFJSnDO?ds##(u%&^Hfizl!`ki8TyB7460=q@9dfH;W$G1k zHo7u=HxqhrBrDS57qQPK64$uf)GO+cJ-e@Sx(>c@%P?RqC-_w@MMIaz3P2I8K-w@ zZA8mIip5 ztwjnAKE8Z-HWP%9tn%=(ZpTi0lUZ9s0}3*=Ti3R(M|XCPC2A06YMJdb7N^#f0^U4N z#V@m3xoZi^?=}QFej%Yz1_D0_wRWX3l4=nlXsHDD$1=hbJ6cihFZ(1E<{e=5;(phr zr(0u5qTW2<)w^96rk!R!e_%lkGOHTi@!CkX0iA#0hn&0Gx^>t?;O_dFd#2y{uZK5q zf75Z_1ro+Wo!YO6ljX^8Z*YD(w(Nz2YJ&o8w?l#_nGIPL?WL3)tixddp?oYPnir$E8(1!8AS)_#8jcR;MmMzYlg<^Zb3Uqb_f z#U{=nSpIB#<90qp72+;zrLDoR%hSm3A#y9{#}6Wqo;Ksq!2RFv0pE<9SZ%l`N0Zj2 zZX0V&%}W5uQy~nIKBs!+Sf!ekWH1=V+@m@m;-~1@Q1A_TN=>G)6e+9RL`DTd>#2C= zHgtE^`CWj4l zy`Q5+r{kDyX}8YSO^PqmnWeOJ*ws)xr zIUY16cXC%uXRsOlf712%F+LbzJJ~Sy))+c}T(bP0>N|)Su=A#7z1h#b zcG(<>xg9vEJq&j7Q}FX9aul^B-aQYVx2eBK+eNLhEA6eAT-jV$wcNiwu046g zxfy2Q`%XwpagOJ#P6IcX7C;f)o#w>^hMZeEd^=C{Mnj{9zvjJUVHw~CaazHHX2@R? z@>n8WSH(n;?~)Nw(hy3ius%TwS1XI3Tv;LkB zL{IF~ULEa0to?A8C7SOk@TEw(bbYIrQEAI>hfBhO)46`LF|^+OC?Vv&KUMq=5S(x&E)xZ!ih&CgI6nBf(jDO z9U$1AcHvAc~W?p!9sp7e>Z3PI2_5ug~ zy0M4~vks9Wn4{vYa^RK8lZu3F?fQPsm6F+8%%OND&2_G>moj)Xk6j-L$3u~Ep^(?y zkOD(;#q(?zF~WD;ZRF}Qu`RB9y@9u4q}bQ4PXP?Ub|{TWJmL!leT-J=^}$-}4#p@2 zO-EU4`_}!&hPM+;Knl0^QJ2D&E8vGvtYAL_VdIn`^~(wy}YWn zXKyKlf!z7Y8oehRpQ71f-=K0qJh(lm-E%OG>)CC8WDUK$@XD z&*J}{_nfa>%r$%V-fOL2J@@kpI&+rc-6cb(hHBBe{yVBZakp{%NZ_RkuHr zKiQB=o7d_*O$&O1QlTPr6}e)hkdS^LSCb5;+<2Qg1pZdi=mjrElct@&G}%<#o@EeTpAgy#(&`^NKL!E zsdAf2bnN&btPXxv>qQee0tE_6RWi2mA8t4?Nt*&UWHvICGljKUIUrJ;K#Bdc(r}~k(fzk%3*V`;g5gafUCxBa*b~*GIb+2dBFhx5-y3X@Hsh z@3~9Kl~)n&@=a}KC#n;**ccl@>!Opf{2e8Um;g;S)yFTzf5r8U=UVd=bccEqy_M-r z3keI7)+4KSc$OTL`{N}AdBTv4#=O%pYOFj3w?e%qAKZkUa{ zp!n!F;E24tKl@}QmA##?{BT>R$=H3^4*P7Dq4sA$5XN8XjvNyde>7EcF`@qpzt>j z4-Iz$kFSe|uS=EAhB1kzr7}I7vN%8pDH>}^yKa!uH=NA7-=-c`Ym08dti8DM3g?58 zMMye4e{=iT+iTbMv1xTKn()~s*$oDA=LeBQNwYITmOsLf|0gBdnm?9^S$%N$CC`vV zC&zI-nc4vytOqYrjca?}KY46_<+^8gdw1!h2c<~K z!&?C`Lupfn;A-RX;WahGSAc==_}U%SMOsf`H@Au)%@u136aP=XW1Uwdbk$^cCVaVs zM1dvhLZ4ikMuQq_oOv!Hd#CdC&My3?A-VUUq;^w>vmIIFVw$a2XfB&!sD=@%-Mn;~ z0dSp2sIZZp*gj6HV2p*~N`8H1yWp=LwiKET#1HN9@dOv(HN8{!DsyOH#jwv>bj}+v zZAT1>8LdXPzi5GW8t{Wbg=l)>;t(SAANtJTa$rbWKI$Y5RGe=gEi!h?mZ}|hMSTl- zAca(4?e1be1)#opItu zYUNJ@?yc^TLGtN6adoOM-KY}VruG`m)m;3C8&V(7G|dt1xTKHSiZae@m)MdQwVCGp zqtz{?F{iDhk=G2muAc4b>{Vs?#EW%B{8eKX$VXF_-`(`IafPB9!81?uryT63giyrl z=nK(UnY6vlSCG*Wicq6Q|2vt-%bD*r;&NRDABIIP|g<@eLl2?8$ zUN^R-d!9cyEWR?iI};ycK!JS`qLYNVBSKsyMGRF5qUOx`Gy$X-qKxp810}N&sr=^r zYA83be?LKWqLL?DRl_*t6|-f7jYoYc<3z1mNO(5=)5d~I^fM|}8PrXnc?SA@|2=jb z35uU=-_r3wZJ3QZd6xIZef8pDqb*OJK8~D^i0X2J)l8-_!;Pz?f5__0h%2mh0d)73 zoIBg*ChEbGH@xc66&oA86>53=mbjETDV0Ls*q;!*++Y6 zNLxv(B=&v2EO_qHG^9CFo3jhJB&2*LnB`=*xw+<%q_@3+^+CX!A59Y$P9&?XMZ{ncmI z^v#Pzga2jkJ*tm6XG8Lkrs(Y?eE23ObGNog5=)j5=6!sFw?q6{FS9vf^~xx%&X*Hr z`(m`zS{32fv}nlNx~hdr^LFTX8MUpp3X|Fq zqTZa`GJIN^?WxjcMV{Hp`V5WhV6IqP*E@m*s2*wV&o&=zLwsH>s5Fw2n~6&avc(x| zS9o<6L6o18=P$JIkk?>%6#F{}ot#`}c*QUr#a|^yeK2|X#7wG-=T9EyXZ%3I@a-c8 z8pPkZ&CTM~E;ZyrFvEca4cc?H1rD80A0pt8U=e1Vg|852Jid%Fb~R}&<{z@@lVTT? z#N>Strj2a$JAVUJtW4;W2GRJ3*%Bf(9S+;Owkp*j>5QW}) zwf(SX5)fQ^VGxVFJ>geL?ilp<@!b9>XAmRtuf`5fCxIuU8wiwfG7Y~%YgOZ^<<=4( z1xV%pcJ1*#9^$S%fPO?t$=Y>JplD$OqXIz@A0d`xOHApLCZ^H9=vvy!?jva^B@0^2F*9WG{~r zXVK}LIP=RtV5Jx$iNDcSG8lp@oBxDb^#}LaY+#gC5nIR%P+kE?;O4A%M=hk%-T(}T z2s7Xw*jJ-R9dB7%QC zadGWa>wdXbWgZDrU5!&Ym*Doihmo2#&V=jx_*jJkOqAe;0o!Vif9d`uc7R@l_Bte32q z7o*15@PPLDvTMm){O@oUUf55C16a9T)GVkeIX$nM>xhCAfLMkFwJ-*uhFuU$QvDu} zI+&;zc#09IqvM>*VAn`UzK{TMAQ!B_1gwT~*J8WQ>n(Awg0mf$Hl6sYq#r8`*6GO} zYakbIejAq*CzCpBM{F(_`r5=>zVb6KFHzmrcCpSRB}M9ci4j+ZW!}28*{&2QGIt*} zxtiHBd-wKjrRG9f2hbSSvjMA1>D2+qt;S}?d#gIS3NtUtBWVGjFqn zX^VXMIMXv-QQJFJ6bg3_&3_rm`iY3QM}xs5>C-FDfzxeno|)gg?)3}QO?J$_q(^6# zAobtI$$tfJtk)}#LX2okDNN7|PBCSm!bZViQA-?;pOzChELg!cn@!fL5UBOxwfcVI zL)*egO&664w_v-t|MP9J!z#3NnXu7d#!)15dyp@v-(#(lFcdOk zctvQ62+_ca!HM_0)XS&mtlGtC?rS>N(he$KkiVIQ{zeSr4zn5!5) zn8}=dJTn|+@TTW2IA(qasLOWK2&DQYP9E&VU+M#@bc=l`6zzuxPE>I}nG;y^kSc5` zGj-FZ6|C-ej-H~_8nV%-JE{M2-QE5ceCH_4n;zgJjjmhjAaYwow{uR! zm@D!dKM<%^;tfU#j=ZYA)jGJ<&~PTlX6R19Ok5CFU;BA(x|_&Z^78ij8jTht=3*V$ z5`)Z(Rv$A{KlOoA&5d|m@wKP`<|i;DIEcY#2DP93u4B%#G1GlJcNwVBruBoW0RQeA z4qS~n&L=m1{G(Y?KX>)t3`VbuyW|OF8U0m{JQ1t#yG{ks>seOi7QytI>)}gD%Dsrc zz^M|{Ay(xLUeI84{Pg!vJDt}A^c1n7aNQ&J=2*ST3d9#I|8PDB%#fY z5js5qM8Oq$1ZMhFCDtZ)v}n7xD%NRy%}y&}(_ACAde}LBlN(YL&^q~*CXTZTqfCBf zj9XK~OwF%w(XGyjcjkXaTg2*WdN^7t(ZE)cto-t& zraP}Pw#~N*onXddwB|3g44fqiowh7<7jvaVZ0~oWaD}c;1y3z?%t%rkcINq$X{k~U zlcBMkAqi(WV`~rIh9BOr{Z!8!AmAOP~`0gM`01Ne|x{w~ulyeNZ^U6~M00t3=7KI?GUTP+up_|aQRN-A0-8h48gCzKpc z{%kvJKEsCNW0;Fb4A@tQSS1EYRga~6aBFBMuY;30hhP*&Sw(~b4o0Z!P5xcY0TdDy zpZ8U+x0yLav2mElpxtQ&pTYG;zi$OuA)`{=yn-y;m9N4RYZu9{3!7-!$YD>q!6g6Q zs692&2wSPKc4S7!g{ZXZ!M*ovz@CbjtGVx9U+T5xBCHy15+iCq$=?=n?|Z5Pz~o#vtfekz$JW*DdEfi?vk z8E8{Pb5y5XSgu*n(R2Ig6V>XMsdqhy!_h0;nor}^Z*-T5_(2PjI4bVCax-6V+)mfe%kux_!4!hDR9T~fwfwz$gl)cg zWxUnwuG?-DyScj0x&?K<-PQCsVexj~tW>S-p1f-Yh!yB>fb~WHmC{{K1z%$|OZ@_d#y{e)IzW zD&wuPhm4xQ;OHuqnU~+DxN!!_bmGmZJ@{<90GgDt1tlos%f#{A zbR@^OFV&r?|5m_$9`JZSUaxGR;g&^)$CkmMz<$UlU3mwsulgkz$pK*TKK)=cbBFgJ z$nbxsz;i&0CaWqPNeSh|_=Iunj3}4S{|8jIWu#jn>+_85rCV!3*f)+OYUIUwyHG5< zKuzU+?H{{6RMZ7}H~hjMI8+BJ@x0nwRazKnHf!KIW08bsRhroC7&lr|Hoa3JVrpj= zNcI_%g|WHea%MQ(M6cyWn9EGfb!B_c;N z`>e6ueI^d}#Dvn8qX~-1Q)7w!XVw-15V1egL;2W~ppt$nVPnI^qNBTWoP)Z}`%zH( zDXgWEWpLx;kh5A^2yp~91Ew2To}J4gpT``;Jz$Y=`v&^M+Vs>nK!(wHBQTWUSjBut zndZV6#$7i!*K2iEEG{JlUuT-UX*f`wSM$xgt4wQh92T;p(Os{}9i2Ju#RwRJNB%RU zwB_D;o^i)Q$}_d0dgXQxk4E~qKJLI|)DF71FoqEhBN`3+$7Y12U|$3M%*RllLBCZM zVklF7&kNj@sQ6aLp53MmTg}+KTuSq{1HL=W*)E@6Y8SPkrPP-4^lWju*#fa6-!Mah zMo^>2oAVV0ABON2_k4Eu2SE#GyCe9fj?d27$oPYnFp%BrI+8uq_H;iycD#~?6L|+J zP3^@L*iQrxmR%WrN*ICBr*y|aClP)20$;@78p=C0s{!P`mS)h6>9{0#dcfMQithDu zi)!FoGoj*-s{kPpUQ@T%*HLev|6S^a(p~DBvo%@#0h2J7~^``H>s|?zs zQpFHnd!^|{$Ot`7=)ea}fQe@?F|G({N(<)?@(?Ynv#!c#k~b_S-%hZQP-n7!1JirU z!UB*y@>e$C`IG+jv4)npQDV4J$)Now!$Jx^6()x@H4T#-$bw7T_5|THI59A~A4{o3 z>~oGM%QY%hzp5bxYrn7l{`tfFM%P_$nXX7steJ=1k1I8|np(~M(1I`D0!{g?eqUFw z8iZY#tZkB3vsuooRsBrzRu7jJ5yp~zFdeg4q>4j&p23Eb*x3oQfI8)fvY85c!Aj>w z9A=o*A73V>bLJ^^e%0S%+z9aFr9yuy+=gm#k7RPuL{lB+d-Bc4<_&k7}!Nl&c=QYB~fVnR!TvNQo?Ch*g z2n$N$mE)&jQu<;k-ud`h+=9xbQaJ6F{bZgpRXMK$@qAG&7OX`+^3-8qnoWg&AOH-w zB=H&O8MaAbJ(+f%%?9h*9GJr|U;)P&w`8xQq4NesgzB~(FW#$z!R*EfRHrE8I)d5C z-~%vAw5MSq&mOjA1AM492%6>SW(*?0_7WEI-vk06>nA%ivehe?>X8kn{?QB}Xe}Pz zZUp$@QMQ#3wVW4pk-Jd-p9FzW9kxQ0Uc!%^KP-Wibk4jUrs4kL0^;cO9$T)WII|vE zB5m3TuU8?STJjTvqULh^U!Btj*Y)F(%%M^EjH>{jgHr4$TljyrqbhV)9;K6Wu3Ih_ z*Xy_tJud(aK-;Cst#x6lvyHM4hm@l@*I#D`Rh)*Msv4P?Id{&efVt0*Xm+e!88MF! z_Z{@H$CTG+687tLXx>8JTnD%7)2H=gkJ3!Jb2B+*oBlkt`LpxLBbtX_v*`c*`^A1a z9|0Np(eA&w1UdI`;0$R;4hwNUGU45Qje^NZRIahw{epk92}f;#L3QZV-T5KP;Cm*l zlD?K}23oIIq`tC#LVWLPlnMJ!*2rN87jN#QV(7jqfuxtLUbSa=rXM{-2_iN;o83=f zZ9V$=>Wzge^;Dj`bCXE7*ac~`Q^)Roo@uROix7vx;=d`b0fwGe9v=^WpP|)-ECH|M zxaQH2ls@6EycDtSCKaq-4bl1<|ZM>VdiVdKFL<{zEAnG>h6;R`si-v_%+5J**v zftWs#?e-T5E402NHH6I#4bU(k1|(4Bx;uAjE&`b*$1ZV{9HYa7E)jxmS%TooYV|&U zd4z%^Xi4eqh0_Wb6MH^Gc+u<)rep}<) zf+KmP(S_f4vkX2+W%sxP$rnB76ebi`#YSHbgk7O<=JH|bGKd(+&H%)gU8bpJwpn8F zA=~#OyMZAqe)+zvetSNbbc}E%qV`XIHKSS0{tdcbF6p>lNDK+A{Fn*K6=G30{Hn`39 zOl7vx3IdvbsxbriM-TLf&mZgbv%Bc)jCUt0D(QdynNBa6RSEsAj7`x0 z^Gv@m`JI3K`|W?Z4$Y)#TSV&@f?4!$imC|aMfAS9OXJ-_gT0i!LCbfS4-J{DZ(t6z zv0*}iYeyD+!%;tNn86pnA*Uw^;6*0c`=0^K7er*YPrnECO|O0gS2|U*L7sv(&Af(e z$e<(W)t=!FhYc3ZZI3GD-Od8T7qZ;lP>(qhBg!Vh$YiaCI1*@NWx$n!3EdDw(#QH# z$K|>x235dV$(lS}8QFk^bSVG23I2HMxadl7JjuBnA*y)vu#rdN!$`b*UOIOBVBJpY zEf_p<>>coTh^l**cfQ2eM6dVcLe+>iJclA7-Me3i2~G^>n5b`BF}O-UTjf`*Sr7B7 zt#w*t`gDpCHo;It@tE-Zvb)}$4|G+Oeaxg6I6CIg0mtbsOxvzC_m_iMD8b~21GN?W zy;iouKS$J&RKs8aCkRZ$TZSiC^X@aQ2EL_!y+0YAgaoJj8%#Ai80oXL!j|%9-Y{@; z*A-!8t56ulwJ6C3M<*8q&Bki_INYs8J>#;>t*8%p-tbM$6cjf0 zYq;wi8TIfum}aKSn~yb*SPsV-A|&wa;c4KXqLa~wt6=QOn`DLzt!)Ey0T5j$rK7aB zmA<7dIzwPkC#ZIie8r9bxq*N0bq}7V$B9Aq)I})V;}YSbU2Op_%WA+f^MVE=MLBbE zqZQ;#mss7(Z#(LABpjFJQz9vyd0e^%RtK__y67Biaji|=cT|42Iefev_Bjl;7fKK;>0LzKfuU-tGudQ1bY301Z!4S1@XarQk`KnO8QndDUxeE~ zD4o4e5vEA1G1IsH<-^?!>|dP|B4IR59rrF>FnIk^D>SxpudibPNkiuIW+!_N*NGct zDf&~Q?{W&%0o2Q2ct$d5TI%F8-;$KKUJ^=26Mr3KZI1NH)ar5~RhaR6?Qa{x20QiId}ZuaH^VeI>)lI>CQ#Kc6vv%1FoU{vU=@6Es1z-)RH!oe5XZFC(v z8R*sf=dw8ZcMqTz!QEHAqq;k+8k$a)D;?EU?)zW`xv#KxH8Htr)T0;HN>awD(imVi zXd?!cCASfej>Dt(qgKVWb~7_zmhR|!ulULQsCnZb?Z_3GMTBel_2E+U%Q*zJUPDiw z>8T%MupAUY7%3DK)Q059Yp4tf*y50HknA9A@s1#@{?Q-Co*X@fDK-@jB8%J2XNh{H z>Jvxk%b(v($PgyB3W|YFNLcD=={ZMbUHNjpV%Q!et0n>QAjdpE$u}dT0)P+c=HwEOKie7Ww9qY7C|% zM|jZ3u(NGls|}vSemUtI~l{OM>i!jTJM)mlglu3pj$a zN<=-5;L+wT*y(BJ6LRVI-|Wn3?$ea^f8OWFpWB7DJiPmS-yW(wRIxpy;G1wVmsy(R zFZJIVg)MCJ=PbK?P>Jj5Vali_=3Q4gn`gQcund4#`Ce34mz0`Me@8KTyQXKNf7R1y z{>@|3WjHK$zi)@7M+qaZ=7CTDHpAD4g7nL=0mLuWw&{-xtH)wqqAQqq?*_BLfqJg< z@OZlR{&`Vvs4#!=;srop?s;IBu#@asPG|vmM}gRLC%X(GF$<(I=)&h189-%@1563aY=J_0 z_zQF)zeYMZ4%(2YUXBh5zQRgnp}7uU7wuxiL_xWzP+K| z8uprNG0Q~#_z_67j@E!Osxy?#2T5Xe;>nDbn=YPp0Yxb^>$gCm(Shr3~ibmUa?icl~m;8{wJ0yrOhF)BE{>L8TCV9N@C%CH~bE3pw&YB^9;xVZE z`x*zk{P|fe4&4@leA(wmn;CzwrfLUcR~}@nS!$9uAdD^4adfttvxkwlITlW5Saeql zxm=t#(h0c2J63SAM^4P#)4SQG1*=vcOsgqbeH>pxT&Z8S9^s(CBbGXd#Zz+wW#7)f z^nKUlRaKicS(CPw@c6T-#Ox@nnk9^&c|-Qgr*Om>JE*X#R7N<=yv;93UIPJK41#&fZosdiWp{U1(0tS)PJfP8 zEfPhS`31}Cpb`-ft9N(m0g`d`evA2;1GjTpHX@u<2beB82IX{xwV zQaU*QyGqzQ5|MHho`kh`PKv|55-*nd6!%z_OQvW!WD|dPqkovJFdf8Y)~-Q7P}r zcLn#^G4901fpv4)BPQ@^9WQS(KIwp|kqSvgGkxSB9!LT!#4jk8REN{a$m)80RzUIIY?u*iR`% zkw&+73pSWGE17kcWv6Qa`6rPoGw7c5n<|OEh0HXpvim3E6Q%ZB#Eeij?$rk zaE&$hl_X(|MgNbp%Ha-cdJ~kGIOd(Vk~zPnrwf>6eU-5*EhcHtRu^rXUqZ;ghvYtd zST;3Fr{$a5#gg36x+8p5rBK}Z)Ksfj2Gmi9a{#rSJUVfe8kcI}DtLJZmErI?jo{2V z0n)0C()yL&4*wSpOAdWnw|uh$e><-MPA=5@a-PV=rT3~QJCrnJb5MkbYOLimpvc*| zPaOL`xkYg#siCYsUJ^vF0=}Zmm&Y4fv-3Y7kY>8StZ(ne(a`k=p$wew)7`4WVdJ*I zVJdJN0E{=Vj_m$|uk1Q(K018=nr+8PTW>NUL$;Hl(M9R zgamzXN?zt+m}ID^Zmp^Njj0Ci1cO?@Ihopf0&IO~_45uuz_86!JRNJ*^r7AZ1J2hc zFN^X8d=+zQaHZC0McL#T<9jTp$dc^=0^}kmAH)j<9iv-kgo|oBaeh&?WU~Hu&cGpO ze{RPCFjnwM6Y~as0HQd%n^-r??yq|mUIyW85uvoUu`4v+sGrOMp{kh__TQ=KjHjzfUhbM0#tR9xG4Fe@SM?3DwMPnH&b?TD zICQH7bP6E=JG zp!Vc4K5uSeC1gX9|E1SKE8Et?e&%(>c3Ys?i1=9(&yT06Z}c@K@DMZ-KIBj`!PodI z-3Ot)0B#tLr|*tAJWLcjxLMi#YHv5hf$SELa$}?-R(_}OSdA|A zo&;Gbp__l-+|P!k=o3Fh5^&)Ht{^t;hm_yMs!B99%!BUNMeiyUo3Q*U0E`G zZY`K~>WG8I`GDd(D%~|Dk}CQQEJ)CD=p9Q_t4AM6U01zRtVC?RB2tH55pV`Iwl|bM z04gNhrsKQ|w0T_e_&Mf(c1&1UtWTLOYuei`EG>m0H$0fD1pH|x*zqHT50m650{UK1 z9_BBP{Acs;ztOJsQ<)nAv>on3%e5k<-;Vd-NN%mw>^i_K^_gduLLl#jP3774TS460+q1u3J6Mo3q8dZu z{?PrbgobpCzatywPLcn?-g$zyybY9>811*y<e62&B@cfH0)!r?ygNqU%$O-!#Dtw{NGT)F_sW+dD+UAfoW4e!Dluf@g9 zQGJtu$Ml2OpV0pC(NJakTz=T8z{@X)6uq&wFLyNlP5FRErV`OjI|H_gSEg_LM8h2f zjDHu8pd>e|mBHNTe)v4fkrE>qY3WR#R4&+(lePoBup@gdryod2+I5htkKzi6o{AF_ z9c*1j`rCI@8Te8P*sUC#Ikc<>7!Dx$nwA^*3SGD#I4S02o=>QMf_L|+nTRX;RMv9Y zPR=k=COkFwZE-=cq0d%>KJJazRD^sB}vI&n~-C3`q$bHlz;XB@B0@ReQ=lCEH#Xm z48r0M2;Qzjo@1cof#;_E@^X|MTPrr8?6Rp%Yg28X4zl%uNBW$~be~(&x<5sImnW>Z zT@EiUyKQXz{-OiMOR;;7s$e7WyM9NYZ!@}$4fH>^hsU#GWk#)0>+bS9ktN@!2eZ9` zk2*~dv5p@ICZ#2sm9XY%CpOQ;rC0jXQ&WGRV>{O7g9%Zmi!M~QbwD}k4aZ^n(H?`( z0qeHxLU$?XjsPe2PEWl9Fq<3zUnKd;{p9=D1h`(Zhn>C#a6S^>KMb?^L<6o90r*70 zAA{302uTE9dxTcCnGv8*CT+8}-#h9$j|UePIpp_%&R;r>1hU>y^+Wr^NOQw+iJj4= zd8T%mf$Kamps-p`&1i48lxR^oMOHq9>A1dq`jQyz3poaim+pJ8z!baxtTyAnM%Zpx z5ZSY$TgPi}$7!}vRX+Sy7+zM!dJj~Yy{%U(p^o#7=8?qwSZ7iRJrV}wIkt=i(Usn~`jvQL9-@4Fw=X>8h z;I)u7G2cQ+XFp~ans0G^J~?VXqh$^XP)ke1pI37{rWO`L0DNemvOw2kQ}NkB%H&S? zvajFMssZsq$m(bT?6Z|%5A?##i*bs1yleG!AYnqS;Q$5*V7I^PwO9&nq6(JK%h8 zCo_#5=CVdD-aE_wSDFBe9cHb-+C}wC6^7dGMRmd3=W1(Mdx7$g%jFCpcvQlsNucCb zsBkxc%AgA#Jj`$3JYFq=ZCedqzdjl5Ju=xr1rEFLE@KDdF5%rbVl0b#|{x< zAGDCqHlR`i3nBfbw4m{dM4{3?wI zKQD5mvPT&E4@1H^SMj}^LZ`bk#^V3r>tEMkesz23Lk7U3Hl~OGN+bC;>ugLUjln>~ zqmm?!vZPYO7;W^>j|q=8zRWM2Nv_h~wLq*_bP=WhJ)%{X_sY#nGL zbN7#ypRSSpP(faV1)Y@{s3f<5YNLd#J1M<2osgQ`LExi_RyvTJZTAhGNIiQH)EJ}} zy1^? zh2HyhYBt|G30;nfj8(cYjR>6y{@}GC1n7Bm2tJ=JF?gHRpMe$bAWp128xhJ5J`QmDVViRVSJw7 zSO5ZySi;&s(kFYk& zcY-KAAc}}@ZNaJ>9&*2@(JpRc@&aPGGnV~>{Y_|qagFbia^hvo4+rc^!M_mhsmkqBnp)dl)X}##al_)-e%-GGmr8J|u2#H2>_`Aq@xK>p_bWdAxznupOHiWMjO%49 zg7b~mc3-d@JV&lOd=|IkVOs?wK|{*A(lmt_D|Ng)N*6Z)G@ z=NZ#>F{5)lCV{>Io^D?PgW$&0QQ8a0=RK`@pI+<*{F4#EBg1|5jnr5L#jVSXN<&t)l|mSS8?Qn2ncpC zt>B0N>qEOdZ%ZCTJptt?Sb?$uqJKZuozfTnB3Z_$^`2Z5OssPEy7(w5?Tc2Ox2-i9 z2|H61*-lxq5Np~<7y1^1dak1^>?hdaHi+5EC*pbBoIM!@Hg(5bKO?0Wvk;SY74}!V zLdS9alvRUsU9CR3sXYZ`s;i;rL-Kc7{9LHYoUq%>zO($IbH1Lxmd z|CKntKO5ZunQ}L<$=m##o0kWVNaBn^>U+-t=nTesF`WyCIBnpspms}@;p6!hR7g>@ zc&=Rm7sv2c3~?RodEL*2&hPr3M=AXsP`+u^4 z_B4OJxcV4PWY>ZTs+#dq1Qz8m;u(|NFF{e1T)pHPt;;VPN84*J38d0|>i`Sfd4F@} z0IhBBXz@Vjy1pasjdVEuo?(*G9`p8=K3d2zWU~$YV#H7q*M{Kse>Y$wLt|lKc}7%a zcpYlVX9&%69&SY65$e0I!Ou3*Ts8gKHvK-_34BJ-h`EtQ-=3iX0utR>a^}N9K^t_^ z#Oq(2$D7lBDi&|RnfCyVGzk_HBGAffyJ|f$+MUSDR{>CXbEEUm=Ho~EzBY-&53kk$vZZ=$%$UgM8p( zR_Ei{9=u@8@+Vxs;vw~ay;)yN9(*6q$oe-Ox6S}E4AXVlh0oc|ZKJIDF#zBuij%qL z#4ZCs&NTP^%4jF%^q|71+y3!p+XP&i|Me>90zMNifXPJX-_Ae=Le?GtZZTd~)y*il zDr2u0eXZqmE;N6&b2bWcTfl63c3^o<2v-FL*$X2{V^4S*7O52VJ$Gp~=Erk~x$SL_ zR6)jt6kz2kpnhjNi^OZQ5cCZ1g4p*Gn3JbpYIN&8?Hd>`H@S?d!2=RV*E0^|x#Zmp zt?Tx_T8;p(V$=@t@6VTO-3|7KA>RLjt;9{q}NqK2jCPDvi=I!idaR!b+3uCEN8m$*Jj* z{0;u@k7dxIlp;-an@iTZ91FsP_(=J>9Ec^95e5SWzs2s?@TG3Cp_2xtR}8zQ$d@s+ z1X)f2HJCuFsF;W(aeEr8W)|iaic}ZfRQ{q@&y_|>yf{_UPt_Kas&f-by3&E4tobkE z*R-6U?rgSdUFU^=WlNOd*M#+PpJbVd&GLMfMK;b)XW$B%eL#SSzWkWr1l@Y|InOqd zc*#$V7g2B4{^al_uy}3vF)Q*gaC;w^-`{U3_4B8Y$~-t5>ki-e51(jr&y zp-S_7NPc2APBX<;cz)F<@8Y#ux2K{-qh9C{3`9lv1}j6r&3x;Yi1YV#g$Hze8)oN? zmng1`DU4VLE0tB&z;a=6u@qx!AQ9{8>TzB3+Gf<&8yV`k=}{Z|CB>MBi_w=o(6Wsb z-9^O6&)8n!<5krslu zfXN6?qXxK46_}#FS~*rCsI0}fyj;y^mjaIutXAN2T7?)BQr4fpapW-AdDNXYkMlk- z-v~ku-n9LhHqpnHD(lC-_3;{}tdF5kuMvx__~Iqo@m<6kp|DU?8Cw~^QGxR~m-$lN)90r5{LVybV4dERn)t17NpikK#m z+Ce?o~Ac} z2TCJ~!R{%Fkn1}4(%j(Kd+4Ft&S4+w_BZTX8@%4^c$Ed!BPe`8AU^J!z%b2|a>n{& zmR0YaE{4kWG}%*RpxY|JuU|QuGZkjVrFC)27KFwF><2*>x^9D+4rAQ8-dfWYWfEe_?TZqL(<`x>QZ`csg}d@}*b+0}JJv3?Rt z7S5{WPaaNxna0w0&_LAY^>}qskm?(uW83}$(HE-rFOQqo#-cmZONBArsZ$DAtZB{7 z6^UXL=Jz;ox`8JXkAA?p`Jq01aR{NY-BW`s2q*>=FvWh;tT5Y~PBzgsb zF5Qvwu>FX%<#vH+%V$DsB92m7k3>3UBtIpx@lhtFYIu3YoG|4vf~}Q!@u;r;PecNF z7RUE4;5t=l*vSQIUN0Wc3k2=^(h=Hc%T^8F-hIow$G#)8kIjSD{OE=HHg%s^(*&ou z7*Gd8HR^Gg$CbH3#=#rG0neUr=zHBw#_qHw@1&SK5T;UjqzAWetvsYp?#0cKk0H?u zy66#;3W7LptHUd~@cv2fHWY7fl1P&&7!HsQj!5DyBg%?4aFW+A9xa$GKMz&^aUi=5lNOHZ(vbnRxTUkU)3{bIF<9+Szxe#~c-9`Hroomc!X z<3$SnEtHrCI1bOP;nowXhHE}ro-l|Nuz)|Sz-y09jbP}~s|RjnhL zup5uOa3D5b-1SvZ#wjee-uV=^lq84U$Cr36{B?k8A2c1?S};%WYo;>)&T;@C<}|4ipw&DX}%DNZr2F zpM_6f1wbJ=8WDRxUr8LlUO@EMS5dFqf#)ci0$%{&BK64OUX(ec+fP{cvGftp>Vbxm zR5l$hut>GQB+ah=?8pI(8*`?s58f-=Gt7(dnA}}*Om!fgF}yc<~OwuJ4qKp5K;gpL+3hF zf^k_MmzMMMA&f`sR3u9s@wbTBp~KAAcb(oRP%Rn=*69p!MUEb$C>8>V|9MejFJfmb zsC6k-n_$yUR`fMeW^!l#YcF)xkYn-R!q)<t8p@?VVyT2$2kW zxDtdU>wgH9neLZR{kDEs@3{H4W+ZDo^My$Wg&)U1+1NGN@HZmVpP4|J9emHzi@7!Q z*$xaNT>gdy?sl=eaBCazp;wnD4sV`18$K^x$xOa++JJ`ShCIB$CEI+v^K&^ebxqdE z)5zijtDLFeg~$u;GKm1awBR={uzCjBG!jT=-(nWAwfeae!}qgF^NSuM=ih(meK9p- zZp&_NcjaPyTZyBlez_flAL#YIPv3t#B$>0>OAS7eCR!xh$ZMCAogE;D40_s59v^CY zMX$wieC*qYsv~dk54*&#URe+qy-y;sxF>w$FSu+dkjdF66(klt*XMlp%_cyEa5Hp9 z_ZT&hJAB9EI8`*n88`Q>(b;#&Bu0-*XGq6UUrfz#1c&;s72;XHhn#qiz0)e?m7c#n1f=nR63&Z}Dy&h(ME3;C#kmbDBb3$aPzBhJy|ACZu4{b!I%r&Gxm$s+%cynJqc&$8g0}70Ijej-H;d< z7BRs=bYh%r&$k9rnES)%mZ$VTw2r>RVEO0?W~+`TXf2UmGnMS%7kpDu(F?HBV7VJ} zUA|YuLO85Hz!`clps2-*5Wio0F4zI?`Q4JRgT1{|yJYayNLit0zPGqYra>JdihWsn zFOOeuG+cQ$w_xj${^SezO%l`m(!}Vhhm#5{~*>K(pENGab*2}ZY?}*VihhU-?h5KJAe&-{wO++J#bM$JZaU5q;tyz9g zZ%Mzq@g7~y-(ZB42pE3t)eCORFmgvD=lSK7-ka5%v7nT&DNLG(>ao}<@hh0}`OF6$ znw80-Ozlg>@65jwANAem=UX!S7gN>ghb<(dNRH`3Rr3$GJZ1u1N#j?SpsG?te7s@z z6W%DOrhV0Y7K)md4vWmu#Q*ZIEbrRX3>|(X%ThZ?l6w26<{K-mSlt`7O7ieTa%gDo zMiMow9$Y(83b%1uu0Sb4tzsjxA0~n*RB zJO8|(Y2DL742O!T1o$bRKTdMI*?WBsJG-HE_DpDHFOGZrF3hB+{WoSrjPPu1I;MeX z-lPK|wvCTXW^%suAv!`egXnvDsXVqgSZSfln9%@-;+^sK_{sM`-2b2 zn9Ij`6*UuSHPPD8VCzU9B&{yP)uENu zasQH{@zO~v!~W)rywLHB*YMZmeEBan5*XqTzip%Ut^Kk!!3#C@o~5P+=W)@KcL?(g znfEu-G!*Ubu7p-hl*oI3{(SJ4&7%u9`dJ-I^J-|j+T!Sg4h$qXdOz(T5g~J6Rt8~l zcTdxC|E@EDYk|*W3Y$grjUHAQl;CUo@^VOgW)GW=pS5-W2|u=y(t95#cq>6eQ@hMk zMJP^rk?|;@@`gkrDX+_NSge)9DpB*pZ-q_8pMgya6c-{xrTF(B>p7v;TEXo(W8|}Q zp2An%UhP^i9HeT>YD_GHxc#c!V-zEYd=v+nZsDK0|1)T| zw6*;u{Vsj)fL(b#-RSEK(4l}Co)WBr;4(ZmUM%@sjFcnc0vK$?oBCJv%M)Z7 z2u#iZN%vX=gkEq?6}d=n6J@(NIMUF9@!c{tM&kjb*GCvqE4~>3TfN)Aw&g8Ke|jv~|C3xXbJg zV<#0=8HZGSs~~37iI!d4r<(qD360&py7_>k;u%CD#vY>dRcl>1mb)JRw(Th{y& zRjK8zZl(hwj9^PPIp3<$Z=DEvpWUxf@|8QeF*;Q#H%bos6>IGLZu^2Vg|d-hCLg;* zB1gmSGEgD+pw961EOW{}?)CF%CuJw#dYH}B74X6R=$UIuRB+5Y87o7=RT=FVxl0?G zd-amfO8T6zf4si~(!hL3?>5}|Q!QtmXalF^qj zVfOn)Ou~Q6r1;n0rq$P`UZE99=lm?_{gJ@BXJg0VS5Z+RsjW_l43*Bwtd8ZQmyYkr z(vv35(aq3pa>}!kYkJ&iYP{@nCOPs4>!f^PYP6;61H>-F2FJx&Wt)Sf3_}DqR&jXm zHr$bvp~%9tBX_r2(@#;;ycq_Wc6N!};2|<4VvXzTh02a2Oqsz^1mEiM|M^2qa^y}= z8e}tReis$hTLZ$A@=fMY;SU1w)r)Jy^TDCOeUG>&z$)7>CUm-tBp2iigS-ra%fM@48Oyd>elEo1F!c4%a_Q;uba zh~pi1ek+xE@c$kgl$k7MN*#D2b?r}-z5Fb}e-Dqr1p(2d2m+7SuL_I9fr>h^5#)+} zJY5s&O=vDfrcCZ?L+Ad{Sm2m{F?gPBMYQ8o_v7g#8!|rO?`@Uz@?q7`x4a$d(v(Bu z|KlTGX+JuVZEKQPX+kLq<3rm7E$;%>JriW|hR6OEhEsbJQ84-4omHfUegGGwq*fZ{?FO%kO#*6HD$3bn$QL~*ZCY<6n5ZI%iV!RW%UEHDxk z>p^h9-4G%iS-{j>>4g=9Xv0ui_iOR3otEQI!C|s--4K(zTZFD5eGce)z_sVPM1d|{ zHP;zd!*;drw-E`i0N(V3>`ZD24aMYJdE6bj4_g4c-iC?v+SL|g_?ZgbP`f_u#PHlF zRwD2_vsaFlAWS0{V+kJMm3OOM)nK2%;@DByd8{^(#l{k|<*lrKeq8%4F1#2l$xm7? zFhGZlmv!|?nG#APa#G~~Qr%p8n~JV||1(-KpUdPCJhHt1lm$A%=D&@ei1*$otpb4Z zS65(5UZ~G`O7QDOxVQM~9$x_l7&FKatvYL01~oXJlT^=teq(a`?f^<=$=vm|z;26& z_22as!tE!xv(A8(x7u`n!b;)S6p zHlZ#di|^m5p3&7HM%$&LH5YNzyewN;&lXXa-YEjrzBm0l()QK#Gqu$p@@ahQwCzl=7xZPX;a&)x~x+N7==-JNobO5rJI!Kfy+ zq)!;vpm83+s`xSPlpR6H_vlIZ&tQi^js&lzr)87j_F+HOt^5u!HuEO?HU>o-)>jN? z+{~+r<(W}SB=gehk^<}ZUE3FR{))``QZjUnSxQcMHeL^KNJLQ+cb#_G95sy(YzR9H zcVC`q%aLdHn9A}}hQCHY#tr(ab!G@nEMuTn`Q%;He43MRnk)hbZm0+g*o&nfO7n!8|HmM;`*#q~de2@;%U>OsJE*DBZQeklQ zm|#;Ku*uzkefXlDq~5y2xVj%mj1_$D-zp^NKzK#!u6>*nMT#2~`j-=r6>2Y`Hg>e? zKB^z>%^}|n6^}Tb(arTWjn>^Byb$i)K}==r2rK6@5-hO&u%lu*STMq_4p`71=%I`@ zq!+pRSHIRM`B!#1AA1}AQF2(J743h`hwe#V?lOW zetv#@0Mq3&2#9|&0wj^R{quvEHHWKZbO=E2I1j3I32!7+Loq)&SR z*|WSIDeL?SZ~wx+O5)1|5qE~cJ1%rMN7#UDvN$6G+`cVFAce$-`?y{0ziB@Ed(>Z+ z4Oz7d#yJ_5OF2X(J5D>ILafE@enoMuv?}iiPqOcY+iqQ5HHNzLe_{Fbw;p>4y-00P zP34t%n|}5%YrcmkteLLJnmAT)BU0*e$J~dw(9Y;l0!n4s9Jb!aq!Wz~^WIboiT5JILC@)vCKX#_omZ{h@^*B1F&os1A60r)ANE7q#TI&KmmNnDGG^D#PV;?az>Z z80L!itkY4txWh(srr`?+*W>kO)yRj%l=kWA>GdDZkuG}9yDQOvGdO!MKgkXdaqNN& zFCq~57(2Sl;tfOGIVt6SGn+5R5}e24O^0tSw4L8f1bgYSWE((eA5bo{RbVCVUiHC`mN`#ImL|08Sp(d1olIlM#L zW3l)LkY&TzHYGI;0|@gvkO57}P0D2*@PPL6v$j%NK!)Z_bCFoMrBO(HhpstXH*yD6 zP3<}|Gt8i(mG;ky$U>Be;DL53V}%9WyMGTWq!19=xyniP8~F+K)N}hz2@m6^rVA9) z_xK#?IuK7yqhbNq&a&vA$!---*a-tHVCV6jFWzPE@%YJoPT^KOHqjv+APq;qoNVWJ zdITU6>MT$L!a*m%?fN27_y7)1GOoJiH*o@@5ddqG_tRfeJYU}xnWvVwO(g!DUf5_% z0w>UVOQNCS=BS8@-@oBsn)o}W!+l|o1go4sk@x$-BCEO4CdPIaGr1)TVOuZqlUI+Gv7kiZf%GIwI4#?tm9DbuOS;X9j{Q9)ha?Kd z&x+AS9Zu~lag4;j4fW`@sD90kj?oi^UsmdEdE!;!?YVFjRK7oF?12t3%I@<85$^%O zPN3H66~JRrV9W00Q?tZ>(SLLAl1%1!fOX+DV8?#CJM!;(f^UItEdXH2&$~m=#l7z| zN)6k|WRn?tfbi@r%E1TzDKncXO5XFpLaS+)ziruV(RjL8ool~;e{MG*h}rm_0Fdck zGPPa525YS_Y^qNeBmR)F&FlG`$~M5RJKcz2dLnA~UbVSn09J{hNyD;R_e>cn$OEvb zbUn}j_7YEYi_N|mPO5^IKT?N{o6%f5DW*?3wi^pJ<(;A+lpLtfG!tY{o`7d|PEJ!0 zw$lZuK^;DimE9Pd8z(@C7l0{vA3LxbrMpmRW?D0_s;a<^2vG6r`KNrS2t?G3P5`w= zepBEf1h)hJTg4S@fUy1oczYXMQk)hh+y6DINlr2EaYzRV7AD;+QFg4`<85byxquk%SAdX*V_!v{A2_F?&KXe6XLS zsn>y$%F0oEXRr^tC}KBBW{+g=(B1NM^iZ;xeM8~C*I;uH?FrRBW?h1@Nzlg!G4}Ub zOkcF?nE2%RoC=lXAs6PeI3jKkKz^?^_!XLl zmiiXc_=2(UY}x@s>u<5g;&BDCF5Aw;zNsmT(zAOTCYrKoEqZ7CItX**WvclUKK2(G zI}B*$#lfg7c0YUmWy=WyZOnb7L(b{c@@K12z`nh-UAeq``%sG$>7p2PK2vKTj3&r; z@2VOQ{GAwjy)-xRwZ7a|u@wpC3#Bo6q7$6>UFqsM_9SXGIZ*@X$ZEwIFmSbTzBMd4 zp@Iv;@culx{V8hc3}`rJC(P**yPg6-z~d7ap`m^BwOUAV{Wt@9`OSSEy&Y(q96J&5 z*Yb8Vdt|5Bo`OsP1Z>ZNAG^N1xIS|5 zj{q!fJcz&War@(6fgi3j?q}DVF`M0(g6uiCK*vaE;e+}_l&^e%DPURx;uuxB+*(6SOnlZSxDts=ud&g5 zP|miW|1Y=m)~+J}P75R_HmeHK$MHw5RkFA{+V;r+=dxy*#H)T{^#ecb zAMPtn!z>HJ6OPT?0+8P|8YWWGu@s$}q(FG{naKkO*FYwwG8|P6Pjj@BB7me*OrR`KthMn7UKky02Va!5GFdBFby=+7IJ0{|b&IPwpQn$~oXcMH z*|md|T^VgS3M83}T&tPayO*waiNC4b_pFpSD}?t)ZkBw_uYP=9>MRgRKF9L9 z-Tge<*7*%%xPt$76_g|9Ck84+h%vcl$^CYJ^+}o(eN>-jq-=D6m34)~Bmk7KAfvJY zT2wTU8?&A&XYk2HO*TUF-g=09>kfM{THJIIwP2#-hmTI(Q-krt{(hH6=Ddg~sai_Ws#ybX zlm1vr8oY>|do0jSzK8sWM}Q~KPV;f%0z0G~!t8VX2^=+3t!5Y6r!4pc8EPQ0L;~lS z9E`A_zI`3fku)Aoslr>(y68;JT<;FujFTkq3*XLctQWMNk^+Ru-{TlJ-#F^w;E9+4 zC3Ch4om_wfbO`h_b$I%H3cbuWGYo(=Yvpn95Oa1WNuI4YXMN5nXiLn_{=qoy^KX6? zECe4=r}}{E=Pl!k!An*FK9cZU)yyi98@8+5C16@5M8q93rD)D=Ik#~6Cl zu)3mZZZVQ-R!+}3H!INL;6vH+d~MJh@8L;BG-CeLkD4shF7sVYfe#bteUGmHR`^E& zjZl>L;m1`Tph0eey%*{Ojl0<~Eueh{P#}W<_4IUa6hyC+>L7y#gS${o;*Sh1-j|!B zEnlKWh?SC6NKLm{emQ_Ykn46Y=j^G)0m2r1u+!n&@pazOdm4Oa-5^xcOMBxv`d63r z^&j+{NB+Sf)-ZufU)C0yI)lCtRTS z76#Y>$n;y_sHR6NUKdLHZ}IjKIi0l3_uIR-!tIj!t&WBB3o9?Bv9)P%p5>J|!qXa& zO0qaVV_Mi#6EkFQq$RtP{W1t(N>=R`rbr3O|I?YN(FN-#|OW zPn+R$!SXvMVC9)ChX3iRx{!7!d7C-OkL8_kF=rCDZGQ-4#NB^+TwtAA*{OOf{T@k8 zW#bl$x=CUf! zACJsvu)dm3&O%@M`}c45%DnSV1z?Lcnw;GBpFCJo4&L7!amuEC9uQ?O)t!y3F3zs6 zM*{cwI&{EZsj9ub398x5tK@TRNHg3%%KM;67j++J@!I{gNJzf5hXC7EzHeM1;Sp%e znlTWUtFv9wnCFI%>KDR^h@8&`W+$Q$<%kCTPR3iDi7{D1Vxi{&9M9g{V>QpX7Y9%@4_ZvY zat2ptpeVvv%1md3@t)hYOg|#ExB#*0NujfrpsLUb3W0lQl>l8aw}Qr^-yJFDivK@f zPzN%e5*7vSPvOgnybZ9S(-{^sx3EAc&75Q3fL5y)$BNH|JPb;aXYLl_ISZ_tLbM0%~5z)PR`Sh2!4r{Y=;6Ep_pd`Hc+ zB7rnAUfq*WCL9NQXY=+I{BTBP9LjCFK?sjbP&7l8H`XhieI8`?y{5>t(NgY6X9H~~ z8syipMYW`6<(3tW$~Yro7~W!XHZ~=jekWlN!e*V2SeZ}-i*2L^3vKrCIJ#e0CeKH8 z-ub>1w<2h^LJka5`>;z#p4+bSBG&2tp)rcO-(ZN-kQ z9a_feDLRBZI_J?h%A1sR6zrQJ1c;EuNW!O@B#V7(c$B?CT)Q9oQ?wso&&AuMgOujs zH}!NKPj$gYC3WT-;rwdh7B20cgrPrN+(Bd1z%+?PO1tpxPF2KPJ{&~{2 z$RQJKKglG*u5jCOyW!!0M1II?P03I`hXJiMZq+k$aGvB<0%G88#?LdXtY z%hSUniqTWw~fA7{c`RFVHVXOEg%-Yd-I%M(}al9!sM|n3vT?+sLE6XdALsRhic2J~N#wd_{jdD-H1;0j_+q_&XAy=*peNv%ocE+;+$4Rl z;4riQ6^7k-RIlqeA)TnBI-R7?TK?qda$QXM-vJQw{5HBuDm=NDv{@#p}Tln^-XcR}+&8#T-ytV=NQ~Vpk=oW1YUj$C>!_h%u zydz);P+vcD)kYBpy>|EhfFggnR5^w0k_n~hC0fq{th2-)Arv5lfo_1mdRlNhIJj|( z4`UMxmCxy=c#_9HhRJZD#*}Pbq=sX5>bQVozm*<#xcy*;64I+=ZxSqotnzU-eg&mq zeXO2HU)vV!msK7Ue>4yku>>TPfegKEp%8ZMxqwPP$*b~G}+IZoE#X*Yto)pp1JAG185 zKih0t{EC_zV)fssI{n|Mx&in!8|M$9PyarH{Gk@b<@SydZv6M~0U{)+nDJ-c8?EbJ z9c7R@OEwu~{QXXe@9RnovC#BZ>};_k{zLM>KCRI*DlvsaImBvB6!QIf04!el;I%3n z4%G7UZnIAS!DvdIQs#S_!Rj-8RmbUdZcPVQUNYT<7myV0A`$#-&)e%#Y7ybV?e}Db zX^CDm{U%wmkIT#`bQ-tR&yt*m%a>jp#4YJs*?QQBH&R9ic74&j{P7bwl(4Y|xAp#6 zjH$OK)cc)Ymwey=6L|O--?njPB2@H}YEC!UUXU!1p^9m`DcX&gm@oicd zbh>k;L7x+~U5tbbB`oTil+ZG`qrK2FRkBllUhtcTlvgo5#BtHq#(;P+rKh>>^BD-flUm|%oZ5J zT}nc{v!U6vmDTXq;RTv`TcR*a8k1>bY_ZLg-1=zJm$>}$w~w|F**?o|S(B7SYDg}7 zu@i}{-F`7N;Td(dZ444YmA!TK(5l}=&QKd+4<4e3T8Dm*mOR(Zj2I4GT%lQ(=Fg;6 zcr!-i`01D4@j_S#L;k|Xne_6LSxM>?M*;@TE9$X+UFTWt`R?*h8&vD3D|{~*fS$nF zMCp=&h$yYm;B?E3aVI>WxcAn{LYxM0Bj!~9W$hK0qD*G;KfGMd)cD|^@v-dWRqYsH%B&ft5MXwO7QpF9K{^MrhTTdgRL9k31k?x{9Wq>X>;_E_T8G)Vb|+(61uT*{d&P98kLqUKI>Y(v9d+@^ZbVffQ@mN zdTso#LZrVzMIm2Smf2b&Y57d&P6FUrTt^>BNB^Xn=h~o#F=LBNLzI!dXmO7A>pdbJ zw!HkK(auGeOq~>4shB74|Kej&SLMS&Q&4)6nVK0Yv0O;TgnS z!J6&lr(>_w#R=3hcl2h!FX$Qls+KYeZt6TX3c@eqcBNAdhLP{>7{IH;`;Sa+PlnzF zF4GNYoC8uz9ivX~z^$VL4TPk`2?ww|kAJg1u±fr}hyeNLB&%~INzbC$&tG;9&J z%U9&e;`RxB%7j=iWx3sYpkcy3b7eVBiGZKZv|#lU=o*fz<{+)QKIldi3RB))Xgh+>tiuaO_k)l>yrHqZwf9iYMXwpsWr&}`F z>O_O{H?75a`7mG0BH3cFiip@BI?Ojrq{_2-|FJ=`084;F`gcFZcs&jROH?$}(JiNS zzg0*wAM5%or%da3Kc@@L`yniAdwrcelyZ5*yId23$`^b^gcWw*_*C*QrYz`w$1$C&93r%=~~VnRe)1xzA>%p2_WN~N`DT^Lpd$9e&Jaf2Zku^!)KZozthd_) zN7QZgt;7KVEa-%*?3$3uCSIM-4tP<1Kyex7rY)+va$}pFtMgc9W2Pu7qepf}hL*%J z2Z4h?q{W2Q6#im_L6OB@Wmz8&`mVHwbyf3_*TY<^QHp?phohi)@U zv>E@WLXLX4&xlGgGN3J1#RioPYr* zb;FZeqkUXU9HO#pg^-Z@V_OjtZTUH?mxpt&ImD?5LB6D_wulX<*q9>x^%uEIb}qNe zB=`2TAGgzKMx9mlQ!oF<#hI?YdKEY=$7(5k(9_{;h)x9;0dwK8%je6jkQ82GHQN<; zINEP1>TTKEE4jjCIt)?Z39Qs9TD`xv{@2Px72ATd>4!r7(Bb1zK-WjOsf6C+?hWq6>cKpP53~4dC?a61nrAZ zhKYnI_+xKzd(-8YiT7waD$q9pE}-U#KPQh?&%V{z^U$pD63)~p%y?0Eg7DNs2w&4Z z+mFgut8ia%)v=KHo|Ef1UT)BJquN&(TU#QMoP`xkI#s<^J!H2|?OlFHKlm9@qjofC zHFZmACB1DD5_*M63a19gPp%O2@#pqQfJz{;FwVcWSbaz zF+N19mDY`fdo%$--ms_Pgzdfb#Ok=s-~DMioTjjS2^jCzh39HN*=s=#^loq#$Zn1q z^Yx<2SS&S~BOn{#kvQSUML>GH~n&-Upp%T-HdN!YMb!OQ-b&uHfw|H)%%=oHiH9x81VPjsql7=|L zl@RbSLB){7)-`u$>`VW|^;AV2TNKQQ5X#@|Wr{L1+pzJl%ZKVAy@_UIr7h_4Tl^ev zMs|3iLsJw26(rZ&Cslhv$nhsIAUbkwe^LH1u$(qxJ%noT*9PJ0o4*8oE~6GaA4!|b zd+fbw>CGhD4|s4aQ;qk)o@p&kSJ?}gD|(4bw9XKm#RA7gla8E-^ozh;hWvGn1hP*t zfPd+Ko?rD+l@Nk*Mb=~UAs)`6*n(4~WxnFt$xxB(O z@`_B&2VX5%_~8{7H|J}iqAoag(<72)=amgvC1fo^ZR}n!vlmOC>pyO+8RKx~`YxoH zXhJJrt6IX#)vakpE=5~=5WyeIC@uH?tURjT`>-pDkChCKPa>C}gvi6cJ}quY_~OS( zo{(fsU(7M0Gb#SRh0dq^H#yk{A##X=uDE#oh|R)f-kJ*tZ~^@Zig;miSHZ96hY0lE zEBbWiF+MNB#Cvalr=fSJop2>W^3aE9z_XFdFQfYr{gPAQ{B`6#VTju1$rg!d!))eg zj3KhvKye)ScKn`a*W}dUEBL%<3o9()#A$h2{Z6#=iv&*%7>QG1#oR6EB}C2<&>q4etW$;yZpcdG5Wh?%6GC$$`p}Vy7xIEQx6{K zZ@aTa8z;B1o}?mC&EcX+p4bI^*kj|YpsYIIXPiijRnW3l>tNFk>}e*VHoDwuEnR%Z z-WF})Ra~`#klh$X#JQF4G7C7Xq*HZSfT5L;TT}?@Wm(ws!Fg%6eFct%N*kt8ErH+G z00%Zb^a0_eBP*t{p5$MJc+5=5{#*K54u&q0OvLej8BN2@F_1|ySrNuc*1O5ukB~{@ z&)d?(wZFif5K-QagQBy@E2XwRyDZB`d}I4G>u^9=R#<}F1jc~Ru{2EUEb zNL%(>in9sMRg*$VRdIZ(s{$5&tfMhQHEMBq8z%94QeK+1&A1Unf5`dslxPB;IX0r8 zDZ(2cReU!hGcoy!10{%l-&krm%VAB@Gm0=d2mkiRpot~zi6+#P#k=!(Pj5JL1DznE zx1N1+dwJ+|K$xLFnRQU?K%AN{tkE+jzOQ}#Ngw}hOsN^nY*kao>$NnX@Y8!iMe4J+ zzLcM-qN{z>TZ^`J`i9%ygnFk(a9JNm`}-;o9|grEah~Rrw_@46H5Wr{SSiJkxNF;* z_LY7(M0ibC)t6F|&lGg;B^qTRP@4W-YBXY=IY-O`pM_bkoJWkV%q8R5zyIxM*)*R? zZBcb5V@iLK3Dm+ai2n);JT3XwdA{I!ae0Cs^x2Sn8cG?h#>l?vxqpbBF<}(NW~&eR z*&wgB?K{O?>22|Rqon*!qTy6x3*lUZjy|Ks-DapHN!DR!JLx}7uWBBfH6_IpXuuO#;~|El13c2gp;<~aLO5j{ZF}{RtJO3( z;lPhFPlft@R&t*uj)OVps3j3hDv1kV3@&M$0x`R871;&HVQJ63rkt%)KFG!FlC?F7 z^Ur6(K`tY^581vn*f@CW-3{{PK-fOxm-$KX-DXgjzHOgG!Ueh-)PVfnz6P$D*DU$S zcuj@&#R_CvANhnH**LY=)7kE>3^=y<ic0CKbMq0)DkC^DR~Iuvh}85QrEwt`#-X4>(^gd@e~>e>a5tOG@iZQ9;;Oe8!#Sh z_NZwAz5t?NGQKyOy zAUX2K#Lz-d;t!AynI0bbRz&qh^q2RoCJc7N1=+h|qb}U`9veA-=RAe@zpZe6{7#W* zgV0g_b7LksnQxqR|LxuL$ckK&`{m&Q^M3{$B9S-0fi0;l4)XrRfoh1(>j-xF8nJ7v zAXpgDL1?xP0}KDzK%>9&mgV#JZ@dm8qi{FV(%*IzhKgNR<%@r&D4Fw2sXreOjpLEVTOk%QoLaqu1mun!TKx zP+>#?!;&<|LQ86?2}M}7MDd?qOr*$!X32`2&XC$$=cz?lH21=~d#_WX|P@Tk(`%@8X}aNs@-)ao~pA*|lo#(Hd~85hf?RJVbQ zN-|pLL)u^Y_p0&^f3m{4+j(AvJd@^}^o)~&&z+Y@T!9c)l<()`vsGBwJJ2!=uzjz! zf7HbhFHE&w5w*saR#{M+L9cL?_=HyKcs1oltcoLwC5=doV{8Tl{H50P#`i*glFw$Y z10%oh5XT*Hn0TZwle3wJz}0n91L*`A|oD( z$Z#jci4=m{V*71@BOg#+5Q_(gj6@1+{Z*1#81j6g&4^RLO9*Q+qs|dZfzhqtaG3ob z5s^yukcI_VsYaSWl-WB~?s{y>kfP+)g@0VIbriH=ISNrY5xNMXgj0;mH<#zA(HVD# zm_vB7;5j&cBVWnA4N{})ae0)9?Z5;{JwLWr|9SuUDgb>W9O&`11VsqU6bUXzTrWW0 zP?nJp#u=R{G-&bjy~N*cCzne2lG&id(br8+ER1Vwh@faE*tU=(&rW${GnaKE=KtyL ztb0g{4tWJ!8oE-t2Y8nS+-LFcQcoNQ+cb6ULq70uB+Tj*?GRDdz6i8#VjN<REB;V!o-)+n^hJKJxgKsWR!`N8A zr94n^&P!}F_^QTY9irxJ%P_X74tA^bRznW__7_78fr`MXW zktwsZ=uGV>w`JP(J}m4>>H4Fxy%i**Esn0o@#p{1IHs>w)lC0`@*U06^^7*_Z&{lL zo>g_2g^RvqGSt-d`JI7R>-E`0VLN=nPux)0@_JyF@LzLtug+Sl9O zJz9{DiWM9#oVRcyj^~x5*u@6H5`W~Mk7)37ne-0Q)NRH{7~0Mt`ose;W?x@Y)<~k% za40aA0GWl%lbLUU$91RX3u*L3gmRPwY1%31gKM_K#;aMgw=+>ttxURc4l5QGu_u+|uol-FR)aQHp8@6zH5=y{|%kIq|?Av=c z63Y*0=!g)Mg&0W>eXGjL!s4R1Wd(CFd1LAaL(V=5l}`x9)xpJclDGCMM7qE1_Sl8V zW%r%O7xc;w+17Q`JX#_LR!PpC9E`LwgjJOhxg-t}vC@9xZA9R}IIrL5C|`+9w+^&8 zS@R6hyV}4#aR2AY3Z0vv#SVIo`B9<}0A(=PQ$)tovO@3sl+syX1qka||hLzI4wOsY4_GqQJ92FZWQHKJntS-QUd-Z7&Pj@Zi%ICXX4 z;%%ubS*kQdhodjdlPE!#UaUt?FIB+0!(oV!$I<_wzR^ujZr8gt5TfjG1vh_vcbb(p zVA}goqKme@HkJ{%oavy#TuS*|e-fmBuzp8Fz*ID^E8VDhFHVuzyW%>(b}`6J7KlqM zOR%z+zSE=GR!k$C{gIcVTnag!=pk>zh*M@ppbLFRDH>= z?zS?EddF^*NDp~*B4N|F))vv#Vr&$FItA_I3b1pzrETXXPr z@~E(?`shg~@CF`o*Z2>oQhyQ>B0Q6JG}c@=%_pwfWjj`sz#($nv50~9Q^qt1&#i37 zRx`~gH>e~kUWzKgASW;Cjb!+5bzh~0M&d@Y+XyxdhkT_YPVX(he0|8sk{KgkJwt*R z?Kog==cj3U4s?~hJ4S;vTCT)*JYdx(%Rs8z{yP|Om>sGTT7T8=RhB(Rh^-u^ckE@^ z!x@GM0vP+9sELpZPnR4b?S1LOpQ)>Z4}HSwO0)Ul20*-#g|tx6(M__0f!=(T_=Ebj zrJ)1N)$89K66R>i1Tssagw>hIe<>qbwH}nee4~T&d2-;8--DjNn&Ia}A8G@C@RS^- z^LDR=6aTjXaDz#D(E#$<#KbdkWk6#(?dX+kXIF!*1`*&C-;~!$mVbp3(6I*B_WZN} zQd>!{?#d9VE-=L>rg>;C!o;f43+U>79;~A6~pOTv-B^C0)vs zcZR~(csQjGG>ox-%YY`-btO3IP+Nmt&PxW{~v4#t2&|afZ4zO*W@u}Vl0_l zKPuWEQ5r6U2x6sI#O_JkWtHil9H@!X97bOr5zj|Z*ZJI=Bms)Nc z?4|hl9IW=G4tj**>U^UfOS=J0P&jZ!%tT2Vi+~{y-J$g}4U?ftam0NLdtOsk>qRUDQ`w<|$Qfa9IqhOxOfN|Ja zDrZ-+a{aK4-?$Kr<(778fZw1OiwL$WIi09PqxNS5*n0q@F&-R-pyxXLH&zA!NiVC%C&yf(8lhu0a}iLU0HWB)EHU4Fq?ByEhgbg1a>q+&%ay-f!G- z&p7;Re$Y)-*WPQdHRm(W8yXr2Bq`W*F0G1t#%$R<4wPcWbzSGQ(mt%oA?+n2P5`!v z(^a(<>D3R%XZk^0I-MnbjXKePg=Dy3n|uTZQ*{=mz*3g7f(0@NB?fR8EhIY+F~vTe z$V^JppM?dJN7Ag_a}Jd?%apUV%_^3Osq+b_3z$PIKj->mb4%yQhJbn3!pVso^v8VoV=>e`}o_-%rFH;uTwkGxesAW7*HJzx?v8@l$u z8)a%m#=g^zXfvsQb1%{K`gMVS^Z4Jg6}l!{F#j}Ls)U6%E0o$>F|)bx#s{m(Z+%ge z>DUzWka{rx>|4-U*!$D4uYTx)JDM&7RiIhc;mzlrNhz_5QyZt09<2U+4fE(e=A}@2 zKwx`dUXJ>~K|DX;5%1tHg^_7cI&4Qt$#*Z7o9mqhWo|=RoBAeL7H%6FuGO|EeFNKI zo{4{>f@Pu?r+^d=pzNA)&bi9mZtwimzz~8++KZtz%say=? zDwt)!;1j16Ns}ZN`iZZ{SRs@?h#!v`jKHNCY@{fT@JDT(IeReU*U3HRACqizq?N1Mimzw z!nE}rrHZgL)ERMaWzja?$P52<1EFk%NbV%O=R*4Dl{M*7qgD+%N`f^i0SEpbbATeF z8U+2)e3a0E?w z!`r91?*_X@$IUc!fok(FA2<4bvb0a z*Da(mQ3B!N zNvTp^xFOixV`q0A8{XCT`F3FXd|`?eGywSEh<(gu&HKE6S9bcFpRZtD&|k!FQF$O$Ym_P;L9Jdd%*xkPtB2m8uy%&j3Y?M<;UHc?ybATuQ7M)@GF)^nj<>m-1~F)1PtlD|3o&`W zFL>xw3$Ko|se12Y-|1wKo+$0FUx zmB5o>-Ny!*0uO86f7-DqDoQJpq^DUiWA}?tgbf71KF%;n38X0?a6j$%GsFxI_+)rk zKG>Qeu+a?>uL92{KZ8d24@tkIYmHTgN-{(=P{ES@z+>x><%$>FI0e>#2)5Rbc06t& z6ggejqR$tS6Z|RWrzA4lQNcj@pKE$Tua&Rx-%%hpqvS+nTw;dbe8^g8li1hzUO;;=U5i+rm#{j2;({57~fB(sOXxw~hC0dm_a8fB6mYOF-u_BhxF`>MFi zhahuea<$()=MR@KL|}1FMhURQ)#<%We4|8LRqW~q8|y67(7}Jy&(ZX}n%5kjuf~+j zx%*!E-;nCcDWRd;Nv_ERB~5vU9P(y4;Kb!rs~qoaS;_E|nv~lvPQgbeb_k>KBDFvL zm?ih^;+@BC^g%>!Y&MRYxnaz!&KzmBuoE-9@4MGl7SYMP*P(!Q+_K?SF%IbXT9%4)QyOWQ z{Y&(2CD{#ko072!~QU5H!Gm^iK#XS;XhzkA(9G!qBg>=c}FfxCtGP zY`G|Pu#3v|IIesNfV*Sp=uqUeKi-R@6>5y);)NDvR_V2^pd22q`&%(EYeW-dWTw~&Hb zRfoKbB~jIrC`BjvCGS$#whXCb4h^>eP?y|E)5`rFyE0@-UY8erZBgrI*6uXE9n@28 zja=cp8tR92k(Qa~3A`hBo(9vF-9u5G-~Y0ZQpAGS%{CchZvS4*?WsfMY$U-OT)V^H z_^{^CL#$3R>etr(;py)XRF({}-G&6gBVT6^uMT_|dEF32|C%yo{6C$~9y(ir%`J8> z?OUtNB%2BRdQ8!+$VnA!YeQ#-Y8a4tUoxh>OehqwDyq5g@`BSou5o{mLW?PaFT&`e zOGK+vuIO`sXs=*p5)*joM!jDXbvCNJri$@lZAj?Cmp$ZMW_;VJXn+&7luB%Lhw21y zO8Fbj`%HUlk_u7$E?-;0vWUqarnIdLWcqx>xBR2S90+`Q9o zD~=aG?nQVUgHqM9SK-5Y!rcqUy9(=rb3M!V=4O@0@uX1>`BT>Iy~+oq>X=n6*|$s>^b-XS`MssT+-^tYDS zJQp4cbXQvtc0_fFZ~9VZ$!n5)&Yc)vJTv2WU|5rJL3mB`Ts&hP^%Q+e`-h;DC;Q{z zXs098f9zagUJa5W(--nG^(?C*#bx`#6#Ec)|0g5T|;% z>nx0ckIum^sSz~_2qqL(Ob1Xr$JQ+(_V(Rmua=}nVImIW2P9x*Zr=>J9xsO6iW`e zJay3#y7L}MwRh5IgquIO>soXJKySz?%21AyO_6Va*vJ-z* zS0P-Ck!d-&OH6ayn=J~{lZ~kJhHrK=i*}Iu<@;76oxjH@MJ^stzGuYC826Y}G@(9f zTd;k<37ywN39g!}>1z+Qa+k(hOfJ!%LCjg^e*s7l?@|BtQyRm9yj@ns1NP4n8QWk0 z&5Z+#J<(!ZgZnqY@gL-SDPYw6#Q`ub;Pk|sm`TokERpKmT^U;j?-U!pPaj{JPD7M7 zs3oCmhnkLdPPWKtalgc0ysmhg>z_Z{V!EC0MF z#CDvMRm24k4o(`O9lq;j3*&uy2SeP%5TWW-DSZb*J78SWhjNHRt+_p#Mow&L{QVQq5s($F~GcsG+6OAb7Iaa^Z3hg1IIY8SrkVINuq1mPb)&20J z)OSQQ6bk8P%0Q8Uyf%{4O{iWhcuqCxn>y=m;;$T9)IANIU?aVFvH}O)AMt9;`^mL; zdIfdLJ()#GWT*UUs+mRO&$m=A0IWq^bJ&9ProKoVTT82bY}I$N0b&QMmG9PK5Lz^d zH8pAbQ^Pj;X3hvPm{O9=w5q>@afGwAO(D*9HNR8NPdqR5w||)IjYG&Y+MIWv<1^Z_ z3_wu>HVxo_L5}^|1p7A4)%nLkSzY5=8tsE6tg@7x+bcZgIqY>Bg0(sLFs8OSSnJ2( zBET79ANxpr^%Y;%UfV1CMzWHSygp1fsfXw=XDeBttP1b_RVSs#_}QmMFA0y3OQiEb z6jF>+Lg96K2Z)enYg4m~^`b}rCAEaF`UsO)ALvRe{3dTU6O=(=!nJ~?j+?Hw^3**0 z=WF6_2y^rU2HCP8q7*@=wcb!_wgEIGD#`g8)}`zkKqMpP_%D5A1r>je8mk&JE&r>j zV5at1+j9|7K|sK^ymHyJj?DzI4Qhx5JyJIL3hiXVWx44y&4=TGz+XO8)j)s3|4$%Y zw^kiD)t!I4`BxfCLy=>z_$qqJ&4I>hCYp3p;rr*8QRaxhZo4NTcBD)XUmV;QYx6$J z8w>eV**3+ z#Ko{d{x3)CB;%Ky!j|upU9(Ujc|@yA*wT|LvA-*gd4DF-z2OX$MZZRLS3?O4@xkbs zk-CB#_uv}59ezg|U6;zmx&oa8;(C#yh}7?cHFw8s?$1$<@NrQ}u5$SY>{NA)E!zIih==a<8NUrJ2m7EM2% zU4h;n6K6jB%t4UcZ%||1CT8wilfU(b7tZWc-#z)#k`m=1S%VwM0{7{1RKx*N%QC_V zpdikg+;}&S+6VZ|I6kI|g4A_;f^rzac{!EuA$ciocJAj{mMU-RI$&OAf>@;F=?Qsb zGBv>W7ab`c+cj*^O9++M<{$ZQn&A5zggp*v@wjuOXkW*kLbs0ZyOkofyT*k^KjGS@ z5x%q2pDMSW0^xVjcj$Aw(=;c3o#_4>1ym6-dUfEu^Ji21@^>M}y5O(isG^$>uSpJo zHZo==z@+m{;7Ly}@dj|F0R0|>wh+-nHC!}DdU05 z7ixs}wBPkvfCZ(BbbJ=tP?{xf!k(!R1`^^PnyolBSg$Nv=@X`^md+xcZ@HU`^NEvz zXO-12G$fAP|Kw(eA}pOXtE4o@-`63|4=EXlK0d$MHB9BHbj#E0u-Z*~Re^=$z6wUD zG!BEfFeIf``lFyMc$Dvy{Qb-%@v45?^eFWQo(Yr(ydTWXn&{C#OD|%Ln73-ci?xsN z4$Pr+F{)MQr$Sy0$U=>@@xiK7d)#MqRlus^n64hzf14=IlpIJWpKNe|W~vRmdYc%a zb>hzwHe&sv33Iyiby@b&owRsejqPLo?I<-oRC-SM7ne;h z6=B3ashIvH;O;J)K@Mc6)CZVU$+lv8Zl0K}lRMkorPd@+T~^a3j_$e1Qh--@w=7~I zDvr5sCj7l4TqO*i78Wx{IFKvge-rARr z0@;OVV=n6}hNJ*R#g^LX0mAsdu53ifgb8t*jjiqE#Du@YX?)~$&(1$D(ck%Ni`ZiH zPM?9Eu0=1_dD01bf!c;M)@5b%^(sGze_Aku&;?J=nL5laQ`ge~;DMo39Ksz>KS^6F zv!6=RuA&{sudH=%ew%0ykJ(&LuCW=G9Pw^mon5@vJD{3iEweE;>Szarz%56tHaY`exDf!+l?Ynw1({z=>@04Y3;B$kcH7ra-pny