From af2e516955883d46078bb90b1282fcac66a64c55 Mon Sep 17 00:00:00 2001 From: martin-rdz Date: Sat, 28 Mar 2026 17:42:25 +0100 Subject: [PATCH 1/3] implement the cal constants reading from database breaking changes in database format (tables with unique identifier to allow for updates) breaking changes on how the constants are stored inside picassopy (but it's now more consistent between LC and pol_cali) unified function for selecting the value with the lowest std adresses #38 --- ppcpy/calibration/lidarconstant.py | 53 +--- ppcpy/calibration/polarization.py | 34 +-- .../config/json2nc-mapper_quasi_results.json | 4 +- .../json2nc-mapper_quasi_results_V2.json | 4 +- ppcpy/config/json2nc-mapper_vol_depol.json | 6 +- ppcpy/interface/picassoProc.py | 61 +++- ppcpy/io/sql_interaction.py | 260 ++++++++++-------- ppcpy/io/write2nc.py | 2 +- ppcpy/misc/helper.py | 12 + ppcpy/qc/transCor.py | 4 +- ppcpy/retrievals/depolarization.py | 6 +- ppcpy/retrievals/highres.py | 2 +- ppcpy/retrievals/quasi.py | 2 +- 13 files changed, 256 insertions(+), 194 deletions(-) diff --git a/ppcpy/calibration/lidarconstant.py b/ppcpy/calibration/lidarconstant.py index b97e441..db6d903 100644 --- a/ppcpy/calibration/lidarconstant.py +++ b/ppcpy/calibration/lidarconstant.py @@ -4,8 +4,10 @@ contains :py:func:`get_best_LC` """ +from collections import defaultdict import numpy as np from ppcpy.misc.helper import mean_stable +from ppcpy.misc.helper import default_to_regular elastic2raman = {355: 387, 532: 607} @@ -40,9 +42,10 @@ def lc_for_cldFreeGrps(data_cube, retrieval:str) -> list: print('LCMeanWindow', config_dict['LCMeanWindow'], 'LCMeanMinIndx', config_dict['LCMeanMinIndx'], 'LCMeanMaxIndx', config_dict['LCMeanMaxIndx']) - LCs = [{} for i in range(len(data_cube.clFreeGrps))] + LCs = defaultdict(list) for i, cldFree in enumerate(data_cube.clFreeGrps): + cldFreeTime = np.array(data_cube.retrievals_highres['time'])[cldFree] cldFree = cldFree[0], cldFree[1] + 1 profiles = data_cube.retrievals_profile[retrieval][i] @@ -94,7 +97,10 @@ def lc_for_cldFreeGrps(data_cube, retrieval:str) -> list: print(f'Can not find a stable LC value, skipping {channel} {cldFree}') continue - LCs[i][channel] = {'LC': LC_stable, 'LCStd': LC_stable * LCStd} + LCs[channel].append({ + 'LC': LC_stable, 'LCStd': LC_stable * LCStd, + 'time_start': int(cldFreeTime[0]), 'time_end': int(cldFreeTime[1]) + }) if retrieval == 'raman' and int(wv) in elastic2raman.keys(): wv_r = elastic2raman[int(wv)] @@ -124,43 +130,10 @@ def lc_for_cldFreeGrps(data_cube, retrieval:str) -> list: print(f'Can not find a stable LC value, skipping {channel} {cldFree}') continue - LCs[i][f"{wv_r}_{t}_{tel}"] = {'LC': LC_r_stable, 'LCStd': LC_r_stable * LCStd_r} - - return LCs - - -def get_best_LC(LCs:list) -> dict: - """Get lidar constant with the lowest standard deviation. - - Parameters - ---------- - LCs : list - Lidar constant for each channel per cloud free period. - - Returns - ------- - LCused : dict - Lidar constants with lowest standard deviation per channel. - - Notes - ----- - - Since ``LC = LC_stable`` and ``LCStd = LC_stable * LC_Std`` so will any negative LC also have - a negative LCStd, and thus be chosen as the best LC. - - **History** - - - 2026-02-16: Added additional checks to hinder negative LCs to be chosen. - - """ - # list comprehension for nested list - all_channels = set([k for e in LCs for k in e.keys()]) - - LCused = {} - for channel in all_channels: - lcs = np.array([e[channel]['LC'] for e in LCs if channel in e and e[channel]['LC'] >= 0]) - lcsstd = np.array([e[channel]['LCStd'] for e in LCs if channel in e and e[channel]['LCStd'] >= 0]) - - LCused[channel] = lcs[np.argmin(lcsstd)] - return LCused + LCs[f"{wv_r}_{t}_{tel}"].append({ + 'LC': LC_stable, 'LCStd': LC_stable * LCStd, + 'time_start': int(cldFreeTime[0]), 'time_end': int(cldFreeTime[1]) + }) + return default_to_regular(LCs) diff --git a/ppcpy/calibration/polarization.py b/ppcpy/calibration/polarization.py index 25ac814..1cc1d5b 100644 --- a/ppcpy/calibration/polarization.py +++ b/ppcpy/calibration/polarization.py @@ -6,6 +6,7 @@ from ppcpy.misc.helper import uniform_filter from ppcpy.retrievals.collection import calc_snr +from ppcpy.misc.helper import default_to_regular # Helper functions @@ -176,7 +177,8 @@ def calibrateGHK(data_cube): data_cube.polly_config_dict[f'rel_std_dminus_{wv}'], data_cube.polly_config_dict[f'depol_cal_segmentLen_{wv}'], data_cube.polly_config_dict[f'depol_cal_smoothWin_{wv}'], collect_debug=False) - logging.info(f'pol_cali_{wv}', pol_cali[f'{wv}_{tel}']) + print(pol_cali[f'{wv}_{tel}']) + logging.info(f"pol_cali_{wv} {pol_cali[f'{wv}_{tel}']}") else: logging.warning(f'calibrateGHK no {wv} channel') @@ -371,8 +373,10 @@ def depol_cali_ghk(signal_t, bg_t, signal_x, bg_x, time, pol_cali_pang_start_tim pol_cali_eta_std = [float(0.5 * (dp * std_dm + dm * std_dp) / np.sqrt(dp * dm)) for dp, std_dp, dm, std_dm in zip(mean_dplus, std_dplus, mean_dminus, std_dminus)] - results = {'eta': pol_cali_eta, 'eta_std': pol_cali_eta_std, 'time_start': pol_cali_start_time, 'time_end': pol_cali_stop_time, 'status': 1} - results['eta_best'] = pol_cali_eta[np.argmin(pol_cali_eta_std)] + results = [ + {'eta': e[0], 'eta_std': e[1], 'time_start': e[2], 'time_end': e[3], 'status': 1} + for e in zip(pol_cali_eta, pol_cali_eta_std, pol_cali_nang_start_time, pol_cali_stop_time)] + if collect_debug: results['global_attri'] = dict(global_attri) return results @@ -422,17 +426,13 @@ def analyze_segments(dplus, dminus, segment_len, rel_std_dplus, rel_std_dminus): data.polCaliEta532=data.polCali532Attri.polCaliEta(index_min); """ -def default_to_regular(d): - """ """ - if isinstance(d, defaultdict): - d = {k: default_to_regular(v) for k, v in d.items()} - return d - - def calibrateMol(data_cube): """calibrate the polarization with the molecular signal Converted from the matlab code to the best knowledge, but not cross-validated yet + + .. todo:: + had to calculate TR_t, TR_c again when calling depol_cali_mol() """ @@ -469,22 +469,18 @@ def calibrateMol(data_cube): background_t=bg_total, signal_c=sigBGCor_cross[:, slice(*refHInd)], background_c=bg_cross, - TR_t=np.squeeze(data_cube.polly_config_dict['TR'][data_cube.gf(wv, t, tel)]), + TR_t=onemx_onepx(np.squeeze(data_cube.polly_config_dict['H'][data_cube.gf(wv, t, tel)])), TR_t_std=0, - TR_c=np.squeeze(data_cube.polly_config_dict['TR'][data_cube.gf(wv, 'cross', tel)]), + TR_c=onemx_onepx(np.squeeze(data_cube.polly_config_dict['H'][data_cube.gf(wv, 'cross', tel)])), TR_c_std=0, minSNR=10, mdr=config_dict[f'molDepol{wv}'], mdrStd=config_dict[f'molDepolStd{wv}'], ) + ret['time_start'] = int(cldFreeTime[0]) + ret['time_end'] = int(cldFreeTime[1]) if not ret['status'] == 0: - pol_cali[f'{wv}_{tel}']['eta'].append(ret['eta']) - pol_cali[f'{wv}_{tel}']['eta_std'].append(ret['eta_std']) - pol_cali[f'{wv}_{tel}']['fac'].append(ret['fac']) - pol_cali[f'{wv}_{tel}']['fac_std'].append(ret['fac_std']) - pol_cali[f'{wv}_{tel}']['time_start'].append(int(cldFreeTime[0])) - pol_cali[f'{wv}_{tel}']['time_end'].append(int(cldFreeTime[1])) - pol_cali[f'{wv}_{tel}']['status'] = 1 + pol_cali[f'{wv}_{tel}'].append(ret) return default_to_regular(pol_cali) diff --git a/ppcpy/config/json2nc-mapper_quasi_results.json b/ppcpy/config/json2nc-mapper_quasi_results.json index e004a44..283f7c1 100644 --- a/ppcpy/config/json2nc-mapper_quasi_results.json +++ b/ppcpy/config/json2nc-mapper_quasi_results.json @@ -214,7 +214,7 @@ "standard_name": "quasi_voldepol_532", "retrieving_info": { "Fixed lidar ratio": {"value": "__polly_config_dict[LR532]", "unit": "[Sr]"}, - "eta": {"value": "__pol_cali[532_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[532_FR]", "unit": ""} }, "comment": "This parameter is retrieved by the method demonstrated in (Holger, ATM, 2017). The retrieved results are dependent on the lidar constants and the AOD below the current bin." } @@ -233,7 +233,7 @@ "standard_name": "quasi_pardepol_532", "retrieving_info": { "Fixed lidar ratio": {"value": "__polly_config_dict[LR532]", "unit": "[Sr]"}, - "eta": {"value": "__pol_cali[532_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[532_FR]", "unit": ""} }, "comment": "This parameter is retrieved by the method demonstrated in (Holger, ATM, 2017). The retrieved results are dependent on the lidar constants and the AOD below the current bin." } diff --git a/ppcpy/config/json2nc-mapper_quasi_results_V2.json b/ppcpy/config/json2nc-mapper_quasi_results_V2.json index e004a44..283f7c1 100644 --- a/ppcpy/config/json2nc-mapper_quasi_results_V2.json +++ b/ppcpy/config/json2nc-mapper_quasi_results_V2.json @@ -214,7 +214,7 @@ "standard_name": "quasi_voldepol_532", "retrieving_info": { "Fixed lidar ratio": {"value": "__polly_config_dict[LR532]", "unit": "[Sr]"}, - "eta": {"value": "__pol_cali[532_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[532_FR]", "unit": ""} }, "comment": "This parameter is retrieved by the method demonstrated in (Holger, ATM, 2017). The retrieved results are dependent on the lidar constants and the AOD below the current bin." } @@ -233,7 +233,7 @@ "standard_name": "quasi_pardepol_532", "retrieving_info": { "Fixed lidar ratio": {"value": "__polly_config_dict[LR532]", "unit": "[Sr]"}, - "eta": {"value": "__pol_cali[532_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[532_FR]", "unit": ""} }, "comment": "This parameter is retrieved by the method demonstrated in (Holger, ATM, 2017). The retrieved results are dependent on the lidar constants and the AOD below the current bin." } diff --git a/ppcpy/config/json2nc-mapper_vol_depol.json b/ppcpy/config/json2nc-mapper_vol_depol.json index 532b00e..daed019 100644 --- a/ppcpy/config/json2nc-mapper_vol_depol.json +++ b/ppcpy/config/json2nc-mapper_vol_depol.json @@ -105,7 +105,7 @@ "plot_range": "0.0, 0.3", "plot_scale": "linear", "retrieving_info": { - "eta": {"value": "__pol_cali[355_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[355_FR]", "unit": ""} }, "comment": "The depolarization ratio was calibrated with \\Delta 90\\circ method." } @@ -126,7 +126,7 @@ "plot_range": "0.0, 0.3", "plot_scale": "linear", "retrieving_info": { - "eta": {"value": "__pol_cali[532_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[532_FR]", "unit": ""} }, "comment": "The depolarization ratio was calibrated with \\Delta 90\\circ method." } @@ -147,7 +147,7 @@ "plot_range": "0.0, 0.3", "plot_scale": "linear", "retrieving_info": { - "eta": {"value": "__pol_cali[1064_FR][eta_best]", "unit": ""} + "eta": {"value": "__etaused[1064_FR]", "unit": ""} }, "comment": "The depolarization ratio was calibrated with \\Delta 90\\circ method." } diff --git a/ppcpy/interface/picassoProc.py b/ppcpy/interface/picassoProc.py index 0fa54c2..bfa1e36 100644 --- a/ppcpy/interface/picassoProc.py +++ b/ppcpy/interface/picassoProc.py @@ -12,6 +12,7 @@ import ppcpy.qc.qualityMask as qualityMask +import ppcpy.calibration.select as select import ppcpy.calibration.polarization as polarization import ppcpy.cloudmask.cloudscreen as cloudscreen import ppcpy.cloudmask.profilesegment as profilesegment @@ -67,6 +68,9 @@ def __init__(self, rawdata_dict, polly_config_dict, picasso_config_dict): self.retrievals_profile = {} self.retrievals_profile['avail_optical_profiles'] = [] + self.pol_cali = {} + self.LC = {} + def mdate_filename(self): """get the date from filename in YYYYMMDD""" @@ -349,7 +353,8 @@ def polarizationCaliD90(self): """ polarization.loadGHK(self) - self.pol_cali = polarization.calibrateGHK(self) + self.pol_cali['D90'] = polarization.calibrateGHK(self) + self.etaused = select.single_best(self.pol_cali['D90'], 'eta', 'eta_std') def cloudScreen(self): @@ -459,7 +464,9 @@ def polarizationCaliMol(self): logging.warning(f'not checked against the matlab code') if self.polly_config_dict['flagMolDepolCali']: - self.pol_cali_mol = polarization.calibrateMol(self) + self.pol_cali['mol'] = polarization.calibrateMol(self) + else: + logging.warning("'flagMolDepolCali' set to False") def transCor(self): @@ -481,7 +488,7 @@ def transCor(self): if self.polly_config_dict['flagTransCor']: logging.warning('transmission correction') self.retrievals_highres['sigTCor'], self.retrievals_highres['BGTCor'] = \ - transCor.transCorGHK_cube(self) + transCor.transCorGHK_cube(self) else: logging.warning('NO transmission correction') self.retrievals_highres['sigTCor'], self.retrievals_highres['BGTCor'] = \ @@ -614,10 +621,8 @@ def Angstroem(self): def LidarCalibration(self): """calculate the lidar constant - .. TODO:: Add option to read constants from database. .. TODO:: Find out how we prioritise raman, klett, and database retrieved LC... """ - self.LC = {} self.LC['klett'] = lidarconstant.lc_for_cldFreeGrps( self, 'klett') self.LC['raman'] = lidarconstant.lc_for_cldFreeGrps( @@ -625,7 +630,7 @@ def LidarCalibration(self): logging.warning('reading calibration constant from database not working yet') # Prioritise Raman retrieved LCs but use Klett retrieved ones when no Raman retrieval exists. - self.LCused = lidarconstant.get_best_LC(self.LC['klett']) | lidarconstant.get_best_LC(self.LC['raman']) + self.LCused = select.single_best(self.LC['klett'], 'LC', 'LCStd') | select.single_best(self.LC['raman'], 'LC', 'LCStd') def attBsc_volDepol(self): @@ -666,37 +671,71 @@ def quasiV2(self): quasi.quasi_angstrom(self, version='V2') quasi.target_cat(self, version='V2') - def write_2_sql_db(self, db_path:str, parameter:str, method:str|None=None): + def write_2_sql_db(self, parameter:str, db_path:str|None=None, method:str|None=None): """ write LC or eta to sqlite db table Paramters - parameter (str): can be LC (Lidar-calibration-constant) or DC (Depol-calibration-constant) - method (str): 'raman' or 'klett' - db_path (str): location of the sqlite db-file + --------- + parameter : str + can be LC (Lidar-calibration-constant) or DC (Depol-calibration-constant) + method : str + 'raman' or 'klett' + db_path : str + location of the sqlite db-file + + + Notes + ----- + The unique columns are needed that new entries overwrite old ones, + otherwise they are just added to the table with same timestamps. """ + if db_path == None: + db_path = self.polly_config_dict['calibrationDB'] + logging.info(f"read db_path from polly_config_dict {db_path}") + if parameter == 'LC': table_name = 'lidar_calibration_constant' column_names = [ 'cali_start_time', 'cali_stop_time', 'liconst', 'uncertainty_liconst', 'used_for_processing', 'wavelength', 'nc_zip_file', 'polly_type', 'cali_method', 'telescope'] data_types = ['text', 'text', 'real', 'real', 'integer', 'text', 'text', 'text', 'text', 'text'] + unique=', UNIQUE(cali_start_time, cali_stop_time, wavelength, polly_type, telescope, cali_method)' elif parameter == 'DC': table_name = 'depol_calibration_constant' column_names = [ 'cali_start_time', 'cali_stop_time', 'depol_const', 'uncertainty_depol_const', 'used_for_processing', 'wavelength', 'telescope', 'nc_zip_file', 'polly_type'] data_types = ['text', 'text', 'real', 'real', 'integer', 'text', 'text', 'text', 'text'] + unique=', UNIQUE(cali_start_time, cali_stop_time, wavelength, polly_type, telescope)' assert len(column_names) == len(data_types), 'column names do not match data types' logging.info(f'writing to sqlite-db: {db_path}') logging.info(f'writing {parameter} to table: {table_name}') - sql_db.setup_empty(db_path, table_name, column_names, data_types) + sql_db.setup_empty(db_path, table_name, column_names, data_types, unique=unique) rows_to_insert = sql_db.prepare_for_sql_db_writing(self, parameter, method) sql_db.write_rows_to_sql_db(db_path, table_name, column_names, rows_to_insert) + def read_calibration_db(self, db_path:str|None=None): + """read the calibration constants from database + + time interval includes 24h before and after the actual date + + """ + if db_path == None: + db_path = self.polly_config_dict['calibrationDB'] + logging.info(f"read db_path from polly_config_dict {db_path}") + + ts_interval = self.retrievals_highres['time'][0], self.retrievals_highres['time'][-1] + table_name = 'lidar_calibration_constant' + self.LC.update(sql_db.get_from_sql_db(db_path, table_name, ts_interval)) + + table_name = 'depol_calibration_constant' + self.pol_cali.update(sql_db.get_from_sql_db(db_path, table_name, ts_interval)) + + def adding_retrieving_infos_2_polly_config_dict(self): """ some infos from the polly_config_dict should have there own keys, e.g. reference_search_range diff --git a/ppcpy/io/sql_interaction.py b/ppcpy/io/sql_interaction.py index 7ab1935..57b08aa 100644 --- a/ppcpy/io/sql_interaction.py +++ b/ppcpy/io/sql_interaction.py @@ -1,66 +1,102 @@ import sqlite3 import logging -from datetime import datetime, timezone -import ppcpy.misc.helper as helper +from datetime import datetime, timezone, timedelta +from collections import defaultdict +import pandas as pd -def get_LC_from_sql_db(db_path:str, table_name:str, wavelength:int|str, method:str, telescope:str, timestamp:str) -> dict: - """ - Accesses the sqlite db table and returns LC for all cloud-free-regions )profiles) - - Parameters: - - db_path (str): name of the specific sqlite db file. - - table_name (str): default 'lidar_calibration_constant' - - wavelength (int or str): the wavelength - - method (str): Klett or Raman - - telescope (str): NR or FR - - timestamp (str): the date or timestamp to look for - Output: - - LC (dict): containing all profiles as list +import ppcpy.misc.helper as helper +from ppcpy.misc.helper import default_to_regular + +mapping = {'far_range': 'FR', 'near_range': 'NR', 'dfov': 'DFOV'} +mapping_inverse = {y: x for x, y in mapping.items()} + +def string_to_ts(s): + """string of format %Y-%m-%d %H:%M:%S to timestamp (timezone-aware)""" + return datetime.strptime(s, "%Y-%m-%d %H:%M:%S").replace(tzinfo=timezone.utc).timestamp() + +def get_from_sql_db(db_path:str, table_name:str, ts_interval:list[str]) -> dict: + """read lidar calibration constant or depol calibration from database + + Parameters + ---------- + db_path : str + name of the specific sqlite db file. + table_name : str + default 'lidar_calibration_constant' + ts_interval : str + the date or timestamp to look for + + + Returns + ------- + dict + in calibration storage format """ - timestamp = datetime.strptime(timestamp, "%Y%m%d").strftime("%Y-%m-%d") - - if telescope == 'FR': - telescope_db = 'far_range' - elif telescope == 'NR': - telescope_db = 'near_range' - - conn = sqlite3.connect(db_path) - cursor = conn.cursor() - - #cursor.execute("SELECT name FROM sqlite_master WHERE type='table';") - #cursor.execute(f"PRAGMA table_info({table_name});") - #print(cursor.fetchall()) - - query = f""" - SELECT cali_start_time,cali_stop_time,liconst - FROM {table_name} - WHERE wavelength LIKE ? AND - cali_method LIKE ? AND - telescope LIKE ? AND - cali_start_time LIKE ? - """ - - params = (f'%{wavelength}%', f'%{method}%', f'%{telescope_db}%',f'%{timestamp}%') - cursor.execute(query, params) - rows = cursor.fetchall() - LC = {} - LC[f'{wavelength}_total_{telescope}'] = [] - for row in rows: - LC[f'{wavelength}_total_{telescope}'].append(row) - + delta = timedelta(hours=24) + start = ( + datetime.fromtimestamp(ts_interval[0], timezone.utc) - delta + ).strftime("%Y-%m-%d %H:%M") + end = ( + datetime.fromtimestamp(ts_interval[0], timezone.utc) + delta + ).strftime("%Y-%m-%d %H:%M") + with sqlite3.connect(db_path) as conn: + df = pd.read_sql_query( + f'SELECT * FROM {table_name} WHERE cali_start_time BETWEEN ? AND ?;', + conn, params=(start, end)) conn.close() - return LC + + ret = {} + + if table_name == 'lidar_calibration_constant': + d = defaultdict(list) + for index, row in df[df.cali_method == 'Raman_Method'].iterrows(): + k = f"{row['wavelength']}_total_{mapping[row['telescope']]}" + d[k].append({ + 'LC': row['liconst'], 'LCStd': row['uncertainty_liconst'], + 'time_start': int(string_to_ts(row['cali_start_time'])), + 'time_end': int(string_to_ts(row['cali_stop_time'])), + }) + ret['raman_db'] = default_to_regular(d) + + d = defaultdict(list) + for index, row in df[df.cali_method == 'Klett_Method'].iterrows(): + k = f"{row['wavelength']}_total_{mapping[row['telescope']]}" + d[k].append({ + 'LC': row['liconst'], 'LCStd': row['uncertainty_liconst'], + 'time_start': int(string_to_ts(row['cali_start_time'])), + 'time_end': int(string_to_ts(row['cali_stop_time'])), + }) + + ret['klett_db'] = default_to_regular(d) + + if table_name == 'depol_calibration_constant': + d = defaultdict(list) + for index, row in df.iterrows(): + k = f"{row['wavelength']}_{mapping[row['telescope']]}" + d[k].append({ + 'eta': row['depol_const'], 'eta_std': row['uncertainty_depol_const'], + 'time_start': int(string_to_ts(row['cali_start_time'])), + 'time_end': int(string_to_ts(row['cali_stop_time'])), + }) + ret['D90_db'] = default_to_regular(d) + + return ret def prepare_for_sql_db_writing(data_cube, parameter:str, method:str) -> list[tuple]: """ Collect all necessary variable and save it to a list of tuples for inserting into a SQLite table. - Parameters: - - data_cube (object) - - parameter (str): LC or DC - - method (str): klett or raman - Output: - - rows_to_insert (list of tuples) + Parameters + ---------- + data_cube : object + parameter :str + LC or DC + method : str + klett or raman + + Returns + ------- + rows_to_insert : list of tuples """ rows_to_insert = [] @@ -70,60 +106,57 @@ def prepare_for_sql_db_writing(data_cube, parameter:str, method:str) -> list[tup method_db = 'Klett_Method' if parameter == 'LC': - for n in range(0, len(data_cube.clFreeGrps)): - starttime = data_cube.retrievals_highres['time64'][data_cube.clFreeGrps[n][0]] - stoptime = data_cube.retrievals_highres['time64'][data_cube.clFreeGrps[n][1]] - start = starttime.astype('datetime64[ms]').item().strftime("%Y-%m-%d %H:%M:%S") - stop = stoptime.astype('datetime64[ms]').item().strftime("%Y-%m-%d %H:%M:%S") - for ch in data_cube.LC[method][n].keys(): - wv, pol, tel = helper.get_wv_pol_telescope_from_dictkeyname(ch) - if tel == 'FR': - tel_db = 'far_range' - elif tel == 'NR': - tel_db = 'near_range' - LC = data_cube.LC[method][n][ch]['LC'] - LC_std = data_cube.LC[method][n][ch]['LCStd'] - LC_is_used = True if LC == data_cube.LCused[ch] else False - rows_to_insert.append( - (str(start), str(stop), float(LC), float(LC_std), LC_is_used, - wv, str(data_cube.rawfile), data_cube.device, method_db, tel_db)) + for e in data_cube.LC[method].keys(): + wv, pol, tel = helper.get_wv_pol_telescope_from_dictkeyname(e) + tel_db = mapping_inverse[tel] + for line in data_cube.LC[method][e]: + LC = line['LC'] + LC_std = line['LCStd'] + LC_is_used = True if LC == data_cube.LCused[e] else False + start_unix = line['time_start'] + stop_unix = line['time_end'] + start = datetime.fromtimestamp(start_unix, timezone.utc).strftime("%Y-%m-%d %H:%M:%S") + stop = datetime.fromtimestamp(stop_unix, timezone.utc).strftime("%Y-%m-%d %H:%M:%S") + rows_to_insert.append(( + str(start), str(stop), float(LC), float(LC_std), LC_is_used, + wv, str(data_cube.rawfile), data_cube.device, method_db, tel_db)) + elif parameter == 'DC': - for e in data_cube.pol_cali.keys(): + for e in data_cube.pol_cali['D90'].keys(): wv, tel = e.split('_') - if tel == 'FR': - tel_db = 'far_range' - elif tel == 'NR': - tel_db = 'near_range' - elif tel == 'DFOV': - tel_db = 'dfov' - for i in range(len(data_cube.pol_cali[e]['eta'])): - eta = data_cube.pol_cali[e]['eta'][i] - eta_std = data_cube.pol_cali[e]['eta_std'][i] - eta_is_used = True if eta == data_cube.pol_cali[e]['eta_best'] else False - start_unix = data_cube.pol_cali[e]['time_start'][i] - stop_unix = data_cube.pol_cali[e]['time_end'][i] + tel_db = mapping_inverse[tel] + for line in data_cube.pol_cali['D90'][e]: + eta = line['eta'] + eta_std = line['eta_std'] + eta_is_used = True if eta == data_cube.etaused[e] else False + start_unix = line['time_start'] + stop_unix = line['time_end'] start = datetime.fromtimestamp(start_unix, timezone.utc).strftime("%Y-%m-%d %H:%M:%S") stop = datetime.fromtimestamp(stop_unix, timezone.utc).strftime("%Y-%m-%d %H:%M:%S") - rows_to_insert.append( - (str(start), str(stop), float(eta), float(eta_std), eta_is_used, - wv, tel_db, str(data_cube.rawfile), data_cube.device)) + rows_to_insert.append(( + str(start), str(stop), float(eta), float(eta_std), eta_is_used, + wv, tel_db, str(data_cube.rawfile), data_cube.device)) return rows_to_insert -def setup_empty(db_path:str, table_name:str, column_names:list[str], data_types:list[str]): - """ - Create/Initialise an empty database. - - Parameters: - - db_path (str): Path to the SQLite database file. - - table_name (str): Name of the target table. - - column_names (list of str): List of column names to insert values into (e.g. ['col1', 'col2']). - - data_types (list of str): List of SQLite data types for each respective columns (e.g. ['text', 'real']) +def setup_empty(db_path:str, table_name:str, column_names:list[str], data_types:list[str], unique:str=''): + """Create/Initialise an empty database. + + Parameters + ---------- + db_path : str + Path to the SQLite database file. + table_name : str + Name of the target table. + column_names : list of str + List of column names to insert values into (e.g. ['col1', 'col2']). + data_types : list of str + List of SQLite data types for each respective columns (e.g. ['text', 'real']) """ column_names = ['id'] + column_names data_types = ['INTEGER PRIMARY KEY'] + data_types columns = ', '.join([f"{c} {d}" for c, d in zip(column_names, data_types)]) - sql = f"CREATE TABLE IF NOT EXISTS {table_name} ({columns}) " + sql = f"CREATE TABLE IF NOT EXISTS {table_name} ({columns}{unique}) " with sqlite3.connect(db_path) as conn: cursor = conn.cursor() cursor.execute(sql) @@ -132,23 +165,32 @@ def setup_empty(db_path:str, table_name:str, column_names:list[str], data_types: def write_rows_to_sql_db(db_path:str, table_name:str, column_names:list[str], rows_to_insert:list[str]): - """ - Insert multiple rows into a SQLite table. - - Parameters: - - db_path (str): Path to the SQLite database file. - - table_name (str): Name of the target table. - - column_names (list of str): List of column names to insert values into (e.g. ['col1', 'col2']). - - rows_to_insert (list of tuples): Data to insert, e.g. [('a', 'b'), ('c', 'd')]. + """Insert multiple rows into a SQLite table. + + Parameters + ---------- + db_path : str + Path to the SQLite database file. + table_name : str + Name of the target table. + column_names : list of str + List of column names to insert values into (e.g. ['col1', 'col2']). + rows_to_insert : list of tuples + Data to insert, e.g. [('a', 'b'), ('c', 'd')]. + + + Notes + ----- + The IGNORE syntax somehow did not work. + With the UNIQUE colums defined and INSERT OR REPLACE at least the new values are updated. + Though they are given a new ID. + """ placeholders = ', '.join(['?'] * len(column_names)) columns = ', '.join(column_names) - sql = f"INSERT OR IGNORE INTO {table_name} ({columns}) VALUES ({placeholders})" - ## IGNORE means: skipping rows with - ## identical entries for 'cali_start_time' & 'cali_stop_time' & 'wavelength' & - ## 'polly_type' & 'cali_method' & 'telescope' which are already in the db - # MR: not sure if newer values should be actually overwritten + #sql = f"INSERT INTO {table_name} ({columns}) VALUES ({placeholders}) ON CONFLICT(cali_start_time, cali_stop_time, wavelength, polly_type, telescope) DO UPDATE SET data = excluded.data" + sql = f"INSERT OR REPLACE INTO {table_name} ({columns}) VALUES ({placeholders})" try: with sqlite3.connect(db_path) as conn: diff --git a/ppcpy/io/write2nc.py b/ppcpy/io/write2nc.py index 017b4f3..5351e99 100644 --- a/ppcpy/io/write2nc.py +++ b/ppcpy/io/write2nc.py @@ -80,7 +80,7 @@ def write_channelwise_2_nc_file(data_cube, root_dir=root_dir, prod_ls=[]): ## update variable attribute if "eta" in json_nc_mapping_dict['variables'][v]['attributes'].keys(): wv, t, tel = re.findall(r"(\d{3,4})_(\w+)_(\w+)", v)[0] - json_nc_mapping_dict['variables'][v]['attributes']['eta'] = data_cube.pol_cali[f'{wv}_{tel}']['eta_best'] + json_nc_mapping_dict['variables'][v]['attributes']['eta'] = data_cube.etaused[f'{wv}_{tel}'] if "Lidar_calibration_constant_used" in json_nc_mapping_dict['variables'][v]['attributes'].keys(): LC_used_key = v.split("attBsc_")[-1] json_nc_mapping_dict['variables'][v]['attributes']['Lidar_calibration_constant_used'] = data_cube.LCused[LC_used_key] diff --git a/ppcpy/misc/helper.py b/ppcpy/misc/helper.py index 7d3113b..d848b84 100644 --- a/ppcpy/misc/helper.py +++ b/ppcpy/misc/helper.py @@ -1,3 +1,6 @@ + + +from collections import defaultdict import sys import os import re @@ -913,3 +916,12 @@ def idx2time(cldFreeIdx:np.ndarray[int, int], nIdx:int, nHour:int) -> str: cldFreeMin[i] = '0' + cldFreeMin[i] out = cldFreeHour[0] + cldFreeMin[0] + '_' + cldFreeHour[1] + cldFreeMin[1] return out + + + + +def default_to_regular(d): + """defaultdict to regular dict """ + if isinstance(d, defaultdict): + d = {k: default_to_regular(v) for k, v in d.items()} + return d \ No newline at end of file diff --git a/ppcpy/qc/transCor.py b/ppcpy/qc/transCor.py index f06cf3c..46ccf1d 100644 --- a/ppcpy/qc/transCor.py +++ b/ppcpy/qc/transCor.py @@ -30,14 +30,14 @@ def transCorGHK_cube(data_cube, signal='BGCor'): print('G', config_dict['G'][flagt], config_dict['G'][flagc]) print('H', config_dict['H'][flagt], config_dict['H'][flagc]) - print('polCaliEta', data_cube.pol_cali[f'{wv}_{tel}']['eta_best']) + print('polCaliEta', data_cube.etaused[f'{wv}_{tel}']) # similar to voldepol_2d vdr, vdrStd = depolarization.calc_profile_vdr( sigBGCor_total, sigBGCor_cross, config_dict['G'][flagt], config_dict['G'][flagc], config_dict['H'][flagt], config_dict['H'][flagc], - data_cube.pol_cali[f'{wv}_{tel}']['eta_best'], config_dict[f'voldepol_error_{wv}'], + data_cube.etaused[f'{wv}_{tel}'], config_dict[f'voldepol_error_{wv}'], ) sigTCor_total, bgTCor_total = transCor_E16_channel( diff --git a/ppcpy/retrievals/depolarization.py b/ppcpy/retrievals/depolarization.py index de37aab..fbff880 100644 --- a/ppcpy/retrievals/depolarization.py +++ b/ppcpy/retrievals/depolarization.py @@ -63,12 +63,12 @@ def voldepol_cldFreeGrps(data_cube, ret_prof_name): # data_cube.retrievals_highres[f'BG{signal}'][slice(*cldFree),data_cube.gf(wv, 'total', tel)]), axis=0) sigc = np.squeeze(data_cube.retrievals_profile[f'sig{signal}'][i,:,flagc]) - print(channel, data_cube.pol_cali[f'{wv}_{tel}']['eta_best']) + print(channel, data_cube.etaused[f'{wv}_{tel}']) vdr, vdrStd = calc_profile_vdr( sigt, sigc, config_dict['G'][flagt], config_dict['G'][flagc], config_dict['H'][flagt], config_dict['H'][flagc], - data_cube.pol_cali[f'{wv}_{tel}']['eta_best'], config_dict[f'voldepol_error_{wv}'], + data_cube.etaused[f'{wv}_{tel}'], config_dict[f'voldepol_error_{wv}'], window=config_dict[f'smoothWin_{retrieval}_{wv}'] ) opt_profiles[i][channel]['vdr'] = vdr @@ -90,7 +90,7 @@ def voldepol_cldFreeGrps(data_cube, ret_prof_name): vdr, vdrStd = calc_profile_vdr( sigt, sigc, config_dict['G'][flagt], config_dict['G'][flagc], config_dict['H'][flagt], config_dict['H'][flagc], - data_cube.pol_cali[f'{wv}_{tel}']['eta_best'], config_dict[f'voldepol_error_{wv}'], + data_cube.etaused[f'{wv}_{tel}'], config_dict[f'voldepol_error_{wv}'], window=1 ) mdr, mdrStd, flgaDeftMdr = get_MDR( diff --git a/ppcpy/retrievals/highres.py b/ppcpy/retrievals/highres.py index cd71cc5..42d83be 100644 --- a/ppcpy/retrievals/highres.py +++ b/ppcpy/retrievals/highres.py @@ -109,7 +109,7 @@ def voldepol_2d(data_cube): vdr, vdrStd = depolarization.calc_profile_vdr( sigt, sigc, config_dict['G'][flagt], config_dict['G'][flagc], config_dict['H'][flagt], config_dict['H'][flagc], - data_cube.pol_cali[f'{wv}_{tel}']['eta_best'], config_dict[f'voldepol_error_{wv}'], + data_cube.etaused[f'{wv}_{tel}'], config_dict[f'voldepol_error_{wv}'], window=1) vdr[data_cube.retrievals_highres['depCalMask'], :] = np.nan data_cube.retrievals_highres[f"voldepol_{wv}_total_{tel}"] = vdr diff --git a/ppcpy/retrievals/quasi.py b/ppcpy/retrievals/quasi.py index 89d82dd..4491b97 100644 --- a/ppcpy/retrievals/quasi.py +++ b/ppcpy/retrievals/quasi.py @@ -42,7 +42,7 @@ def quasi_pdr(data_cube, wvs=[532], version='V1'): vdr, _ = depolarization.calc_profile_vdr( sigt, sigc, config_dict['G'][flagt], config_dict['G'][flagc], config_dict['H'][flagt], config_dict['H'][flagc], - data_cube.pol_cali[f'{wv}_{tel}']['eta_best'], config_dict[f'voldepol_error_{wv}'], + data_cube.etaused[f'{wv}_{tel}'], config_dict[f'voldepol_error_{wv}'], window=1) if f"quasiBsc{version}_{wv}_{t}_{tel}" in data_cube.retrievals_highres.keys(): From 91a530738944b54ed1b2fb7540bb1a645b97a8d8 Mon Sep 17 00:00:00 2001 From: martin-rdz Date: Sat, 28 Mar 2026 18:41:10 +0100 Subject: [PATCH 2/3] add documentation add forgotten file fix typo --- docs/source/userguide/cal_constants.rst | 105 ++++++++++++++++++++++++ docs/source/userguide/index.rst | 2 + ppcpy/calibration/select.py | 43 ++++++++++ ppcpy/interface/picassoProc.py | 4 +- 4 files changed, 152 insertions(+), 2 deletions(-) create mode 100644 docs/source/userguide/cal_constants.rst create mode 100644 ppcpy/calibration/select.py diff --git a/docs/source/userguide/cal_constants.rst b/docs/source/userguide/cal_constants.rst new file mode 100644 index 0000000..8fbe9ba --- /dev/null +++ b/docs/source/userguide/cal_constants.rst @@ -0,0 +1,105 @@ + +********************* +Calibration Constants +********************* + +Polarization Calibration +========================= + +The delta-90 polarization calibration for the detected time periods is done with :py:meth:`~ppcpy.interface.picassoProc.PicassoProc.polarizationCaliD90`. +Resulting etas are stored in ``data_cube.pol_cali`` for all the time intervals and ``data_cube.etaused`` for the one with the lowest standard deviation. + + +.. code-block:: python + + >>> data_cube.polarizationCaliD90() + >>> data_cube.pol_cali + {'D90': + {'355_FR': [ + {'eta': 10.602253607032777, + 'eta_std': 0.36947564349830514, + 'time_start': 1755916590, + 'time_end': 1755916800, + 'status': 1}, + {'eta': 10.45194101666947, + 'eta_std': 0.27254858255926373, + 'time_start': 1755968190, + 'time_end': 1755968400, + 'status': 1}, + ... ] + >>> data_cube.etaused + {'355_FR': np.float64(10.754878969562773), + '532_FR': np.float64(6.852638549685832)} + +A method for polarization calibration at the reference height is also implemented, but not extensively tested or used operationally: :py:meth:`~ppcpy.interface.picassoProc.PicassoProc.polarizationCaliMol` + + +Lidar (absolute) Calibration +============================== + +After retrieval of the optical profiles, the Lidar calibration constants can be calculated with :py:meth:`~ppcpy.interface.picassoProc.PicassoProc.LidarCalibration`. +Similarly to the polarization calibration the values are stored in ``data_cube.LC`` with values with the lowest standard deviation are selected into ``data_cube.LCused``. + + +.. code-block:: python + + >>> data_cube.LidarCalibration() + >>> data_cube.LC + {'klett': { + '532_total_FR': [ + {'LC': np.float64(76651241245552.4), + 'LCStd': np.float64(38384910061.7245), + 'time_start': 1755907200, 'time_end': 1755910770}, + {'LC': np.float64(79019754266917.53), + 'LCStd': np.float64(20060545795.042507), + 'time_start': 1755910800, 'time_end': 1755916200}, + {'LC': np.float64(77383458139465.62), + 'LCStd': np.float64(30826575129.720245), + 'time_start': 1755916830, 'time_end': 1755920400}, + ... ], + '355_total_FR': [ + {'LC': np.float64(61221519590893.43), + 'LCStd': np.float64(15266830837.91351), + 'time_start': 1755907200, 'time_end': 1755910770}, + ... ]}, + 'raman': { + '355_total_FR': [ + {'LC': np.float64(57637355636083.1), + 'LCStd': np.float64(42772775199.508766), + 'time_start': 1755907200, 'time_end': 1755910770}, + {'LC': np.float64(56121416303705.17), + 'LCStd': np.float64(65519957166.51127), + 'time_start': 1755910800, 'time_end': 1755916200}, + ... ]}, + } + >>> data_cube.LCused + {'532_total_FR': np.float64(47109585570330.64), + '355_total_FR': np.float64(51991322600216.1), + '1064_total_FR': np.float64(418368070699337.2), + '532_total_NR': np.float64(6649087258385.09), + '355_total_NR': np.float64(9550053707341.717), + '387_total_FR': np.float64(51991322600216.1), + '607_total_FR': np.float64(47109585570330.64), + '607_total_NR': np.float64(6649087258385.09), + '387_total_NR': np.float64(9550053707341.717)} + + +Writing/Reading Database +========================== + +A sqlite database for storage of the retrieved calibration data is defined in ``data_cube.polly_config_dict['calibrationDB']`` (originally in the respective config file): :py:meth:`~ppcpy.interface.picassoProc.PicassoProc.write_2_sql_db` + +.. code-block:: python + + data_cube.write_2_sql_db(parameter='LC', method='raman') + data_cube.write_2_sql_db(parameter='LC', method='klett') + data_cube.write_2_sql_db(parameter='DC') + +Storing multiple times will update the calibration constants, if exactly the same time interval is used for the optical profile. +The calibration constants can also be read back into picassopy: :py:meth:`~ppcpy.interface.picassoProc.PicassoProc.read_calibration_db` + +.. code-block:: python + + data_cube.read_calibration_db() + +The values are then stored into ``data_cube.LC['raman_db']``, ``data_cube.LC['klett_db']`` and ``data_cube.pol_cali['D90_db']``, respectively. \ No newline at end of file diff --git a/docs/source/userguide/index.rst b/docs/source/userguide/index.rst index bc049c7..76c9993 100644 --- a/docs/source/userguide/index.rst +++ b/docs/source/userguide/index.rst @@ -12,5 +12,7 @@ User Guide configfiles.rst preprocess.rst + cal_constants.rst + meteodata.rst todo.rst diff --git a/ppcpy/calibration/select.py b/ppcpy/calibration/select.py new file mode 100644 index 0000000..183e5ca --- /dev/null +++ b/ppcpy/calibration/select.py @@ -0,0 +1,43 @@ + +import numpy as np + + +def single_best(d, name_val, name_min): + """select the best calibration constant + + + generalization of lidarconstant.get_best_LC + + Parameters + ---------- + d : dict + dict per channel of list of calibration constants + name_val : str + designator for the actual value + name_min : str + designator of the minimum + + + Returns + ------- + best : dict + Lidar constants/Etas with lowest standard deviation per channel. + + Notes + ----- + Since ``LC = LC_stable`` and ``LCStd = LC_stable * LC_Std`` so will any negative LC also have + a negative LCStd, and thus be chosen as the best LC. + + **History** + + - 2026-02-16: Added additional checks to hinder negative LCs to be chosen. + - 2026-03-27: generalized to also hold for depolarization calibration + """ + + best = {} + for k, l in d.items(): + val = np.array([e[name_val] for e in l if e[name_val] >= 0]) + min = np.array([e[name_min] for e in l if e[name_val] >= 0]) + best[k] = val[np.argmin(min)] + + return best \ No newline at end of file diff --git a/ppcpy/interface/picassoProc.py b/ppcpy/interface/picassoProc.py index bfa1e36..7e1181b 100644 --- a/ppcpy/interface/picassoProc.py +++ b/ppcpy/interface/picassoProc.py @@ -674,8 +674,8 @@ def quasiV2(self): def write_2_sql_db(self, parameter:str, db_path:str|None=None, method:str|None=None): """ write LC or eta to sqlite db table - Paramters - --------- + Parameters + ---------- parameter : str can be LC (Lidar-calibration-constant) or DC (Depol-calibration-constant) method : str From d066eb9daca4bfa08527d19ed79e64643ae16d7f Mon Sep 17 00:00:00 2001 From: martin-rdz Date: Sat, 28 Mar 2026 20:19:32 +0100 Subject: [PATCH 3/3] add plotting function add plot to doc fix some typos adressing #38 --- docs/source/api_doc/calibration.rst | 7 ++ docs/source/userguide/cal_constants.rst | 14 +++- .../userguide/img/plot_cals_example.png | Bin 0 -> 67052 bytes docs/source/userguide/preprocess.rst | 12 ++-- ppcpy/calibration/select.py | 68 +++++++++++++++++- 5 files changed, 93 insertions(+), 8 deletions(-) create mode 100644 docs/source/userguide/img/plot_cals_example.png diff --git a/docs/source/api_doc/calibration.rst b/docs/source/api_doc/calibration.rst index 16bf139..1583f64 100644 --- a/docs/source/api_doc/calibration.rst +++ b/docs/source/api_doc/calibration.rst @@ -24,3 +24,10 @@ ppcpy.calibration.rayleighfit .. automodule:: ppcpy.calibration.rayleighfit :members: + +************************* +ppcpy.calibration.select +************************* + +.. automodule:: ppcpy.calibration.select + :members: diff --git a/docs/source/userguide/cal_constants.rst b/docs/source/userguide/cal_constants.rst index 8fbe9ba..b95e1a0 100644 --- a/docs/source/userguide/cal_constants.rst +++ b/docs/source/userguide/cal_constants.rst @@ -102,4 +102,16 @@ The calibration constants can also be read back into picassopy: :py:meth:`~ppcpy data_cube.read_calibration_db() -The values are then stored into ``data_cube.LC['raman_db']``, ``data_cube.LC['klett_db']`` and ``data_cube.pol_cali['D90_db']``, respectively. \ No newline at end of file +The values are then stored into ``data_cube.LC['raman_db']``, ``data_cube.LC['klett_db']`` and ``data_cube.pol_cali['D90_db']``, respectively. + +A basic visualization of the calibration constants is also available by :py:func:`~ppcpy.calibration.select.plot_cals` + +.. code-block:: python + + from ppcpy.calibration.select import plot_cals + plot_cals(data_cube.pol_cali, 'eta', used=data_cube.etaused) + plot_cals(data_cube.LC, 'LC', used=data_cube.LCused) + + +.. image:: img/plot_cals_example.png + :width: 55% \ No newline at end of file diff --git a/docs/source/userguide/img/plot_cals_example.png b/docs/source/userguide/img/plot_cals_example.png new file mode 100644 index 0000000000000000000000000000000000000000..383a056ae7f0d87a7425c6fbc7d0b34357f730bc GIT binary patch literal 67052 zcmdS>WmMH|^finkAdPgFbc1wvON(?$mvl>rbW2E=bf~3cprD{o1yvn`{P5u=_#V+gR@)g03hguG-wSr}@(n1c=X^Oy zvG*PZf0x}|h^OaoPLFL@r%T)i_m7t0+^g8&;H#{@(_W~EGkm>P6h{c$lo3a=XF#Mz z7Y{+NimI;eZ2a@((>K%K`!!xk^$jVjWO2Aw9|KUi6Nk_FEJmx{kMxI!GJLk@1>BG5 zQVXmQ$OF;E5y+#%Sk6%Y`?6!f4^SmmHpojRM}Kr0Y7Y~YginG%{NF zVMRqtD}#$69?B66d7>k%L1j{ED&Z?)UQrntc;61P_(4f0C-#lLC`v{~Mm;z11u?G! zRMsfZ9}Zh?qR9e9lc^$AC|Kwh;7@+uT^_(Ez`VgnN%woY-Q3y9og`W86>C^jG%=y{ z7ybP3;2kM#Vj62QF*PNrt4kCX7M52~Ak897Gdw)(y(6&NXo0D%dY`$7Q{}1-pY$O@@-$b}#P+ zV$j7ysj!_a19OcmpFN zeRd5YB9!Fh{w+t%5kbg!;`TK+59g!jlUYv*fB)LX#l>}*Dbiqx+|Q3DgWY4J{i3d} z&ZP4r^ytXW_okZ+6Auqhdo!6yw?5MKV216})!&n5(#p|zQJYk0S64198im(1G{KRP zkuOA`**Q3FN8H=H)A^jWYD^J3Zinr8U3SsI%q&x->9C=Vz=AG$?qm$5^VO~vB9r$F zQbG;@v6!FLCO-2_NTNXiU5G8Ak+i^A_&w4s+?;HQ=xQ)l@lhZQ)1P_VyY6@vil1+qwFWMDz3WF#CvL2OVr-4#9*OOYDn8U!(s#bemjRi5@iNw@t4p zu)@Z}tF@Y>ssFerZDmEjyu3_Y9Mg4h`Od&#q(CukuGy85nwr|h-CYLv{=J!*8RXcN z>o&a8(;JHMy`Vg5Jy#<2J`r13TvS#6HZk$;=_EwNVN~eaVJE}gX;$C;bW8%R+U-z> z_!YdUNAg>BbsYQFQ#kh(PulBccQX6-TN)6W0q|&N2Zqn7ATeZt{ab__zz4`=%Ntg5 za&tFA7q2~wO5W08Q{bZn8EDDgBZiPkw&B z?R-tpH%85~A&yZMFq2+9L*gSoZhRE2at+_p2*Z;od+&hKistZR5csXY<_-_YgVtDGDffq9K{fW%AG=yM|v{9 z2ggRuy9jGbOQYEeecQ!)=7*BlOCr{E(f7vRUwh51#c1VL!2%M(a#!?F~ z5MV*ZUL$o^S1&Ft-S~pcVzkog30XT+Gc%txa$+#Eh+g(=i|6+Hqm8RYds=4ZC~En{ zla|%T=%pn+Qjb3`zI48Jb#-OZZ%*J@@d!8{70Ciw>A)E^1XPxd%Wr%B8tOQevk{Yf zJ!=ksPCcUc&WxzCMnisYb%F1lnh@XuU=nRogIAByo-^8={6&pK`N(_r$U#oYdi-nA zdbeK+nMFD{w6FnU&CShSX6@H)$EUSfkC7muB>rY;iG(Nf{whP}lK_&C@WY7&C`VAD z@84r1<1+mj^?Tv~fj|QJ8U*i_$n!1Q>FKH9*#ISNp3LmdIV3z>_(|NC9N{%exmy2% zWAR!*Ag%4w^*cjDk!*v!c5NAGe<-EG_~__Y*4EZwyeI_hfh8pjOLp$2b)#D#fvn%| zWbF*|tfF+-HLeAp#9)6mdNYt_-#&P;H_5k??%W6_;rVzu+w|uP^38eHGt*!#`nvr2 zH@WPHh=|+UQNNL)6b?CM<((0MW7FrS$Jk!sZ=*XKNgWT=czO;mw6(Q|(s+)mZuN)= z3BehE>I5$z92ubj<)Yypb3ooI)y*)jUQxg&s-+|!E%sf2bKggKgbJ^V5LQTeH z6&B9Cx$yw2?$D$MvoAHl8I_ob>wCYr%Ko7X_T9U81EZsb+1a0`OEo8FX0C!}>y(w0 zUe%fpuY<5g?l?O;gBvROK6{{DW%L4)Fr#=@{&d2VQE_lURx`y%0^v0O?ibYG#OEljSE2Nk0*|%*KOp)Sw7s<&JxZ&tBGo+6qc{ihw%{ zp1unPBzO9O%0JiPCxRw)pY8lUQ*`yV^O%q+5+cC`fFuGoBsx{ptIeRfrDc521twU9IIN^(YGy{<9-MBf zHW`M}*@MA>0a11Jl9Cd*3ru%6w~6_AIRyn2?lk6yr&!^_^71IfG;TJF5n@u1+bY_w zbXb~qDZtT*DJhY8Jau5F5C39N3055)2dB5uaYGAi;qN3?i|eE1*xn%Y^}HLfZ_>H3 zMMZP~48jGFJF(qgAMbzyLIe<<;KNDx!)c74>0H%^C9i|>grub05{B>`qlU$>@NiIY zNlycFC->f@D1c%OQX(ir;o(Riyzm9wY{7c=yqnkh~^2)3?90w}hbw(TJD020O6&A>v3MuSCt*x!TSt}-w zcb7|!{Um~q7n6_-{t_A2Vew;wkGXOJy0ODtqa!7#%&S@*zCv|%bs$>Lp6-v&L4pVZ ziOLBaNasQ7q_KsCe+VQzJRE{5;a?F$O89!yPG(`tf^cC6FrYg`=WOu4QMz%4uH<8yPeUS9lQ zir6?f)yBO@IyyR*A3vJio^3sVq}F3y=X8JAz~TGg0rCnL7Z(Ht00aiU0gw`O3E-GW zNlBIIx8Q=Vz{ZC0cOvr!7zRK@Sa^8-W=Z@vJh<-=eRO2FEEOT2cWhQ>Lfi(U&G2-6ZRBHJeF!!ZvSGeq~hkj z3MC;1T}_TGbVi&w=s3uHZ^%JpBPokqg0foWviI)lcvTrhM3SBpB27H~13}g|Q6dz$ zxcGPiVQ=2>oGP1HmcB>=2_+>ZQ&n&xr$7KXI%=~)QPbA6E*YgI^zrV}4+UN~_QE~AH(A|{j)P+NKd0YxFHaU5qCkzV$e!=QN zhLDp(adv)w3sTr`u+Q?K>Olc~>d|HZ8N>c`>Fv!%%!_0mds2|-A>m=S%rAVmO9{5F z5uk=HG-i;RJLwg(@6C}wV6K1DGyshlo6SrarGXW;cqrBuBOlsD+x3b(*f2D~6DX6u z$nVX~sU6Rc$4rLp>EO_Wta-|EDN0OoIypP@9k-stU5JW2-E2U1AA;|OG6Wh{3$GD^ z#KEAvg~4`1rUuP>5)`xw(y*g>~KCd2n!W zIJvmiuJ1ubfTZ!#>F|;gkb~58bt4uJ`vBfgO??ZZ(yc|P)c_pxyE`uiMn*^ejNMBw zP&x29%uzrq00M)MhzJoDx~H#?-Dy)6yyW!$F_KIqb9dJgOy$${kyPa@hs&;tv9U3Q z=4B|D!of^AB`a&p!^6W15}ZdaZtiy)8q&75e?K@Yl$N{&B^I=C6qJ+!h*&hlS9!9a zaKuvilX3)uvI=mhAEZ@!0)T*YM&9-few7vzj8~8_qyJH5ad8-Ej-O2q7wT-6o4K5- z|Lgmsyv=-b;GtiD{e*=4s}L9vTWV@pQBhI&C~%907#tM+)P+GEB{{`N%%%x zU0xm;BoY@lw_q~L#p@I(@QwGwX+N=`$Jgy`^ODk1ZIEgf>pz-w1tMbDPzR!ioZmn4 z62AfkC+6&;ZK> zcmoJgpyx6!1$t5lkP88L!K7i4Ek%JY{EJfIM)N@hI9CK<#E=F?<=wk35F%C}we*By ztcAR0(yIC_^yvx`vik=Is@T+m@MVLh>-CFvCU+Nmv86PD=rE8V0TJ}9M&`n9GsDEv z@o;v+yWN>rSN9EcA~jZ%1p`taa3g{aDq5oH)e4BulkLd^@tNA$c`!g31t$p#vgsTB zF%X11HT52l1o+&x_)SesR5Ua>bV%S&LHL~LjreT5{Qdj)IA;{-T-hDhB>|Qe2kl3$ zbQBWWYh3xg2rRG0kitY2_ z;WoU;%?m7ArA>gudu_R8aQTJ^GDHP3lpyTDMvTZvtLORuqOV{)1WAvY4_Lq%{xpPm zp}8#=`0sti(-Hpvb6x+($+Q1AJn#-B---~m0uD(RnQzSZ^DBG!)B-77Eb&k=spc4S0mXsCenk)kPQo;nJkM!&5(s^XA@Blp>J#&})uscQ!A$VmIJSZXJtP?#e~HHL+fJKcXa*W~%FC zibEcF(3(CCw!sO*dvUq&g8|7DnO|u)+Okh268AaW< zrCQRpEqpbp=;G<`zUhB!Re69ucvaqAW~K0kmSZVFG_5vPDnCL?weAb9ySynEhze>d zf8hYj)V>y+<9AaUJ!O&FOM-Xib|jw|8=~5x(xjJr(C^W?r|pe(H@KUq7y|J*ITcOh z<}t4t4tE4?Xc`)!mlkG*%0DtM74^iI?9Hh(_-70YP`e?L!@ZQt9(_H0WIi)pVIqD8 z%QDyzh=E8M$}Wvz5xr6NqqdRi8PgVyXfdsRa_bY*{wq`)jN1C|gq9$M>yz z$<8)Vy^7suYTMm}qAH~OsG`V+6LGcm0fG*DnDg$Vjk5T7?sm1Qn9FQI_%y}l@vdKm z$5{A40jc~2tqErT0B8FDOfghx-0?8fhHV$6_5A~;eNy}IOkcTe(c({JDb}hIbDWos zFVoXJS;BK9H6TkTUbuVHS7e|4{q*?gk%^Qo?IXvrX?K}`hcI_)-H5d`Riu%0A4-@$ z-n`ZqccE^|t(O{+n*KE3sB zTt0fi+Of)twMPb%k#j&$!M3SaV4*?(OY`a5aF2 zkio_IuS&=}PGIe>uNyV$Sop0{n-@5f8x%1`UUKQDSIsFdj?`D zArue`OEgdf`l5*X0d+_b5t1@&4rm!!1qILnC4x7lQAqwgZJiIe(;8ktgvV0KIDt+v z(C~a`wBN4be{Yksz#gKDxHvQc0f9%`%E>2y8BRT_bjo3V`Qmo0#;hL8RL;tcciqD`>PL@!7tBPFz^e4)kCUI(%eA2@1l4Zbb3iEFb&!+xwmqW%}M(0m>V~ z`TeaOJ-WvGr>h0toI#ckv^OSwC=lwLMln^=+L}SFP-zQL0>Tgomz2imOauN7K}65& z=436BfCDHuD41AS1R}nIfZwFlH1y#D9q{ARmwx~p-H=pjl9Tytm?LO%&lF%-LD2d@ z$hze1jQ?afLlV$QK%}K=lYXb;wrhY7DI*k`Wdh?rc*ZaJ~g;}3L?gzT@7$PU<*zFYs z|Ee#MrWsC#7l(fMi&=k6{?!LG{ohGgwpL+fpff1slFNpNQS5GVV;*rO`!MpG^7Lv^ zuIbp|-BCwO;G%tZG8gPlYF#6@?N|96FgJ-w2k#s5kr?3|4z7qp4%{l+1stXA{d2utOTY6QyC z*SJXMWtqU{W{T2(b@=b@&)YTYgfiRuXyK&7#avG(8AiX=`99}NDy05V>4@mgx(_hl zkB$cR0d=tu2j1XBQyt6Fx{KS9qM(dPB|#aN0jVIJxXh#e!TJiQZJ z_M(CW;c;Mm!JA@J@_YYJgyeY`buAL9HDW{d#q5YMI#u`zFFNek#wr4Ai zM-oNO<|wKUKt+Zev=Bn_+1Z;I;xL>0SAs3>So=rP$Z+G(V&eMx9&SfnweGBmv6O6E z+yB+`uJ>PKJU3+%^cS0!CI-_{7$tI|P{O#+>TDUE+g$0k^0lO+dy+G7{y_`fI8ij{ za3z_gZ<^W3kv|>}Joyd$V3bucGIPWi41K;`;Tecp)M7aAiDPGQrT~Ub#f5*pqKWWx zI?nS}p$?+5Qp|QrF?aiHUHrwg{CIyart_&3W=kNt1$LNaK%grizmMAR)RbP4u(pa0 z4>zbd|1NiFKV#Tl>y~@!ysnNWA<|X^h=G1FPmT%@JHe4Ct*mXIPoF-GJ&sneIR;3DlXDmWm%7*I;xk@Hjvk~B_K!wPAI0oTi<=Tcdpu`2{6xWX~VHnuP9?b zepPkhAM$LP+ZpN;zSkSd;_#e#k7I8|AI)XsXg4RD{auyR^flnI6ISx(eOHM8?qmk+ zIn#`$maH?=cUF6$jzpz-W&UIUC!cLV7VUyypv{@~khLw+eU$Wl0}+&ux8OM55arb` zMz)*%_rhI{OXwR2uRMRt2x@mq!@_R=y9< z2;5v4&%^fYo0Wd#8`mA(!WP^Ws6R(V2 zo3BJTTmO;s`3VhHfd7rbGwKqWNbudHp;0nF2?**RZ6_4%zkUr^gJ60q3n2IS`)|>j zhA%UsdjyVOA022KXK;oUSxw3|9&#JwGX2!l4su>@X1(m83&dY~2zfPrwd6dt>@E_8 zCL}q@>L>U3{J;yTOrh52BTRL{qsA2 z7@2$P9OF@T(8W=Z;zYB=yv?>ekRw7fIfxE!AUQ%cQX0FKz{HwCbtf~zN|t6N$M_3s zRP1yRx}n%pDj_(q8AsP5pULN8D>?{vhTU_9l7lMC(i1zL4)bWo+b7stgM#RGhSo0a z1c*{$66!zE@12}$L}g_Vn5baIpWxzQStfvvkM|PND8UjDsw~?LQOe}RL`UmIvT?vj zc=3N98KD$DRY5q#$+c<}e!@R)7w4B>IwBzfpH*=s1#HAKY~$=D#+6pA;X>^W^-l zQs?FVfvONGfAFT--irniY=fEg6IBbvpEwz!MY#bgGX0;M!CEsw)>>78Mg^O*fr$BH z5V0s``#Sg5J9{=4KMIm04B7k>aa-@NFv6ou>_=dqCd(|;+Cs`*Su{9Hm{ccJ3r2p{ z`UpJj?+@(LDj5E)f%Vy~>&9|5q5G_gnzMUp^!0zwD;epGjf!S{p*6@zckyRXYVh$N zu?KZ^e2xjWg8Wx{EImerR2&bD4rY$keI!UUFdV_1sFd&AEB=iYNT= zK{vM9YxJ_JO_Y(8Huz{N(XAnYg)A$!IVB2{3N?ByPP-~1>Yr^ zG{(kk+B9p0`_o7Ie} zG?Xhj9G9Ou4oPL3(!uKU3J|-mV9S+-qD+QUy9HW$tTr~KBu}^0u#NFpdiXstkA`(E z{f-QO2w!4@S4rI?vtrs>(^!w#<`kr(zTGT-7uq>!sIPb*!8+dN%~$d}h>a3P!j`EW zhT5g|ls_6xB!TMAnLUrnL{_eIVq(C-9+~pp;Tjb`{jZ)WNf=f5tDKk!fhOlSj7W}} zdeN>El5tnZ>tEWAornF^_1d0cXzX<5%0`2Hjy&gRloRKul^CD)trrZvQYy@TJGi?| zVq|=tBo4B9egT3*Vtv+1;duavFqpU7DP70c=S$9`Y<^GPfDa$f=pjw$6eq(Ap%)

E494mJ5@;o$Ch)5z_4;rTsQ8D*02+!f;lI%sRXOI4+tMRsB=w;a; z{^aDzH_bl+9_SKInA=Nwauw0>BgSHfM}3CxWq-7!C6LqiQ1`rH8tf~hBRg;Y6}r$I zTSPD!*0f8J5Z~x`>FAj3hsa=G={Vq6qB?PeEx@Sj4K@MAJ*E3Gnx zpqdq-e;Ec5dlcoFK;82){C4ri@p5gB8Imppg z%}7adki0TSE%MN7^-ELdZBaSu=w`9E%1-Y58jL86+k5wSwxd6@(v?KPQIP|-u%7IE z&o(H4RmKh|@p84wKpgkB0C3o?gqXKE6uC zprioEG}?dr@=L)Qw>AQONTIH9TP#>NT%_+Wp<^R@U%PumAc5y+8x0MJT2%jGW{`Gm zp+rV65|RWp%g3AejyFYn8(B3Wvv{*Fs+Zdn6QH?-%AY=z;^|C%kto;^L4nN-HOxj_ zz7vRFO?l$itopbAijv_iI)Bqj*P!A@B8Q`W3DUO_fl&An4&zB;0{XsN_R9cEVY1Fd z+;p2+3|p&Cc5BYQA3p|^L0hKQCuq!@IG+1q*}6St)66S8Y>AiSWi?I3uevu)A4Klz zXEJZTv78>boaa+P7~)e8OIzDK&}i@PTLU)yM3LXI#`bS-zs+yXFq9Q6YGxLqN zkVz=OK*Camlxy{i2(GTlNd z?{W(;7!hz;V|UyyA1CuV8k^!ky;nEe8Oeaiji$i;1(a~}k@UgChSd?^vWc;hj;g$1 zA{I{<(@W!-ysGqDV`PxN^!gM`L`8$?bK2XZ@AXU$)9I_IK1ule5?05%cT~Qpm|3*`!g~_DNXW?A!14oZT%>3s4{rcxZn``6 z9oRlX(v|nP?tpabz)QU3GO4KNx%JKa(8h2BXqu9hv$fXKB?Ifg(wUw!`GJ5dF8O@N zIEzmr&I0sR661SSk;Gybx6ak0sq};A1b;y14s=>zs$6$!YJ9pofHb*~N!oWjUKI*Y zmjIO&NFl)a06Yf30s#a*Ia%3kz$(=6{u33z7+`>AAq~i$^a|{kF3p*>X|J3u-@?tm!hR#tnTo&Mwhy)ghg@G3|gKT5hlA*E~CeRcTZ_Ppt|k0^!L zk%G13!O8R#P-d#wW}>XFtIz7{)K0tgz-r?|C6HL)60rWGs3IU+4Ge97 z4g#7^Ttb2*P*H*O3}G~2;9p7NqD4eU_Z?rOCD)G}T(18U74f6`UTj@@2?1-4LIywA z13x&h@E1n^d2mmqGA>~G0daA0;gYUoanq(kCkd3l>S`H3CoMZW7NBPWr&-{tjiqIF zOG{c=Ss7(?)(lwU%exn>|4|@HdGVQ|Xe9gzfHiAZ2?9`f}Q`^V^ovW^1&D+Gn0t1AeZ zo{sJh!~*~^x7oK{;zB$-K!!t3TC6ccKJIuH23rROAzL7{8pP=e_^|Wc+!$br8in|7 z66n?b?|jXp7{kTLS5=_gt+ahQ^$-N289Oin0D~z3OD)ix!9GK@;n0|vkd+kPtH1LS zo}PRqB|vNMzB%0h)~Dm^JsDcGcNqV+By0OzP<9JYm&!Q*<8baNycI;ekR$8*Z+UuB zB>$h;p?RH)8A|#C_FZzk+8nk=7U%q9wKu&HTQep@)!b}bow@ML*g zr~OvI7*}H7oKlf5kI0wdgo!DRBAhCnogZ|iX_&uVzgXmVB1ZUp!g~ppf0}1Sa#r9N z-H8_@0BDG~HA3gVZO;F6PPXuTDOgXyG2nSi9bBQUcvbt{BWMh(81s}`;hyxpSua){ zu=oi`^*_dLxBaX>UGWX66@KtsCulbCMqV_0_V=8kOZ7d+*FSDwoJ4%BL=Ra_i&u)> z?kBG+xNO3j{<2Sy-thsMZ%jtUB!f40NaTU^zSYyE{f+^QV6PWL9ECJBxbg=|e~Eh0 zd3#@rbH`xUb!yXL)FYEeM+WMwH{N7=lS}&*e`f(rm|}#Og^2^g@)ubT-b)5EXz*hQs-Ic%0avwlOGGFTB2d{AlSf zx@in~AZKL~c!6&~SV1Evar@Zkt&B{ls4#6<91?Oc!gQ)J+OGvxKfHA2fgV$is#|Nl7=2H`p8r0;3PI>ytk&n(C+^%mX+zq z<*=pNOja0?`b?Wg2aQX&;D+~K zk1x^FRDQb6lBX^YrweTAg&B&coAQ_9NElP-d=e;FUX(3eP&OZAskyfA*3o`Dbpa6M zGEs*qLO&wb{~foJbhAD1-dbk~DSg&ZX8fdkg=rtQ^CcY_pHzgnWjwmqhyfuV}+nGa&V9Avv%5K(*w=fylEONYAS5W&fGm zp6FiZ3Zt~aNo-)?WTI}4amJ)rPucvI-<&Ivjl&aKm#E<__|?-|sE*Y9NJUKSTx`RrRQWcox`Y#n1*T81+DN~@@d?#ZRZV7GCK9-``E zb_Cuyp+*|+Ws!F1=s1+!eXZYln#StL^i8ez6P5>@wuS8j|8(>U;d(Mq>||?bz$WF_ zl(MHYUGJT_R<~I3RBvn5nAE^CVTq%>0*Jghk+twVisyus_w$IQqy_f|EBJnxXHTa! zP}Ic4v~JQMjw{U6hTT<+7Dq!sXTo7uKH+1SOJUT9q!xqqB=i8Ss{Ej%t>bZkYj}|( zc9ASyg0a){w_Mt94f*ky`6)MMi)j!q@DeQLOo2q*?6d5x^7sON;nklLffHDvDVN!* zqcQ$^PP-3YcP{LI2_UpmApmVljfguTrEGZxaUP51c@Lv8vbR}eu^#;-$;e0!75Wj* z!>JrcPL^M&{SgPT90I*tD?R2#Y+prxI6pyY`<(j0Svx7cdwZr$K{-pQK!^j7%IohZ zIrYf=A~(9x_GII;6?TLJoQ} z%%Jk1R;>q$h_`sBO1A>7UVqdR)tM+M83GtX48oLRmzjHv04Lm6_e6Q&=m?rKE>#YT z9A#NiaM9!AulXK*XVktjCh%yQxp@NZ#~%fBNNN)WZ}?oJ*Ah8N!53ip^!Wb1Wf^I~*ZkRg_G_zE|Bl-%!%RUTWguhY*i^DB6i|$zVn6+M z(fxy^nN0eOMt+<8+~~l533P$_S{f_}Rj)nWu}|=w)mJ~x9np??;4G1(b-6j?ohAv^ zI~s@Cn-(fxL^j!5v}w$u%KjvtgmP;XsmtZHPlOCQ`elbx8)RIWq9&A2;lBp<*O9-% z-oXx5SvV|*9^CXs8lL?80e^FvKz4p|5pvtPNj5Mj6Dp=~e)k$*5)UT`Sc_St*;2;l zbI0G~Odmw)c}vU2yxm#Am^~Yal)0Jzk#)LJ`)Pd^jrRF#EFvYOG{GSb?*ynIRE|d~ z+Y>wBCBnhFz%9`*rgxqmypIL%{SG@w(q=SyxAxsAg#e8%%O@nZfjum&@%ckoYa7J6 zdTAD#srj&u5~PgbU^YHJr&9t6Rcy>M#7Je*QEe*_ZM@^Hap zY;GC(_37_f_E~-3DV)K4>*D-+wlX_bb^+(wdaI`wdmC#(c+`ZOXvga=79K zZT-6D#*xj;!i0biI{!!DxzD9u0*O@cvTXD6KofxmORd}YzYG4Jn+QWMU&c=N;z5MD zbwpC z3NkFrd(SxX6Vk6&B-PZ{E82C|XAui5Wj@{Y6a3tBVF&sp#Iy~NzPrLGH=ejRmuGnR zuE#A-`oYY1eqyNIoSZ7fEa@)h2BS1i`kHWbaqki?E}NUnM8kNZjA2)aZMH$|5D8VS zd4rysX4gms;$D)Cqnoblx1JFZ5G0#!GxI9aP^~b|VX+*uLl~G5)+$<#jDZ%u1kBjU z!`lNCA^&8{-4tbY6gxys6rKp+>EkVSR;9j+V7E_pgRhq!6$ z_d;KnfDmcF(MBsFy)&l*W@tF`_o0dXv+*_YMQ^1<$VrlBb}<7`w}x&pJED7c@|rp}AY z-g9`V$LlTB&C*`=bP6uz9hkMv-1#an00Q$GSu|7YeqbC z{vHK>0Zmc^YyHO`A>i%oCo4RH7s+1R5h*;~K~udBARBvFTwN1--G8}PR#o@4Hn zOMgo*)m6BIUiI{7JyC%1GF!Wbw84DaAn2Encv;^2(HepXeFomzCid}g;A12ad z>!~@wbaF9eRHrl~%X1;hYd~T_3Ga%oSoU}SWO=+ehm9zp!*vRO;yzt{fUDNj6RT

W*dJWPbW2iWZ~Wy%$!3Q~1BTQP~93cvJ57&Vf!pQFk$L`rEF?dAh1}Znizkub9L( z?M#t>=P~`x6YWgtig(3@xW_K~gU%f%Q!wiaOt}>yIDc35f6bK7NeRu>4)14rK zDzL;6zkt7Wyrt5?q{NH2m8YSlfe%>~izuW9hU5T^>3cszhoblY^^)cVKi z7qtFExigx`!9p*kev3z2ZRp%C3JYsvI-8scMJZ{6jtp}psmHff+gnt>xhwdmszy7K@7g!e3-p@KZ5p6X zC!e1K1SQGR_ccAzYF9;Ee;JlAjRG%1Y$T88I~|h2%P*2EdzI)a*xpFH*d%Zn3X}{= ztepks;7@tg*o%HjM@*6#%xs^j;c{-s{uVpj2X+VEY`ZKf8_Xosw#dc)niKAG13DF_xyoUPTomHP7 zBrfzD`=wN6j2;tQ6EGphl@I=PoSsPFCR{)F>hs9S&X3ID7mg2DF@I-EZ@R!iS&Gb7 zMDHa8uJou#>~gvORF-@>X2Sr-g`X*YWC>mavEt35bOfhKulDfB_{7U2vzY0f@#ZLw z$M(iV6j*t`zsCTr5PK*7s3Pe%o04ls*GpkL7hb;EbVzYzWY8RNh{gi#h}i4SxSzn? zNxU!IGjW#x*J@rOsWO7a5U(Z8{oL-hg|nG3B}F*iL*5(2L%U1KY(z&h;c;gG?Eb+T zo&S77N*FmG!Lk42P0Slf%CWJ}4MhkN?=btAZXG?E<$h2(*=;u^2v}=d&s5uQ2DvET z!DYm^caD5;_Ru^XJ93BgWRBa7K&(s06gn*}p!J=eCUbmX|6I-Et8K%OlI<-3|E({+dUN6~ z%rSY=7*xN+lc8A$%Pjd;|-S&vZBXAukv;?*eYIbE|}(yjkP&pttWjm zk`ZVZ`V>0*ke5drmx|fHOU3BS%u+zI5K!J4Sj1NwCcKlbhS|Zu�wd8!o4+J!D(l z<84Vaw)lO`*_eI70VQRq|4>TdZHi>DzWF_G959Q4yWzlt^#Pdw`!6=F!C4Q(S3{Bb zD^Z2d~3yY_J=OZf1u>D+7fzt0-p&Qci1LdHxy2rqm?~Zrx zJIW;Imyd_Ay1J^%N0}iM#1f(+IL1S4C}E^ur0l;huPSDA#48Oe$1#5cZcJz2tJ^%8 zd>K%I<@1~M8XN2jH$?2$S6c8Yjk9&xL@nT-NCXh3W{iIgijfh03cQ*2KQ&@&il*ta zURH3KW~;F#Q+1uz>W`*A3XjL?fVTdu;m-JRh5LbAL?Zx)nnKRrg{P+27wcuBWS?7u zku?eLj@KJE1odLdxpDPAumXdTPnf~)GNeySM**V3IYc4F-f-cEe4s*%$P+IxVG0zG zIF(qXN`sbPL)oU$opwCVV=hTLlr_!bSQu!O{zhsX)tWvTMa62s2!l~n+Fk8NLQ>f( zKJ-Oiwta`x<;tI@7ebldj(Z6TxMSKLrbKmqO%4U}?%!GqlV5`bDo>?MV{@S(@&GLu z$j`6BgbSP~T78k73Go3V<`HG=bl~$3A_{}Jm5IW0z8Wi%uVAuP%cCG=c8M=MKSD2a z)%{4n=ejC<)%ziuVsN_iMj}YEQcGc|7cbUjRA@IB-hlaZc4qy~QV6(75b6P%ge}9O z{gCCnzL6xdSvaFA4}gDuMm{k1f&>()fIqm~Gj1wrmOXp9#4+*5`cE};3`S*xQ*L)} zRY2+Ymi9tCeuk8VFMO?;iQPvd0fKHvCZV{-r@r#AYnEjrRmwzY4KmM&vK0ANa zL9YUgbw4mIt)V$l=;`s*-YE^5jf%-ZemM5DU;{jAk{0aWBf`bc?-BWklY?JTU)9K0 zE2eTt@=t-_&G` z01$9~I8PCv$|Rd9ZHb2jY|gxM5!KbhCBh0hsc4ArJY0farfjwS>;V^|$D%A^UW2j^5{Fu@Ms#wFwyp^XhCKmT17IU2YJ`gnhh zg@Y3`ZYfxh`7)$iOim8*{^~FneDnplAqsA@LS9U{F(84tu&^-KU{40FIU>v!mz7DX zt7BqfQbdqV?{Nl7gmyB`Eky$Ln=>XPFwIVq42#34320)nas}t#G^-=YV(4j!dhyw| za&odG!n=}Fzp!!Pv@QE~DGL{JFr#O0T=MKbjQUlLW-4}-H=w&O*UedsNLs}7wV5IN z-+3u$YvZe_fg(sAO#$8qc>EgqTjoR#Ff-?;3EPkPimZt(og6P85Cd`B!agoFsu6bV zr@OX#o+!F;+nRsq^015gwe!}QKi)Q_NSf&Pfv~ug37`=&i7-CdE0sA;RduV%N$A0i zkXxvI;p9RSoV?;j$Viyhh*cKv0~7e9i#mmM+kg zLYPHT|JNr$*n4jf{rmg@i2WFH>rYwu*3@*o>f>U)BjnR1z6*}!(c#Cb>6Q*}~d@(y6@TW*1O;g2jX4tNTV?UY7_J;B8RepYf|u6xzeQ%r}X%3*FK z`?gQhK&(nTQ|X-VTL+P2?bQsm4Fl5h11_}gAdyhTdSxdAtT za>{B^u-=TuMc16KW>KB%AS@SjI9hdic)#GGCZTQ&IC5BSofB*qT3dYZy*l>v=<62b zd?fc3xNlhUn|Y+=@r3(fRAs?kkcj&QPj3jhMc2~#(xt|luvK$v=;_BA>;~3}=6hfn zihc-HF8WVs%bXzN=r_1Iz}PwsK`)#RRw9`K<9J^`lVIiB{67#ZAgxPJuM z0Ne$?3QdcQrou-7N-JWJIQTq`QR4UMy2zJZWB9B^A4?WGi73NZHfp?5u(JU4ITf1+oNanLLrNS%8^?!-iQ zeVxPVKrtG~{hzOyQ?2(xyE4z#knb|8V zWbe)Ian|*{@89pZkNdd)x{u@fT%W6WpYQWLU(eU`^&F4q<8hK^Q(X6T8T!S3jcd2W zF#vU=jWYI~^zl>tBmLgC<0)#Z4>C$AJQsr#W9tkv`_3uDe24B!+fL9~2-Dm2&T%7{mmw z=6FXoF{?q4T>yU&Ep#x>*kCT&SLi#~Cu%^)^P@ppSklmN45jJ6lctz2?@@GgT;n1; zci^(a?rUH9jSgoHH<;ut9QonAN-^GPKiS;Wv#2NNO}4S6X1tTB!*=H~P0+-`=V1#P z|A0$zM4M3AU%FGa)kF{heoYm|M-{R=8?WWqSOZu6!LhB+}tykQXu-n zFsos*_6pDXlAFayxh2khbS3xf_f=;+OvoOf z`p%i-vAi-@aLzksfYL zvoXJToI;d*)9dUcl2CKux7E8mNhOV0Y&x2?` z_MaXP4wF(&-<}_JzxG6ji4{xVgFSUw85&2p4qqvebOoH7e`AFTRT3}#?{X5(_MSe2 z8(wpq8((;I_Wmf`*n00~;CxBxC9N+!R>S~5Y)=={o*m8`Yg=#5S3e-AqA29$ZGQ*p zVhH6qyV;up^@FqX?3hu;kPyh5?iLf29lW5vDnA*YcxYT-fpU9cca7S{vf_DzmM;XH z1B@TXl3b0=EVA$REZpUw-b2mN4>U-cgeGgrvj}~?uSdVmTjiCmkXTx4N=TOSJ013H zbtsA4+1r92b)HhJU#u-|iIrkrr{ld(d-nBxl~ea7*h{Cn?5~EN^>S4)E)5B!216UZ-=>)F>-bxO+$i^;xsy|& z=b;Zjg7ufjCMLbZzT{*#PX+oW%G|1rp8DmsT|#Yk?+x1%w7wggD<^wcwdVskFM5Y| z8pVOU3P2pU%*^wc<~bAm<3gyE)mU38r}NTHf~$deg&qfv`Q9I)n1GNo^mktaH1}2is5EYRbQc# zyBl6GD(!x%QB`$Sm&GZayKl?SU{U$Kq+o+T~KABNx-J#d`^ATCq zu1AjaqeCSbPJ$n{GV|QE6D;g#^AU(3zkP1YNrM#2k0$eM77vBK_yZX)uo823V6~{B zcA?4|9VPm1-ni)LXKOy~M zzEVtlx!PKMY09@^wO#Exa5P~7HAjb-QFNv9QTu&~y3oSY=VyGFZSaQGI>oA1+tv=g zDS?Y;#x8~B{kqOUXGfbF+j1$Wo46*ms`#?n`sE?#@rapV>R&l4BPD|C56LGAW;~{D z&e)PzTKwXDqvv_EbePV@`Cb;Uh-}gRC!S&VRp|bhO5Hc*rk`iJT-Ys*Nq@}t5llF- zu`I9*Fe(5+8rssLM0HxZL~xbOLu7|%agjJc|a&- zPpGQ&4F@V}s=Fn5qr|-Aj+bPa>Ge~;gVR%1AAL98aJ3g&?q}rP#DA^Fd#U2;#&lUs zhb49D9;%R+bb`V)`+bVZ_3Z2&)|`ZaA_NGJ@*1ct4p2}Oe|6cJJ^lR@0jL8DhMzoH z8aSWU*B1&K+SbiQXwgU(It$e60ExZTa`y+ifs>GP?dw z~#&O4>BntABqOSo8!tn!Pf(9AT76+IHv8>GcyT#2@qZ+CTgTAiQl*p6%ny1 zDMQ!uHAe9?xPxebp43ST0$CN9T|jTq(=P~=8@LFIk4#SM|2JxXA=}{+gRt_#asE7d z2VSRR$4{Oo?;HE(ds(xZr%NV8&b8yNK(8_y3ML}bIfaRyna+hG3z)*F&h_{dNuifj z;@3ZX-=>1pAnEl1fsM=L#<>}Tmv~vfMo2_RGyjL*k_vl4)q8*Cl97rno$z5tfxwA7 zp2=2p+b%xZ^6X0Md8^Ge%MWN?aX$_g>y~)J%y?IAW@&P~Te?iecSDU!z(avTJuPmdyPR}%F;HWeEzzWAHhgZ`X<|8p&kZ_in#np8Xrhegrz zV?AxVp0|Uk76YvWTaeJ75jSep`+zJ?TL}L}pgoRa8+sMidT#Yb4Lhsf<=yf5Y0udX z7u@9v<|1Q|8r2T;Gnx6p)L}&jpoTRWxI{Fh1Ps~ zUd$a9T_|@V`46NQs9G1m5GjQFj)ryp|kPwoTNff)Ye8(2r9`B_s*`JcEBl&=O@E@OQD$&(*HEI5>}w$Ie`b>+O0 zxw+@@755?4f{=OV!bcDLL@#%C5)(tcFZZBulA(8CtVUS6Io4M{DIq%gy}yrkrM`2| z_K>pug_PISwR**ge23`eGaAz+(QA*I$Q^hn+5e|U(Hvz>7j9TW%uevJ7VMIYUV9ff zWV5y}Dpr~m!>L?kbe@0adb={cKUlwxqM0q}wmCi@bk=v+EFK-#6FlDYimH))bInuF z--^sF(*9)C(jZ-`NI1=GJM$pQFg9%Q+`}fWQcii<1k<0jJ;Pr#KAvH++hrB)Si7%% zNZJlqz22Ktf5u~H#;yUBb-ZW3{i7rw9YAB>(s?kj!8G3XP)1O-BbYbZE|0)>^u8pq z6c|!YrztM4GZ0!i=A>fW*E}+^R2r)gc5dH0ou=kpc8+tkRRue zKDoPTst;3Ky(u3UWz>XCviJ9T`7O@5q44%)_6~dEyjgrS>-DBb`?L*5(_<;htBxI= z8aV@NnK~G2cX4rK#PmbANVCArVW*9)#HIfXy`ZS$+Rb-3LLZsOEs5Cj(jiR^fD)sZ z*5bP~N4_F~Ip5fa@A5S-{)R^Bv}4{1vlBg8V~>ceJbyCRHjAxs&l@%c)_IocU71_) z&EC13L=eFiMDY}ckAvrGe5E^X%1^C!+1I%iTUr*4&-&%rhQ(8->Mlo?KREV;k#BqH z{Vy0RvvQ@7_?YnQ{=-u_@ocko+!XHnS)cq~=u^t@Y-3n`@zk9mcFI_~78HoGE9jn1 z|CAN%7JH@^UV&zi29JM7mv{(yQN)T?yjx>DYelnnXXi7alQ(q8JqZwj+9hSSs|$(p zJZpb8)WxM3`?3n({px(zpErAqYY){gr8&ki49peg@PA75xKA&zLsN>h6DU+Q#cuQZ zCZ1XzoT$yOlX2H@Gf%l5ilgOGGId)%M|kzcw~Yn18<)2%hp#=UcFq$Yl9tw3Y%9GW zbe9;{=}k%We1GDY(l;NJT(>3(|C--I5iJ|_k;X!4A3v3kltygyAxB9|rX7Fwd$SM*igW^7XIRqZOa+|i8aP4-i}8AKP+S#>Lnh!b8j|!XpRZP zO$CWJ-X3v$IsWPtrcwHm;v>+a)!%3?JQMRyJ^bv>Ypjfl9-KG!ZPQyD>Ha7cQXV!^ zP<-0W&S-Vw(g!7lnnJOvGSj+)i=;1qb>+zSHtWjF+L0Gn3OxSvKy-cK&id=tZ%K+& zlTGh}912Y-o>g|0yzqbGI1XUm0WD)7jG}ZAuhCF8d-DXT{lcaoexPa4ljsw+HrxZY{qyq4#!)6%JvRO~T-+D1_A3vSV=Iy(G8AayUp5lnBI zh7?UDWm0kE`wJTzMn$%}Cti|qT*u*#Pq}TY0R)@Ylx0p>$bv;Zf z7kJ7n)dOgTmKr|yD&aN1=Vc7617*T2bCK>e)4qiPX5IO3+{tzh=O5_xLi6GD>4-!g z7Oh`GbSNOmFKo4UU!xmAVYHMu5f{;*P37_FJ{FB3$pZOi|C^r17&fiM_l_l_4|^siR@BLL zkrIHBP~uw}*>*24Ug$N>IneygJ6&jr2~1KXmX>|O79Q2ST@f>d#_Cm7CFk|Yb2M#Y z5(qbKXCKjLyvnkEI+)M8GI#Ys{P*HP?Wy)4l8%nd;QqN&T)JBKgD>9c4ZNUFDcP{> zeS|*^+1UTaSaha5EW2y7N&q`-Y-M1<6)S__+oKDePO{39=pyKK&^Z z^Vj^2QzPEICWbBIBDW@KlWBe)6-}YLc}XZWP}OtneHz9lbP$X+`aB46V&(cx1OJQ~ zl+*a{DH}vRxgbYm$4D>dZ_ULXW@3-36Gv}psi5@417k779ta0$1rNy@eF$*?T3tbU zwawRhW8GwxdtgAFp7;T2PJDd%NsSj~+wY>0TKx9yIXQvEDfiVkb2`(SGudYbFMnh4 zE{zdhy`IAA$y9p)A!08LK>}bUm$j!Y;ny_&3@qS}kX5PjDY2 zDwz%Glae)MW>#RbCc|Ae7%6vi z^YeRMm5;s_C)0B8S?T?OCax@QF70Df#c?3&VJ5NIh7j^XD&EI6!)pc)J47X|_X#z+ zck1cdteWU(S058TuO@XB=(`G9liMB2Rhi!miAUP**dl55v~-`Q!f*E7y}JuOp?916 z^t881OLt~C=#lI~!=A_QWqBu|fnqQks($iPyY_jzjP115nK)~~N0)sCPD`@(&zmrmaPtyB;Sg-=KA(EVkl$X#O14h^F=xNms)wxT|6^Lxe%^qX zmK>P_WakV`$>7X^y*qAEZofCNtp^0G;50r%b2YBPt($Dy7O?k$n+80e4t&L!f&rr# z&puQhqS)dSJ;Agxu@J3)$S(hg9j7;lyH=ex9@UWUaKVSRJpS@%0Xdm!3Z}WpRt!M> zo0z!O?&pdS9gR%a;0#Fmb{gn*caIjXy^Y=1;uKv+vyNY~{y8?gWJmhDU`$T_T%Mp< zTNapszX9KPUOxM4N#XtgmmQ5ZV=aB#hycDLK)l7su$%<;7!>*4ga z7a83F*-4_>x;!gQ9E-mmNEtc#lVkMa1vkvBnzS zxj-8!qbH(rI5~7BOSijVtiAfUqx-S;?1?pLfw1EXw%seUE2MMAx~F?rI-0A+_`_C= zb-P!>o2vb~b+g(X?~aX5*?=Cg>}PdX*XL3$gSf=Ro{7_IJ~J|GlEh7%X@*G_yo1km)n)s33o?quMUtNyOf%$!dNqaNYUr^om1_ z)%yFsY-fRm6EJR{O1FPOPoB+eWqXQSPQs4P9}z+HfR|Xc%GB=j=&# zF|Mn-Q#tVZBFBRJ8w%H3Z#SFE9Xdsq)kpdnks?W2WZc43+TRzw@uS?aJuB2R%He^4 zyIqa6-&pHH>zeOrpZaRYcIOx?pNnDVB>OVfQc^$I(aZiAzWYx3%|!EsQ!y{JJfb@q zzFyMZ_pIsaw8l&=kD^@J)aWAINUC1E!y#>U17@ z1B~qMr1ZzmojVwipnG%$1UqaePZkZB9YA}sR3kDq)eo$MGPiC$ZXaV26r{JK2%j5i zmX(!tNlYv$D+9eP8#_Be=PakFNI^$;Q*J_%4D6Q#A9z!f0%|g#8UnR@8t#Sk$VG9K z7r+_cV|Wk5i{KOjtG(Zg7mw=nh@nide4adc5_Q0m($dW`ii%}E^aAlIDXcVmF|Ov? zJ%Q(Q1sN$RJ^xiBBcsKL_5w%yMdP+HD{!uYe3g=#+I;fYCG1}?Xx8SEWf`1bJqfakmhvv7iX5rDGFX zU(^I2rv8Hm8s#q-F?tpv#R_$VnDlhY+)fnTDfvT%N_AMF>6FlPoXU?3B@a%C%+pyDDUH%m5%i;I7mSa5G<{uF9G+Cuc5VHa1? zSW9cwd+^|o{$Ra}!TV0}YagE^Jyk1xzDPUTSrGU%XUFqP$IV4+A}gCqS9W%@rZW4* zUXg^h1*zU3o(Tk^xFL{G|1!D2cH%^2eEcr-v~_f9-tnR8{d6B)d}3l0NGdTZ2%=n! z@{P7<`F{GujkjLB@-Z>d6DDRdRLh991LS4fmdBH(U8qxHm2DV|y7GeXs~Fcs<);Wd z@iR+BWl^#E=pQ}Z_)VC|hY#O|EuOaX?ceng<>SePCY?0k)4*H0Z$c%v^N4EM_wP6G zgKMZDoe$m@{NpLO8(~{7U%jfUufGi9TN`>_gBB6Y>Y=iXHYdn~2L=W}|9eYb-c2=y zEXy3~9$*QLPDmgLl@hyq_b?jxp#s(@S{GrK5;Xgy;dhtYvXhd$z!Uvm;S}h$DXFL+ zZy*|sMos(UD7$$U4Buj{o2Uk7?J^jA?u{_uY zXsO_4F`0$L1c^f`o)*6Av$M8mexa z$KEGCSPYT5cTW-o$SeW+A_C1#O&TD%@$q4y+oim zE^f2M2rEH;-8_=oJ_dr*8#iyBU}uMU=!gjL;_BL!R8TVkNef8RD0Xcllp|oXVDk14 z3ZfO7nrO@L!VK!(3&Q7kcs^nS!VO?HJ7&JxZ`9-hXrcQaP3+Gf{9G_2i;bp`bnp=v zo53?G25$oC4=fqPEZ~Ove|>e?+6IT8ieXF)l7cn#7#|-kBp09wfhrfa4|u+>fem?6 zV`Jkez7>?Xa`Mo+IFQmz=W`J2CT1`Ql3Lg$26i#NDgz&ZErL35d<)z8pN9A69+S}y zl&fo{g~#{4i{UY1iZvBZc&$D<;t~5kHblQT?trfMmAIjzvAdXsHPS_J0keWBmJZIi z6;gnZ5)utEckW!rGJi%3gqJWY(3Ux@{~Wd;$kFidgv%bb=+R7_SpC`7b^~{-kn{?8 zurdZ`30w^RDkd(@-TCThSAH;zrf*^88OX*B&db7MVlI_12z*(4>#)TqGyqT%MdX7e ztNX3rmmn#vZk%4xdZ0uh<*Py#UHGZjRNeb;(0dQD^M8x;DeD{qfw&;CIqK>z zDJW1~O;0Q9V-ylPg4c*a;KfxfoE~G29?dZQmTYrY`8oN+qx|AoZAA_%;>w=ZGD9^n zT$EQbQ@3uKh`)J9Mb$I1^LeW#3l|5M!$qjb$SN!E-Y;a-4qOX{aN(4T#;gKfzkZ5% zz$yvb0bQItGWhA6>gtT(xsHv9`1OWjsD{II2Qgp`?pm_&+7&CgGc;P&A$)r$9KYFv z<@aVF2lM&!R7s78NTVQn!vUQPHsD>jG>ceEf#*oY#t`*Bbj`wo$JEpmyC~J5l^ebs ze7x%D)o=QImjk*kAn7GUHF}Orfw$TWRJKq=Kq5eX! zM#(Nvh=P3mKD<0mmfzK(bV5B&i#_f?8ygQ$QxmcokPpHhd0NPh`DXRmp}JbSeritD z=TdRUp^1=I&A&&o;eyTNQFJIl$$eJ6$aV6%DWGmJ47z}d5`t?S>qqMD(wS>b=&yjd zt+-V;d>?$tRBtI@2<%jQKnquT|4eON-PgD09u~2FFSvh8POb#q_f21*=}?UJn~eM` zG*_T9bOUz;rr@40`FMGFDvZnTr901_=Ytbb3*e;?$N5(H9ddelNhksls~{Bad*HaJ zYkTX++$W3S`pJ~OswxRM9R|nQ1K9sP0kYlv<_T;w6CL5nUQ6RG_)z(mtjG^b*J-0tx96yhRBX3b-(Wvk~j~l#HPo zf-#W`AoG5AGm=-2m5B*ou!u=^v_l4YLG=E*pn+rZ$r54| zPn)N}jvR**4BYbpZ{K?2T|)tec-M<-PE;Is#la&^&}S3hNDzf1rvg#JlwJP&yV?G* zk4-2Q99Zw^4h$LW@83^BaRtT#esbcc(7-^DS`9eC=AogX`!GG65SW>-fR4nb7DxoR zxw(NjA5j9;x@m#=pwxNeE(t_3ga`#73z6B`K^YI5!NPnA+C|Xga_MsXtXFqz6G34Q zHftOJL`MLG^iaG2bHlEUM;s?k^bosIkB8fF*#d9Jr@hmhTMe`VP(1*_bD$7-1LD%s zD1v{8T3MaPUlN9gy^3?|&Sw)w$QZbF1uQ{T3$ZiY9JlRs|K`*!_u(oKDIvBPX-*l$ zS%Ym0Lg>Z21S0?@hRTNOqARMZm`K zgKP%rBTzyTKC4cz|9x7uMt`|q?TSl&diu-zXRuaI7lOJ$j{tH?CO*C($fsccNFm@P zEHN_jaYsjoi2%R~0LkqY!DWlx$;u4QyZ5~YF$D_?%j(t~WEKuiZ-pHBLCs%0Y%vY0 zwlv+EBT-RN7dBRBhG!tmfWG3P^#UBLw^URrAp*b$YIb7H=JyI<2I+1xhw}G)c10J7 zG_70r4KD3j9%@-Ta$|6b=?{6rgCq^{)%^zibnt(_zC7Lvq6rfz#EvbvXHduYg3prF znHB4K2f9d+G9IuS69Iw4q(p<#-roMgheLaJ5I_sq%qB?E?O|hMd(_2HLYxB3%mmeG zZ(ko6?IG7gnbq#Cnxeb9IMu^&Rt0P*@$nvDMt||*g$NDc$AQli zHUq9`*!GEuldQE=U#x1s{T89`)9fq>-o!JnU@i08biV{?eMqZBbdsNl<|FK7u$D{d z=&*p%?on6Q7}&sf{$q4OqJtr_$)$!dE6udlA`uDisBTII= zkf4@>oN4c0>FM{i9@JyaJn;-V0q0*DALTtB`1$h%#6pXZ1VaweODGqCEWWngHXRK5 z2vZ4FI8fJvA%?(BK%5T!Qvd)oM?$q9V5bmC7@|ZUG3n+8S`JXYFD@>cTm}s;v76y6 zd%|QlB{dQf0YpyLfBhCEAXh~LD`f1_l9GfH$?WW`$rXh4@QPyCeHSLNQlfAn6_u68 z!IaZ666oi5fTMeXc2=baxd}l7VPRnheHRQdhKEZd;{WX4VU%?ic`UxCxNv!s!FiR%w3Pxj z)~9#jJh=C+DUOOtK%!;TD2_tN2oHz%TC~!nk?JZO?M%kEzpEMtkBI=TyR=45YqxOI$vAgz zPGfHE{@+WTb2V?A zowajZOfxq#D;l9;U|>)UTH8oe?6ANs6C~E z^OynEE~CMMnxXD$i10P1rewW!8tYiNH=4p}#LUFxnpJ3QPrOP+vCZ%=p=DJ?1&MR`ic|b- zxBD#VsiF2P!k{fj8M2l$-LXC@fdvz35Wq2(?ia#h5%4z4QD>(cSmpQ~B(7VOl65`> zu=fWmi#G{ulBxYTyd_45<(G%)FB02YyFrB^XQ@#Yv49H7uJvd~#p&g_7r_#-ll-Ev zA5Eif-TsAdAwrI2cNG*A7FQA)hKsau=S5tPyH;b3?=2U2Ij)=$?mhf~JE&+R@XZ^d zf#Ra-u!Dglp{~CzM9EUHIAfFlKGILU=b|)DZ`-wBW-AUExAtUb19?j3Ro)eT*CB5b zm>eY~<>Ay&_ur)22H%u=4J+7<+po*&9MKP&($Z$>DGw2^SUh#7BYpeL5O0WMJKS&G!pKKh_%hy5S?&7DN6xT0U*59bQ(Sw>a4uoI0K#>GC=h<`6!hiHR8H>s z^=DzFSnx$%exT&=*nyFi0ct2Ql3JjWf?xnaBQmgjRcUCd`Eq?yV(+1Fga{%|Q%nXw zdq%QT5xEqM!x?U`Z^{NSOpC#JaPQj~ylftUdwsfyv=Yxy;vJj(tM(q_}2hTC2I@4MmTCA4WoZc%uv=P7+AOs@WuwzxS zu$#fDaKo4OFQND)us5f*LUT#1Q%YW*fP+{fDP_yo97fCvsBMeAyV7h%T3G|d-j%jH zBlLv>CKr%zmCUOq#x2+&qCqtXW>_?V^Tt7_`;istwk@AKjZm7P@?t(bw9K(m-yY@| z=ztI3q8zna;mSy>87ywGr_GawpeV|*4%ZcKQG(jTD2q;B-nXTiNFL+a;H{M>lvWG+ z_&7K?I;!z55%VC5Lq!7bl8q8sP(#0-#`QUJ&m6DC^_YpX=hfML#>3{HBe}cJ{V_-9 zoUehMtD406b_#OV7nI?Yq*3?=XwH4T|K&l32G1dqIiKjU_Vx9#j?PX=)VX!i)HCLv zTkqMnpMH?*M$~BKQGS{=-!eN*R1XVGq=Jl$FcYGKAY!kE^f_~ZQvn9z1pf16Q35=%zC|A z({yVS*zsxN!p3{6aj$2OMK-BS$Q-h>zwZ098pRmII75PW7L04kt$aG9ayxxsEtI?6 zHgwcF{+P#Djju@6220EdSJRX3T<+>Q;x}e6ME|EbFamXGuTpX!m9Q?Z2*38$#>VYy zy;S)tc@uZ*3PU^dO-IC2L8VCMp@O1aYb&h~*_`oU!Ayg0Mp~1APuqZ5cAAS~iuvd4?i6S7H$|dTF#usX}FS)W4z0%q>Gws=h1U|aTZbmmV zE!&!HU;dyCih9)~%S1H->RxVb_Scd(_9IB9H`DFn$gQtB88y918bLrk?$x3(I`&eh+g)oZ= zP(0>HTj^)m534=frW~6}(dz|y7>GekPDDAMWij1L=us${v+oGqtX3-)o$Dti>_l1k zXE4QlsXNX|5m zP5h2-^R>yHHG~d`>(Bk@d|@)wh{ks7*|R!f>3k7OAq_2#%tpsTgq?+!);N{`j`xjX zQ`)C+*lA}{qM~X1%6=-B?Z)MmKC$$IvLcqR%_i@yIjA!qv! z&wj59(3O|Hbq>X>{3v1leA>&q4kVOz7TLd=nChAXPz=ERQd2k2y#GcXdU7>XnMm<7evzhabO&+>Ys8NY1u*|3K%+ zVkP%_#yjrw`lE64{15f;%FlJ4le1EMvV)I{ufL>!wu<@(D{nIBA}cV~Kw9Y1JFkIh zcpx5W0(r*Y^foh=F=8P3*11<9Zw{jaTob?eF7YI$*FaEkW>|8sH{FfI(-0^q3Bq(N zA%Bw8v!6|WdOw#Q57M~~5Yyf=ibC*D z+4m}T?l`<2)+zwy4{*6uOGX>>NB0KC>g&zMxOql}a0!zNRRVA&ycQ3l=IGiJIg=v| z(7l+>y5}67mlgt`Ue!(N!jj6sysEK{@cb2@4^hf=7ghzBF{-eO>1FK!{oC966UD{e zQeL~1cPuJCKK|&R3=_knRG#-EcWsVg%U?KjkM@iTkRmr;}j7gw3fKKZtCuw3`IDC^K$y-{^t>zf~K$PeL4in{!z z$N`mAi*A+Y&uK*xIbPL}Iab%bkIwX3Gfa>>sh!(#s$wQNXWj4kb^lpbC0Dm|)J-xv zCpDE^)c1POi@djH=ThaG^%^b-skyACi+Q2hyd{uTM_H%1)r2NCW#0>mL zUHOFbOih1y&+Q)O_EP5S{^M`^qqL+3!&#~YJ!69gO@s>%CXeq{?mXDRF&3--hwH6< zsTuzU5A~jn1xv-@E=CGHeCEwvM-Jw1bMSw2P1lbu_&ww%6RPu~||DH2^T zqro28ao0%Y`B#I8mc_!9q0Wj+)3`sQ*W`tAG9PX~mJh#waK|f8jk%3wS%C*yt2xc( z+m!ud1Xwywhlm7<_~uYSI!bovxtAFo@s+kqi+ zagFQzfvXRWZ|OVfuXnJB^L6^K3*;7o8`Omz+W*%jzjbzn&PR?_)D7-;J|C7pY+rGz z@^~S8_*UAT9>NsTqStqX7V-D={8G>T%yhl}UZ2rK^;Sl*7sfwRtNd%a-bnRlIK5tN zsxkYhs`&c$<+nmNzopo}l!#!yo}(#u@}Qwhd;N*`6tv8=Rnj3fy6!@KX;;Fh9D5gH zhHek-Z+eVxo_O|knfdziicErBI*I<3vX4G@Z{2xgEHQZ2k}EvHYpKh|tupF4X(T;< z{}L|#;_tYv?-fn+Gc)JQ5PaCc8_Vh2Xt&o{*SWnErZmnk%qo1u z$L`HK*V?l=@8Y3fH?m#!WH|h(&yrF*;^p*w)8MA|(ir?;bk!tuH>#OiL@SsCb62C*;{wo!2Z0@=QW|4I-<_q+cN1vI|K_+2Dro8r_j2 zK51!Z(Q8t$Okr{Tf{a@7-<|kd2y{+uYpbd8by{1V!j#JlCak!0s%XH5YD);mS5>7a zlZv1cPfCj7%TM&m;^j`FD)s)40(bBa^dZb=hi+-}JFR{mJ$w4!<8>mf;zGgd0}2I_av=*dVXeWaZH&C2A(_YJj&QBptPK$Cr;%vEqO!HRmb9^~1ZfpC8yQcXx*23W)q4$3<>Bf1DN*J6 zyW2g#emx2f4(86z&5c&r7;ujI_v&ZKtjf|VQEi;YYcfT@5=#fp$fP&W`GJzl?4QOx zk4*c%L3+bDAUJrBNQ!Dc3LJy(3{Lk^mPa+}C2G~bMn-O;?`z=n{e zGg!3CmoNVxS6z(DyLs~_n&FXod50OdDJN6medXonr>3QCJ<7~{Z7I)sU@#N>uE>Yc zdy2Z9fH#S~Gn3ON)}B)*ZL21@8ai9^P9dpIkLQz$@8A3NR9SktKrbCJ5leSH0-r@0 z3B4kD+NCSYgSc=pvH)xAg@-;v$A}fI&%eS_dw!kmanr~g%w6*gjfs_&f^x|_XBKq8 zQ4eDYAVE3k;>C*wN#ewcTlo5sqA zijwPoRBXz}5C8jA?yGIfVLy%_ITejQ;4K=RE~Gf``>H|MxL%J(FGKk-J-w-g1^VIN zg(m(QR;R!H`Ntm&{N77^h_QFM@$PomRV8*R;()=Cwv_n!9;6ekih|lc{m;eK&MxTOoUka)mnF(D zVtNME=l-OW0u4-QeNV?oHWn0OPQRhAb!FFAypN^d!Ny7`%A<$!yCd5?gx}JWP+gyE zRwHf)a!1Kc|2-_{3q0)o%w43a-}?eA<0Te#Em zvw_M$LaKCme$30)*Y`Gu>wP$CLN*gAxX@BRh0x0AS><%SuQkFORqJa+=Trmb-}uG{ zIzQf)5>6jYiA9V+V(iCo=w|{&xWhCFYPRZM9&JXBth7e@*LtomYzY$S-r;6i}oE&IPLy(fRoAa#6@Q_qb9t!Nwy+Qaql@@DIA5^7re3dSD7 z9IxlFR%}~7U;B#o1lOC-%fUli*vRVq3nR!PBfn$ONH|HU~Ean8#O@` zv|Dp5Pl1>EopE=;>{xqXpsL+7bpI;WtNuL_uQ*F>*5fQML8lL-7xXrq&D#>>Z3NmW zevRwNFIeHs<$3VaLpkFCp^5zAZuadl>UkI`x-ZbEA{t^%AMO!d076j}!n!@7)?KNI z&U5>%{$0I(BLGYI>$jemFUwbpWQIpB3f}pml6_~SMbo`Y+Dp4Kmkr!$5Q|Ip*sV#wm1gI%l+!nlFmPBeaOwOcjS-9GN(g@)Vtfa z;-u&c*5*|Z(R{%HQ66+w_x&G>q1wmL-a3f#Dj{|T>1skUjL^Kr3((^^1WjE+Q;q-@ z!9cY!S&o~i%`#Rae!!zo?tSQ-2KZv+g);woUw;cunP<)>XEgwG!l;x^Te`l*cn58= zLek412^Lfcfwe*;rHSAwW3>BmDZvg-f&L24h+7GADD~r4Y9D0a$cJBrb;1NMmFgWv zK?`R#{_2qT4pvrHqM2l3B7(3@oe`sBX|2x@Y;<*XZ9>74t_F8JSQNDB@S-#2|CzHp)|P0`16jifC$eu~hjYrP*rzCkBoZQqDHEf;jpNlPBH}1f z`ktFB?e12)9+oily$5Hhjm-r=KR;JDH~b;#F#E|@yQQV2VH_}}lJ0JXWCbsUqzvWi z%E}Srz>7EkEet95)_-Q{9}|3(NId@C2MPD%*MI-(|DV5JcRz5srTV`wHaZRuB2~0k zRFLN8SB+1n^T)88;-tfWBm6;AgNUo7LM~ha)Q9*>=-^>i&t!h|4j@0T`g3puBF;2l zp4CUB#%DE3Bj$6s2fbg&>S4_;IIUhF{sIP#R_lQ79}9jC)x2`F(fJS71V(N z$0$SV*5*p~V}r*9^sd`)lNIGOtJfV})V}bZX6xYj^^{Rk9PRH_^AC}ckr77;0oDoS zx3s~8xKg^_0mSJn$BsR@=0ydnds1yajF+03absMyDugcpp6(-htIzq>lFmUH1g`@$fs1FcF%GW2`iQ~keE&bW{Jyc2=9D{>{X@L2lXM7ng_kw}~ zLX7#&zgVrl`;|n5vHPdX;&(}Pd7tcm^yn1p{f)zH6{*RNDySJeX#ViyVM@Lz_8l5F zvm?#EIFAlcQeK7r_Hc6wgO!yP8XL-Yq;_3GuuKcx1UQ%|;A!9?!{9NwFlhXzIYl~7 zDpare13f$?f&6=fn5F@C0vsyIzJ2cie&W%qeFWeR)+hFTHbIqzz&i(y$2BBmqwo%V*-=1=PWrFnv&gCu9N%!(|7#C1g~X;Mtrq&# zzdioL+k`FliY6;+PwZIv1O%dQuTAL2#l@Y4AOPySvau5Cc#Esh9)?rj45wHzTNp^! z_~6_j?{5zTPrgdLp=R>xhW6YvWFnl(Kitcufe<(Lv@R4sE-v*4AR0|7QF?d#&KzIL z>KZ$hmxG_ozdd=5xCfsj`|9}f|H|hHeUqt-WXE=Smbgt=j$=?~?oTU(yn}shOmsPDyi8n!G8vS0OuvW+=E^En>UvC z?>PNrv+(7Vlt`&SPuJy?KD2O^lXZZcBqB8eQ>0Zvju=B*WiHm%&aQWKR07(XV!lU! zr}IN9O( zwAHn>1jUpv!^Ys>Jw0DBpAFz(-hA6bLU0}DN8KmwoB#K9%?=N>pCg^UheR6+t3{=y zZirD%aFXlj==i^W4M~j?NLq=7)WygYoZ|2_HY!W+whSSS?d<9{C<=9~*;s)hF=oGh5li)EXef zABY4iP=vYQuy{r^x-OId?%XXuHRmvFgp7OaXemP5W-G!rW>5?O7FUbCiN-d9yBCI*n?gcEY5B~U zkbKF+0#dS{mR4V7F7p&`V1Q~!_2W}JmCV)mB~_R`156s_yFKw$RIV?tAvnNGw)5>PETI-6ynN0LT3=JqSAYkm^S; zukPfpuk0$Bd%Z=mj?pnOSFgube*ezK$yo;K4(=de$QhZLx96++1_lO>`DQ1Pn^P>U zmn_^pnNw1`d?l#CkW^S``gq@4-{NEJRTudK9EVp7sB=%l{7{d;Ab>1zyGVbrqQoBn ztHvLQ7CIG`#u=KI164X_5z`PmyLp7bvhea!W6Ed$6X47V(W8nAG?DLg_9Fte;18Wn zL*fby?TuazSDB_eZ~*|!5Lc$8yw4r?BgRHYckS90fo$nJB|W{jm)9N~)~MWUnuIZVh*jKBD6=GquRKlsftI zvczjT=_AH+xKGFB*`8O5ssze4-JsPU9*;P*ckHDh;wrt4&s>0~sN-;d`D!hNc@-W; zBIJS}P{cFH!W)pJ9y)JOj?hOO1xn$K)w4L#qA$7ZxCjgwg49zLBd{lePDroy6Z1a0 z^ZkBeN7SKvLqz_C4Y-I!{`^7>J?)WB;Y5TFEC2bliT=aOemC*e#k@<*}?ymI*q&vFFwX~h*&ufVpreuzL(S>4SIyi)my+Q&Bid3Vc0^YyI$e-_5AJ+`zt~>a?Uy1ZId8^t%fL@qc}3HA2w{ zr7cF9yb)xYqVD}}!WMN)9rzd>9QnG#_gFu{8K()!wjv*OC-cm&Hoc)Mc|o?86*`3x%9CIjCO zB9)qQ7Z&^H78O=e3`R*Ab0~i_Pa>~X;9TY7=kFaK_X!N7R#jESQQkK+WK*WLEAI}+ zJjVsS!{wZnapf{vDiT2n&#U|;(jS~YLV6^#G|NpiIucIzb52fFonEaTPsl}E|7JDk zP>y>!wRyJYVK%)aaC|o}U3y$EtQ1}MYix{om8ir5dxk`VWbfX8Ick!{*X$fMQ&G_! zKoLfuThM{=u8Xs^SE|B=!!cl!)TE71^!D~5^gubYsJIyO3CLI0)E7MAcra#VInq?x z)pZ2vE}j8l`lsX@*g?DS1z*2zh8+AM9=%66;DO+WM{HLJ{OhHc8y+Z#F7r~1Q%QYR zxxi{J?f)`Z>)d#n$P1i0KI!z^a<@_R%@wY9SwC@icdvJ303>Io z?*{S*)N#{V_Z~dmjGO_53T;6=Hl-8!4eT0 z+r*NO(--rjbiKILU-O>7(DS+fUUUHsq+6>;`_isAzN=5Ok`5VNUz{1p`t#@Y5{^rg z^;KKDxj&_SzQR(09QStZ*>fej?x^d>8z)nwE~fCN+y<5twO1a)wwK`4X-HP#{nMBf z$MijuWF(_WTczf8ZmxiG@@sJBXse=vjjSg%ai_^}Ju^B8FD+FI9D-0CA<%GER+<0x zd}iZNu#CGHijz?C0fou4ieM#R|~7@OZ9CKN*6qtD9vJwdN=Q; z>M+41ro6(9bpd#qv}C&m04x!h9S5ElUyVixQE5YRm65z02Xrv`{{3G(51wj4xfwT4 zly#AFwkceTYv6+8@;S_h?!XpkzJKZL41M#4CW*@6Z|oG)evPHp**_`tRI-c8j74!A zt&;W8(_55lTT-|LqwTq4r)Y9SE%G(j!jmv#B#C66!04vrP49+h5=&|%KXmB9zxcmC z*V-5pO+mon2e3+GR4^3Db$g=o4d5bu1?Ool?_!Y-BNWw$AR8e-Rx&Z#f@oVJ_`ElE zw!gpsWi(81CZ!t6|95Yfnb(JGm%EkzkNQ<>;S!&LKK{}ed=ebmelC}*gA{QZhLjGX#Q=BoVOF=mq-P7ObXW`sgsiQtkr`rm>l;cF zuu6msReyZzmcC|j~R*UV5M-{(*{p-Jy77D+J zS>PA}8tlEd^V)jOQJ9H{h3-3rhz@xP%{jVQU3^-j|-whwJ6EO`O zw|o#fI^#;y?0?tMLg2rnQEQAd~Rsoe8tO8CZh?zKOD!RIIa zja7a>^UL2jhdCLMc5z|K4bXPH+eRO~$cEwaA3xr`e_xe3^a1fa%BF0EvY4hq;|2Au z7KK;iAGytG={Sx1Um2!MPObX}1zzw?Aus67pikjiSv26Fy!#@V?|-%S=3zar-@m_& zZP=!;Gla-c6f$QXGDj3m(qyJdr9@KOtdcpIl7!M*B}pk{XxK7_BpFL6G)hJ4JeK`E z=Q_VX&VRqpb$zbOrryK-y6<(b^?W|pbM5b)H6x(xbJo}B^33b?gG!#ed%9J(P%mER zzaV{u^Yyn^cUeE~)(58K&f-ajXE!_097^eN?m(2(n;JExi(qNXCR-c^3ujCp&NA+> z%Tu$L4pU)F;f=)=KR;J`kU>)HnpCWQM)kfb_ws--YrouBIz;-HnzT{9$4Em5g`WG~ z1&4Dy*AMOfwSFf8g|nA~G9ITN8uBwFXRJzP-`siek^&CwXnEPv2<|HCVzU;<4(L>r zZQ4`jUi$u3lTfo4)ytOa9GO`4_WPF^Icb{@v{+iZzZiwb8}(cEX?u_Am3{B_X~gP~ zFZXil`#jFN)H2(Mb3`LeYinz3QWoeG>y4GhVxqRY`vk|hUEivV95q)=yVde&jZS{h zh>}MmP7ko!Z8c@Xv)1TXaiSZS#5a&nWFlfn~|MY6jhE_N^^YZ`f_@r~}^HckTBbGBff1(uTh@ z0)j?QeUNP(I${4mE0VpY+zX!j{?xYW@l{=M3mh_`ZtLSMO&vNWEujxVp9n2=dV)*p z=(4uD(#a}|L)tO}o}lx1&ylo9K|Z?FRy*UN>ziYLmUO4_pMRD;ojTY{xf5Gdq-3XexX4yh~-yKtObu#61oNZPfN!pZ1lM337V(BMS$0rA>V^f3J z|B6cM`$T-)O(}*1KdvcGD*5#K+89`MTH$gdqB2)5r5A#zR7%2L@n~|JnC16XXsNR^ zo|`xh5RYuF!Rd}+OWGM*U%Hd^ah`quv^x))Pmc_Bh!LSpi0-x>p-nXdd|l)I8mf8h zFWd5An!Z_Y!`%*+%}QJ%n*p2heq-EJ(nq^FU8(mPrm7nL<%_|pxNu8L=kwOx1Dn1u z^eP%wu;O9PK(B?gc>BXlbE>OUp|x#bJjZ9i>AA6KvS9grq~@8+%503H&IckWCj7cxP!lWMRIk)`-d4?8~QqqZ)s8chjRGSogU+I zR)_U@a`FHR6DkXXrTc|XuFwB$8q%uSg#pv@i+4|Ub?EOgj>*rtdAYQu#sULMq6pN` zX^>5WMT-`Zz$!e5%T$oDG+n)3m){D$c2@BCjpESQw2kWfZ-r05Ao}zZl?$96D8aPQ zNs%U@^Ujg~KF`SsF|n?H1bLCq4-RvDd3f($B^Z8|>J(|Pv3ZoSrZO;UFLm`M>>iL8 zbm`LggZxL4_&C@!r{TEMua&+p{jzBJt=a7^w(C3O?|r9cEh;NaO$qkbcsl&x**bma zMtYYrnC(FKsw*yb_U5mSjm( z@SiiY&s=t&cqY@PUskV(&UOnr`5$&0ZE((D)axVZGoYtP?qFk~3RCXfIrrkwJF-Q> z{yEE}5u$#f(dwi8|BgoRQi}GL-N(ONyKe%b9&p;-Xi~B~O--{k@5Jk!xs{pO7P77Z z2}JAG^({PpV^U{wCx{|`;(5vyR2;=?0wv+q?1Qh@zie7qc8jE zJ$C5T{?mk;1BW$EdP_(+F1Lf`8^Gx@K7cj%}1RJ%c<8yF3XT?4SNf?NXym)y!dK zt8eY~nB2+4yOaN~F8Qb3F8+5BR=Aa)nHw~IY;N^}(NkRejK0uWfAhkjr|S&1>vWof zWB(AZneM3R>KpqX>+--M*CsP?o>KR@WxKmE<$h4dzvp>W5Bf8?e@N}U&MyB zMWU+@u|~KfnCz*M=I$TTCjaP;8t=unOQEZ0TW^2b<<`mPlWMifty4VG^UlW=E4TuU z4jB1ZyEV-ma51+&FEEIdCFsnNu_PzCN9zA>vSkzC2IkXk3>78=6j6;R=-HrKm58lI zxnJM4IcHQqYuA7d8FPcCVU$)^?W|W-d}is#IXFa*j8E+TP<7m+grwYozrqSDJ?Kx= zoj9=pR6!W97bNeqXzM+hD|=HLAW0B_0N{jtI%Vk)XRb!HBM(7PW!wS1z^Ho+_z3I< zk@t~{3-%Ed6nSp~&E6O4<-`J+G=~fZnzP z&zie8b+`CxnHfY4Q7yyGrM7{@SvcW7Okg!{@7x9 zix0PSx8D3`MBirdySq7+JpMCjuTpNRV+&rA7CH7u2=AZKSB1)#yoylotfI*|%(I1g#Db z6ZRrdSqp}b-zD}<0(n5PD10Pv6PgzZle`V_wDkM;vxxKn-Ft2AF?#e-o>L|mjGg+Y zcs=rO$?h$`nzUOPknnm#){tMYYjWe{?@S);hsy6PS1@k;C+e@VVs14jq3hvu6I7zPI@$ z*IB@U^cOcpcyewCgDCbe%8Z@iY!|ocg3zh zEZhG(@VADnx268tcYEuG2whZD0}E=j?5}3c!uQ>Y{axZb{^~SW;}PMT@$E$80sHGh zwkl=Rnf%)FGyFel=i|#*ft6Tn)`QDw{`IRGVcx|{W6bD*V#fiouU>gr<-g)*VU>ON zl~F&IUZRG*eY+8HJV&#p4YWF1&$H*xx1?f4h2ikOMR8XG#y_fZZFt{!rYq_yhkf7L zEM2i;$A4&5`SSl>^2o&8#C?sA`g+Cu24n0O_WrenTlA9o*8$#c^JK|Sd#mz837L%( zA5T~}e@6M4Nxy$H!da{Nw{bt)Epju`3^-nQ^T7Xn;s1v+qI6GjBjk!JL2-mu={9C< z?d|dZZ58x>XJ|nN^5^P0(51N1uRU{4#>U1LK9A`I&Q|^7M~2sh-Q^Pgx=ZBEiQ6i>V(moXjJV(BO_bwY0>8@lTqu-dy^9^g*_1gU8k0!#9`08xa_du z_T*56EYi+zML|Le;NYilt1exHsAeKig}q-C0Q7ox>vkZI+*7s%@884zOc2C9r-S7a zWT(JrBuS~4!_7e=gjpL~ac}69Y13XmYIPqC?eps!$IKm9dC(zqxMX#urA~;)7A-Uy zu;u0rOcHK%S$%cWHJ+1moR}NIt%A453`ji7uRgxMaRUT9!vwY|ZKm=bj0-oig|A)4 zn~ZRLxzjxI0|MDg{Hzc+Bx%T3>=qhD0#J#?50u)Xg|7Pi-Wjm*VRp~%{r3@tkGsmW zT%B+4|D`U3bK-ir2S>Vk-Gm26hVCGKV!*(gMq&~m3=g=e4KD<-#GZeAZoWqDDZ3vS zTHnRm0Ud43fQV-#p*czJ+9)XE2SgmYX=8(W@u|X8!~m}hYB{BX(CCn+Z7*sh#y#-t zLSFGgct&)vOmtP#)Xd=H;h7{FS03sM7MBdy$<Xr+OZac2*>{x!Z;og7%k#((OS2_4aom01(4?g^ zSzX%Y#q$Q_T9<6Ie1TE&KYUoe&?2(Bn)&W;9-=f^%4)M^{>IwjB)2$wHgR29&(P_6 zmQxF|6vMD2fTu!u9?NROegTKVJ^aIFeDi2hGK8d1pI6scuzL$x_LKZO7uXvXnC5ob z%iO2=oDSXRb}uzRruui$~MSg0I6ciev;jz8f#QMFz5Oy`zOmOnUthuB+#vente zOLiGGsqg|}l^85y!^$ zs;gOnRAw1_l?(%yhB=hc>x<5XIrjBZ?$ALXGpM93sEasW^25SgfFH`KrnMXsgxyRD zc+tUF3ttv4`->FOvL{D4s1R#*sg8Sb`O{I?(cK1^U$T)XV4qXIJue%!;4sD;6Y8o% z*>>E*K9~p$vANOUf3{tYz0fI=F+jk*DlKt~&Hnoz(%8eCk~dE67PODmZrmzJe)(G@I*2~`a(mPR^xWpB z2S!J0!*4|JZJ*pvi$3b1f9~p?i$UJk3Qh7%0R+lQOZj2*Cy^=9$6jf~yXaGVuN-z% zvm#pO)tB*C0;=7#e!h=RIQvqs_D86lFRTAu=kVA@D+x1WbJl*i_HSn9J~|CRYuTsE zXx~V5?ML8%4*0g+yjS>)qDGAxHE+@4=WUZRJ#vM{O`E>a{$GW zP@z6ze_csTWL%e!C&Vkwh0W1p&z&*j-&VhVk*EFltD_nq8RMxqhPP?2{Ihr8zMoj@ z`1p8!+qn2Zm-WZ{ZS6tAOeG{I&9Y*5&a?Z)?HiCXlKMI(7{|B-zx;GXH+4K>t7{?8 zI@y1}+uILe9HqC;rMcnNe|s7jh@g~U^lioWJ8`PRU%0rs z$~$kb@&QT7xws6Ut|urjy7*e#;`U^zb8)UH6jm3G<(ntHTT*B8-3KuZY7vUM+zyrp zZrZ#Bi3q`qaO?K%uMN(Nh>F+SogRI3t$thGz0@$Lruxaynbqs8Ul(10Ht@Lp(Nq8G zS=gV}T}Mtt=SY5k^<8oS2dQ>y_WHWbcD{TTKp<}8$NVy%MQ4NlmA3Vt26lCy2k`=h zbc#SNayT8&lP~X5L<{vRISD;L$sspfqUlcW-KP)36;@w#dUvaf*J=0YAuzp=BLdK; zby@wg^v?n-V=gDOb^7=ytuETUQCkH31en=JhHPunO6BFFsXwij#h8uU_+yyx9?7-k ze!Hx_Z~XSwdq}uj|2Z;#Y3+|jLS@kOLrBX^axY2hki4|Tw*m9c8ByWQL+lG9tFc20J@kyZA2TCn5GyPEl?Qbc|Ke%x??5NfFCMgqq^Ja{f?Un_kB~DD@V1MAQ^FZ6-#Va64{CB2})Qm}5UJ zIvcl0a$!h3Dwn#__C=T#LT)0LGisS$z(pE&jn_h6%!_hG!ZRR;H_kabeA93G5UV*n z2iP!V#v)F!Vtin~kkh9d)2go2khny(bH2!87O+0(63SPAKUDOPv*xYbmk;Vh`J6E# zZro+N{TcaZT{e&#GYIAU_KbGA3<#M4fwMfebDX(95rN;z{&j4t66<1UV1Pa46e_5q zBHyH>q^!XUkTIq-Byjur-R&Or(=lR~rN=d^^)(O-`ZXK%)el^}wp9tK^t*zB${myQ zI{Qz97UPo|Ebi9aG>)|+$LL6~>PVtIddU zqj|qt`D|g&fdl8z&fq)3;l)lgxl`%*Sr{38mWP&y{Gc6~V%w0A5cw&zvBq@HOp-Sr zacnBd1YT|1SH4P&lQ#z*ZigFX(Bn&6_t*0S*-_oaiL&;zK|z;suyVx+;SEkj#}S(DMWE4PllB@*Z)$@}t?# zkJE=Wdl}q}4}mi`YXMPIv>m7DTI)G-q#0_Ro!M3!^hFfQlBj9!{t8&_OG+fh zj~J%*A&VJEI0HHN`>I)zNT$3KIHhm8$;y@CW=U%!je}^&7)-}r=J^Ff5gYeJKC5lD z?%kjJ5Ixks%EC4AZ}dM>pZDz1VPXNXOcyNS zfR*7!=7)CG_t7JROw=G{-1WW{A0KL!Y0+QiS(#-j5n8WSg~rJ4Y0+o`cajH9<96}; z_um<1B-J}9tISfPsdIdIFqW)v-&`15x8Vc5NV|QR9kDZAb;E3aqeYzg01$mG6j@5j z>95g8XsE7#-_|f`t%Gcd(QCi-MF^xuYs%uJ(bz)HF=++`dZ~7j(wc3HAz`}a%CcSTM;@fMv~{3S(opx zGleo@s*TM_sFKK?d-j+{J=FR9BYs!C0kU~naUPoF=^{W!$iOG;0PbQzdz3WPP!N20 z82k`H10<+q(0FC#X7vL;5y=&DoILt=7LmhIs+V)x*U)hgo{;o|(%~@OuI<%l+aT7V z9Nh)RbkQz*;JxLsUu{bgF-$B++LTGps!U!AHw+cpYC>&kBXo94moG18r7I9zl;0cr z@B+qW(qaoC5S6Rxj5-G`>o}B>N66ffzf25!Lgi+?_u#?PPis&|nUR2rn;}xzJb9FfW$jIW6Q6WD8Hnj zLCltiZtcQ=CbSTYmj#h7FOlS9Tfaxj0_q$TKuZ^P(|iz1SzUo1mlURjjDB0V`Wj!( z^XBnUYk+#T`+gk`@=D??%PS8{T#9L!9P+RF(dKohakp(|c!-xw)F} zr*MuItTHw)rqxTD5xb3ceM|9W!voJq2FhJxTxAO)%`uiPKgdtB zn8kmtWc`4g5cZqe1edtgf!) zN#$j8D5aB^h_#IBTyi=51?49wA{!v85qOAXLW(|-vd~T=LhFNey812@;fR~&ka@iE ztY{_HptRje0twD7hYL{*B4=ZlX!r{es##Z+oZDs1P~W$Rprl)rgsSv~U9>auyzr4; zTKU-^OJd2n@@aV6lB*T1xzv*4)?|(dz24Ci#OxHAcF*zSt+<8c@{Uff(*%a@!`asU zUi?anN+iIAs}U{uRJnOB*lco{?r6s<$}i^8bP)}H4MN}tJmUjcz7xK3`cB`ZS`egst3P+h@2CXamV)SUil^t*Pfx=3 zh^%qzou%!#MONipIYD;e3jp8;OtNglfCU$)&?`O5lX5_}nHoiYuSD0vo{)qQPmw|> zks~Z7{NbkYFGnCK$G*d~P(xgq-j$TB&T2S7v-lP)ii`c3)2DZ$5MMy8ZVei`FjBFf zH>4_PM^)bzDt+GWTy6#0FSV`ffLBtN>=~zy9zDuFn$FUB-@fF`r}RjXW9?m6(LIyCy%^gp&4Q5&FRk21 z15HqMVOVeqH|!K6ql0YhZ7!ZcBb{{##lEZ?CX%BJUoUFR;?H`hKXy$%fZt?T`kM0>DaiRHQHjBA4Y~U{O(NaS+h$>#e$Wl(57|>aH~P{Q(#t?(yU%=uIi5 zD29A+w#C+)zXgE?g&kXZqY(v&X!=V_PZDUkfmxCz6ogimWNeincKj0Qo)Y zH3as&(PP?k&g@2N3OK64r`mGUrgErlLkcSs#kCKC830Y@ktU-qQ`*p~PcN!x%xBh9 z5I3>T)HFX?MGrZ}3XFj9?p(@gX-}lJojZ_hWm&z^Gsq?lxGExhZmqg|u-+M6o_%iU zXtIzG*s@w&_r&x21{umoz6rjSU=j_1hnvWZLupbE|MYP!4}zT%;rUG$SRANPG6&iz z{%ay3Rl}vn;Ib){86fu6<)6R5V$dKVupKor*~6?+4@eLp-`?R~N;v=+Hiq;#tgAJd zyLcLM>79u#CVA{UX6<>G*+%0RC>F*?xtyq}sWGIS2d3YCXDR07M@U1dJ+^XL@!08H zMW;K_!^6YhzR&qx)^0^lV<_jkS0#-03L3xhD5>^iWW{2h%k2~SY`AB&lic^_<|;PjcX%nKV`}ncZnuG15SjeM3WmI z*tgGlR+-|*xVknbH4Ou)-v!9uL^gK}xb`PaB?s-s2B|8=;}0 z?d-YUVa!m?_JAqeJt7a5ikfwmjl$wefsfNYH44w_Nr({SnGl zX)|hEdTwW6-b``TQg!TME+yH2RJ62Bm>d$Q{71t@Y%btJli>ACtWHI`ftnOm0Bxt> z5P*#bM@RUsez9WNvR5yb1ggM@OUo7-;8~o60KdJ7c~cK$%uW1Djqh{g&6_txc^|L| z$gwxMYz(Ez3zIqRFZoz6=H|{Deq|(_rAd$lb+d z-tNqdj&>1a9JeY7gZsOm*Kqe0y#Vb(dhGgD*lx~vZd|k9;7mqEM0wA**bvZ)GY_E^ z{nvf8wcCVdHd9XDO4pb1mv*{lsZ*LZZHn(;z_HkjP96UK`*iOX0owY?9Xo1idu8=D zO^;SoITzFzc1}k6DkeF!A1=mfflGbb9@UMnssoN4Q`}Q5?UGDwSszn&_djZjT_gNm zz;omE15MJRGI~ev2mB?1!52`CMQyOGOG+E~4;}LX0^)moTS-PJy0Du!Z;q|2&Ie$f zJiP3#%^>idm2NK^fKBw=iXO$(#~Kqxll04Am_$BiiYS(qTqZRGr~z1_F+l%_^^r#4 z*L#m0^Ap_r<=ux5A8LZxxCf-qFKCfyYF*V;vxyVOvbh?F#l?(kKuIJ6SI$(sJ?!(S zgWsRXbzin{2A<<{Ilp+J#4E06;J|yDu3>vYq92< zMcn;}h3@*hEG^>NWZAmKIj5~0Mtl|RBOOE`-t~2Y>h5C1JwEm%w-Nn^k|-cCJlT7A z_Lyw@d{eiw>u0-J>)c)(Q%c1y)Q`wXu!j$DQ0dm~J`xC#T^b$t`{0|Y5#@Hm*vdWH z&AvV;K8w&RbsibiqVo55M1VjwAv!dw|39EFMoFWKNJ1=>lcUBv8YxFUvT$RgrE06C z41b6NkL&F+JIZ^!Sxyl0x7=TERcKzONL*89k{?*PdVL!woLAuxfSwAR)s@{n!G%TA zfp8;prIBh!%oO^kR7R-XZU~Qnt5R{D=0evfs|QPm5c4*Ng=rM|4(!IRq;gO+vIu*E z7m-0hs((_XOS6!ScR`xKBBvq}HmqTvB3u_3xv)^cc+w6`5S+#S=m>BO5jEiJN8s0+ zHW!Mv8`2ujyBF(OT38U>jLjQWl$RfL_~DK^!R(sHaOv&hvUx`WYGOYBT*8LdZQ7*c zvT)=h9V@(YK>jq#GPCDdS;|CJac0@=5SqEuDmp%KOnnlP`H#}x!%$Xf;X$QZimtWa zo*Q$IkNiV0V=4PcQc zjG!8MkAUaP26c)^i2{tqV&b6RbXW_*m)@6%uG{Z@Pv|#bKq|_KIP=I$wx*7^^@yY+ zM~_~-dUfjRq&ix6jnbks-tjNlS+54! z=Lb?6ckoM%7;eyPl9$)@iVbD8KX} z=i*JyPmLTXH0UPKYkHLpoNYc76zrqT8zc~)e~y)OQ}qwzZ^#kARTDpKT2Qm@SCNjL z+vslat)U2hNFpg1QBUTJYSxqL2g&s?<})QFE-Tsp zb{WRzA+@~XSu+HJ!XV~87aR~Xpr;aYBZpe7-hbS#4j@xkP6m&n*R5y&{t*tIi$-(R`9m0N(zN+pNDHL6Pz$4NRcd4wiDR%omd!;u8my>q*A@l06> zrir{|zcTKzxjxnbp{t(5BlH&bc zeSWrVr^?D}XR7~?&~PbAs%uscfwyF0!k|HW5qnViPNPqkp)5BsiNIz}^zuNJ)5&lI zmLpq68cn&h2^#{jsOsJOZDHZLZ_VfE=}EH={xGC+%#4u%He4Rt3D2QfwpZx+q0Q_J zauzN|K^`I03u!w9dOmL$`6$Box)GS}t#gYFegX-}>G|e ze;AasR*X96;uCrMJE|eD_`dF*_>-uU7C{4pg!$~*gROlqeY7-9`Jc-xI75Bgx8C~V z(&kVC@9C|jeo4&ebtGw~k>{=+IP(+a#0>~E^L+2Vyf47_An!U)HvOkpNjgya8TnfJ zsV$UDqb%ORrSyam2Gb5txLZdZOUF<;{&=)q3G&F}&_^4IIAC&g_5lVR026W!m^#(A zd-u770V>Jk#KjzCF6GJJEBASE$_r&IcdiNMn2eL3VWFW7$l{GQOipu(#MkV7^zZtwabcSYnBxxz5#=^jU>AsMJ6PXlWMEpO9Y&agsnE@%b z@xNkb$~CrF0CZG?^r5&2sDwNhg(1v75W2)j5X)3RWNCV*O`x%?VJICffL7v6;q|N1 z^r`J|emM;X_`(x0%#Q;{eE|MSRLh6A(CFX4rGBd4O{Pe33tT`=MjMtt)u_wRg-=d6 zktG@8(?G{7n`HiF!FgI%!~of?Eq-~2M>tbozko)I)4+*(ExtiO>I0x51@8g( zx_+z&xwLEq9-6ctaMYwK8ELdDkvr%nn^IJh+#D3!?WNO(kB3K&Y3PGgJunU+X{|Cl zwjlf?f{@jerUyGb0<((60ZLp-SZ?^KP@28&5E?BT?(X!eX>=V8oFMNi^ z#*IHNwM!O5Fs?RPYZA%f%a4&LBQIL_G00tf6vfzU#fn((dN1~2X;qcT=j@rDNMS1x z96_d7jAA0$YsipxybYK1ztB)hy5JT~(Mr%~iULy32tdBH`4VYbb{g~#5DPFNGb0O* zqFh*a>na4E`1Ih0aEqIrY>DqM4&G5@W1uTySjJ1rTeGl>V$7syU6mD~Oc>&6_5JlV z_tL5*S>u+lhtVHe<}4R%lTvggo!`uEGR))pN&}S@be(6ZFLQp&BcQq#fbqMyl`O< zNd>~|!14r{utNkrI&9UcmxaoI|1H-Z`+6G<4#;_c%#79_zsSRYgn@~ZV(n*1a`1%< z7gXGDn2sjmK-luUkoHX>TR#UpI7DS4J{$Zs(g+#SIAg}n#~y}-?8{VSfIOIfcTvd* z*@VFPAEnK3=!y;tHDA9TM@s}Op!6<~-&p|@<1sGsO4FDl~3n;K7?TRvxjyv}xvK+Eqkpc4M0>*D)aNT>(2@s46@WG_i z{3qnnY)~1NN7^9fl3d9PEO$M5{8$WssAdQ?ZJ8*GyMUFc|5C#PWH5*w37Wh=rHdZtzy@=6EhTlF|j^dikeAFIV%uLhJ;J}8PMq9->|GPLw3Ew(^xX8Rm#L^f=r=AyF z{$`4n1*h-eqR&!9=AKWiuDsExjnc}YnyIf4@+YABB^SLs;1q8l%g=4Iz^#n}^A~-H zcFHSt+Q`|{bnZ_iOa+Xl=|yHzat9-ilE!aoN>MTV;K9KQXTiBb3#MJ(KP^-Z8l0>R zH(j}HWVKNZvDy>@@P+%@C39Oy`wKB3j1D4$m_59JdMhKtN5v=qkL3TG&i|vO;Q#*Xs(Z7x;!**VOVc)_ zFn2;*C8d8b0i{06$PH#=h>g!Mm6-V7n;6;71kcGsQOGPXQDkwTqj&wBk(wKf@IV@* zC`JI?q^!<2RkRBL+lP7llw_61U=qy{!o{8dbao!yTrB>!)@s=20A6 zBnGKOb(7qum@ZJA5Ky9jj<}v$YfH0v^K@-(Bg!|m zAP+GgQ>Rvd;zIL8(W%U5s0+Z=e z1UpHibzWl}h&kNZ_f}b{aP|8&ub9$@#SN7Gs-yS&6EjbK#HpP7kgL{kP97=R6EIK$ z1;_+d^j@6Soc}WdadLJJ;oOTJoUSd)HIv0543v_=bR->9;heZ&xV9#7b`JiL;$8pC zoqS|E%Zvz-G*A`S2}DB{p=sVY#fy>vu?39{+i*_DzYU($9Wtq_nkXp55FlgozAZdi7?b2hY5!Qc+v_-t8S=3zSVEH zWRG22vK<=lA_E?mnT?j7=x334;q+F~7_Qq4i|QI~ZgCL7(@c=54js~9_Reo_8*4mj z&f0=bUv-F)@Z!=UByzALZxeNfdB=%49<6}503J*b zuOS}pbhuRpBuM`ycC?$pz4=%0)J9a%^y3~!kWzj0r;PI^xQZPz(G?M|>4SINU?KD< z7rKi1x0tt^#3-<6!KJ9Nwe=g7_DP!hM}|zT={tQ{Z4opDAU+4FF-8T5mLG^EHRu`; zVmg&E!0zC6r8Z<6?)Q`X0Y!=ahl;Up8API^q-FGn_hH|Dj;4Fuyr$t|9)N(_xH1jb z^s>Z+L>?db1!eWYEeEzC*^NqcmO47a7`Mieu{OK>V}u`TS_7n(Tg?%`1~LF?BV~ zGljZFz8M02NY2G!JVINyq$f1+;!i~3h}8QuIb(M1!~+axB#jfD66;Jnp%@k8Et4`k z^ci{*CeoR6z)2hG^mZaA&DN;)C{i_uLu7}F=g*tR)ei%TkF+@_1FlK3 zgg53~>;xl52Zc%U0n#Qu>%oV{6ClV@((dPwOUDp(#Y~EgsP!x2@$^QHGk40A#w1F; z>E5JJa>BGUaqU<9I2#cQ4B@_=rtar3p z_XlT5>MBAisS^Q6Ovk0YN>Cb^Rm$heV0Nvj%^Edn8EzOk5~xp%=r&6>!EI>6d-LFIcAJZugJ~Po%+vVDjO*^DLy}(3pZo`<)r_t8I99_YKMgwMrnmhBRa2L}5Cx;KzOG-iZn?QJy z5}1V$8f6rMxYCmT#ZnD|qggDosMMo%_#-gJSuv#`Q)$|Chpl!y^1lxV=}r?f#lRGsri*+yDdT`n$X`B!)QCuzL|g3d2D8PX zAZuvPn?*(lBOW<8M)4T+$vPHSCc(-k09cewK#z(5Dkci({NsY5x9ofiNhQPDO?Jm{ zG_5OEOM4#hg*d)*5ed2gS?F_T>AbX0mY7GOh!LVfM*2!em}E;5?|9x=p;>*cOX#9G zVgGP;rg(D$+>5QW;8vVoqw5YdL*|I`YzB5@??`~ikETE z(h&zfhVR`T6gz#(KBGa?2AiJ zfd6S8^a)%bW9QZ`{Ns1L=O5&E5rX-1i-_+9TXz~ynsVwbFpG$3gq4R}Jj^0c=lb2Hw%^|yA$Bs?; z=Z{7rOGB5$XZkJ2CYo7`ocY$$AwGz9n&TxO!}1}s=ZDDoerusAly^cf9vhzBpW;YH zwfxQf5|!Xm;Oh{hq`MQL&0rW;k;9WZ+WX1tqn4+>*$LZjY;oOl0b(;+kp+xRZ9jLc zTI=rP`+sA2A~OC@&tpO<15ZLU-a=vqg7g=VkEp760QRPqOHlp8Q4OSDHrY~#8dp|`{rsRg8Xyd6XfOM7<8GZ3)*h^}1Ls+p z5{T)HQ*m)7kViN7OYKSJlPb#&)E1C`@hoUp zJMC-V_IDb#axcjCdDc%k*OgoI*6Y_+V8|C=RAllf;QxFzPOQ^K>Sevt_plD8Y?5yL zel%$F0X(KSIbC4HeeA_ou@^5wDr~LhDNbFqDCGW8yAy2b%s4yKvIk>IFt%@*yBHxL zdMnY{wgZ}O&&Fh57^OtB<$O)8)=%X$GYM@gf;-OisZL|pT@a_)N1rQ2U3`kNIqPaD=xtZJkgmCD3^;Z*p77&5`%3k5d6>@}J$(cb@Wbt!xI6tGAhy=PRY zvKzFomg@TNu^t^X=jx9lg~V(VkwMKx(`YS0P5{zdckaxTHy-6^PxpAiRVKLvD=C+K z=HzrFjJM^w>iC}$bh>o2aPk+VMnib{(8{rs^5WR=Hwi^ger%Ws=}K|IZMSe^^Fmv$ z5)grc`B_=W#ryQ?wFN+!GOU4^_Kh6b31CLhA$%vn9EjIqj`5#AS$~TJrRmetHLaIK z_|TbhFkL&_Hi!>BJvtlJfR%$?V`Xdo3`Q~iczlcu7R_BoV7}(hzvq=DOdv6YIU~85 z1oi>{cYbiBPL_B>zJGtxL=k!2V*PNn#^|YMU|{;S|*xlO)hcKDXwZBb71uDZ-&WT;K(FO@XaKzMTQxwhz;RLa4wg_|!4kYXB1^CHVL6?BJm*wO)T@uD~)V0YEg-=@VANyYBFQ znx3c*cW{1$Y$@yt)7g22JpC2@N6YU#K=Z=G!A zW6dV+g-U=gWW1sF+26X&katlJ82|WubhbLpYE<69jt4k~Jl=?LtD|>YlXgV6OV)hT zIbxSJrvkth2#a`gni>$f;9Vk-?@b#x#AX{d{M4K~&3{@i9tscX0s>&Xma&_SFM4Qk zE+OGCu-GMYhk$CtQ|Z)koJeUPL^DXvnI=Ql*-gq0nbAOGiT&|(Ngtp!50BBo8^3oZufMYVv6+#PC=NhXOvc>8UzGPJD*k~I zDvv|-k)V>%VvVTa5&Ih^4tQ#WW}gowvQ~%zMXN)!l?f?Cq9-Wqiw;^WF{p~D8}>Lh zK)rFn=<4~#kBE3DHhy2g^A@=QACm`k0v^cpXNBEbq_$}wjtawne3;$rCa9Si zO$y=+z&X2Yd@`?zBPgZ7lccM;i^D%n;L9V|+>$jI%U2u}7$DX{tufU!&xQ#Ha^Is* z{m8P0vEyN=!lsR@yVPf^4Ob9K(W7LP=9NJt9-#KG)nLBK6=dEKXW04vFKaV1WHsqo zSy`HE#gq2!+mk6#D6FNW2gpz6p%6aij<4B@EdVW7TU5Fwh3~YMt^?g>WWDy4KKYX2 zqJ=|WBsmg&VV>l*+4h63WL9UP#gwTDf;F(B0E}Ei8^_ zy}5#|8&!cGKgb2w4)kO(@`9|_Cy(4;j#!H738Rz+uZE#+f~`HD5>1)0iC3;b@YSzm zJHeY&&-H$Oyt()$m_LW2l+2xO(vxJez(^#moup(zzh<4EpG}H8Lr?D%RXk@7Y796<22==Bxi-eu zI87UIj%L|i5U%xAJ^J_GB~>(`K&$Mtbg5AFokmjNjk52RU)4XD%v`X!*^ngTCMnIL zh3CKRpfxCrVFaqwaOc18E?dr#Ek)|hLA%9U3x9=V;Jq|WAy+2_%bjYwOIlcjE2h2JWk3%zB z=oBuvn~aRU{1=Rz#P> z=N5DxoQ2T=e7F&D@6DTFavB-(K>jol@PSR)oB2Ow8n2Zf(gnoY;*h!@^V5sqD{S&~ zauq!&IjsEG?o%nCb#Mj8RElR}#+!)#fZ8$%g&!WzC_`#2F{KLLoi64th_J-X8Cp~P z&>x{`J>gN*uTR?u%#3{LPU%S=PU#h@=Vwd#A-!xEp1eJ>;>04Bwzjt7p?6UyCEW`Q zDZ>f(7D=iZ^$)!U{eY1?OK&zEKH3X8NhzIKsEYr^rEtl@*+2*jlGYB83WgHoJ)c(s z8@oC;LWf8pEgW11^ttwtL1$+Itrin5#Vni`8-NcA!(RGnkDb@n8t`H|_D3AcILr5k zj87~fiI+De$1CNl$I9|ll{}9N=M%ccPuK2Cv~@c5Dpd-dlZ-MiL?Nn;YI2JYF1se z5)lp|O&*!FEaxwi)y-o(NoWu?4~d?1jEd@%JmLBEjY(M-ho@5M ziqu?5mj&%^^mDNdHzEyd=Px=$>!?T@0;>Z~KHqGxOAF?owTg|eF+yXo!Uu2PqN=5{|UIc8S+M{b5++~y;qLT-lY{j)1 zqQ8z&xU3waxjU0aL}+~a>cw3J=a+K!m2j-&kg2%f?0`+@YKU~Pf7V;>R8~2U>*6EQ z*K?1L7Jm9P!_4e(`G-0*m3PSiqSikk!dmNiVI~0&aFYm&XiTbh0X#$;06`dp-1N{u zOJA(1%vb9pND-$L=7wm=nZ-um*!&~6s`Q<}sxL=+Ezso4&d5T}KU+9o&~BEwCHhWRb2CNK?;TV76Spr zOGRR}8-UyO+VWX+2hRW{mp(i=AE?3Nxv9|`)Hij|y})2X$uA&t#e3$Umu&#TPteM@ zl%Y4sYsH=kUm0Oxgx#t{uN+ zo!6}D!-{8=G*)XYp$2b28z#01JV(j(#C9Rqg{aY>;8#9A-CPnQ7sGhtQC7C|=It-v zDXKafk0unz)qf{gm>}W(%iF&p`3!eP}Qse z_CUL$$oTV@Etsv$NCWQ00q0-~1JXB03_4@-K75cVIpj{DTbeW7merf_bts#^MwBG- z%0ri=wFIW)0zAzLH2HW+X)~M}vVDD|WWkR;g)}?RTB#&1x97`PR>)i zWSQZF2R%J51se}%RL>8}s7C4Bui%tb=VUaKv&JA3=5}Mnz%+|)BIWqQG&LGZmo#F#dF1uwwRo`L5p1#{)>FrL z{G#yLC#FLWfeBo^hHPOyo!tu<4kGv)DUdWC-vDw$gu=k7YRy*eqE@%KqN}I(7k>is z$M%o$Y+6E(CB3MwT;USEVR@${ZgUwWA+$Pi7o39`IP`vC3)M^w=38DC=6)*x!Zy-f zZV~_kqJh?!S`ee&Y&TFx(N-6_{*MuunEJvbOa6i86;4V8sVG^*0;G;8vhfuZAs+~) zH0pu`o#onMpOv#+1>&VPi5NTnr@~16A#IrBglpiM1gaCky0kxRfwiY$ftl?*A3xp|PYS5S7C0)2zL<;d@mc*Y%U))ZJ&Icz zXB!?1B|!}L5&Jn6Jb>bB`k8r>GQ-!ZnQpxFAzIqQ#$z@nIw7$o29t-E!xz##lJv!Z zkOM+^3mQSXO{ikhNHBBGiXQs&CvZjI;NalLIupPEaOg5Xrib zjJo;>Imm${M|_ZiP^y6=C?OR^orH2mi|6_q-sA!x|9+FzagBH|a;;AXH3j%@iuM{& z{@BZ!*8ejQwAb@<}2!zNk|{$dPUC}u$k3E&KZSE3xf4R1>tOTLu+wxs0Y z{JG;WRB1sN0-EVca}))t{D(!tUefYGpwQGM($h_#&PRHNs>(RmtuXodcD3H;6O z6-=MvQm^(gTlD{QRbs2L6ipEzQ1{Txe}SG?7Zfg9X5WIjghY5y8&b5e<95LHCMOz< zEC06_u*4dOT$%Th5{*>c3;Qu~P@biIjcZ4$JsGCNoeLmQVqQy2@H=GxHr#ku;c#$o zKBZTi@!-L7f;uKBSZgBn5Q!(LjQEg1?#Y)GZ658o1S!phm-N-!OBF|p=OU^saN-6? z4%i~MDM`^nrR3x+=$w7s5PHQwewQ+4ay8Uq5Vd`m}V(+DF-PJPv)Yhr)7p8Wa*4|ooao-+25>K^o zZ@z0t!%Jf|3cCAN-kE%H$H$ph$L~0!pAhn^^2eJhw|eckuWzXG*Yj0&m4A-xd>S$# z*ud0uI#CypG^m>eK~5NPZW0*X6KrqBO*vh$7#!%EjJg|T73bUj`Grkj?-E{Uwr$>I ziwp^Td3Akz{U^j#+d(Psx3hS^j7&}ETU$S*$x-tt@tnKxhGTqa_UEitow)?Ey6A+} zuWqs9cET&i6BE4Wp!DQAw=FA$Oaz5Un{uhmXA;w*yUWT)4Ga`%IBEvnL7h_GIm!GJ z0>^srhh~i%3hRT>&(LARwsPZQEY7ALtJVWwr!EPFoGUAHevz9i2|raXCUD5sZ*V5*@gW9ma~e5&$&y}iadDU8 z;?7MCxxYN0;vIW{om;myHa9n?(ftEB3mCRTM%#N)(#Q-gNysI@ML3yp6Y&K_kB`HY zOnINubKNzoxE{Scmb|Ad<+K()^9&6q z!F!`$jfoj-YiD;(p_u$+SkHaU;FeN9R!9m8Y7WwZYxTaw#L3Ad1-Vq2H*;d-duyd4N-^ zCivVndTu{Ld^$m%Sf<7e|Da3y{B|#Ho~#4Kl^#b9A3hlw>3{U7iWijb;aaSLyRrQZ zYYHlA!q1;y3IP}Ws_)o^#jmf4V7mD>-Go=ghb@hz2OCR*lt+(5y!MtPW~;6m-adCR zCXx1IfJJpive6&HoEVwLB{ZHe{#{!?$;=dtZgIH+K^%1^M1#+io|?(EUf#34^jja* zBZ)!N16b`M>Rl!|uW1olG*aCVK2(Q+wZ#9(Jy+BAuC0oLd_jf{-=$y4QaF4czC*4_8bn6l9^ z9%35lG=5?4q5GabeJbXhBqU>g4O)Jml@5=NHr-ac^lK&N7&zF`#q{7?;G^7JzbA#B z4tp~MP`a?Fk(7i~|7_!lPmSMMbtMDk_W1wb<7g z8F~<6q-w|qTVSgTxyVs`;(aFJ_0fhuIHu04$Ik+#|GIX`h+}BaeRA*qc>yJOT+J`8 z1ah|!m8{v9_Jk=EE%cuZni3zgW^~Zuv%9|I?ZJr~JZR9w-t9z1BTdn9=2^#zs0SNR z2mMM{ob-M%<=dQ`b%7(gw%(EGbGM{QF-v299}SHIfq^?%(_tQ4-mwk`H* zQ&Z{TP+VSaw`$d*rAxO$Lb!eW!RFHu{{Z3yJn49`})|_&Y6FjUmZV=qTp^pPaI0GUms71%Sy1IWX6Al%5nFq_CII^ zf~^_ZwRM0#SLdnN*aLNSUVPRO#|C`(_|ECD&7c|39<$ey@lncitbN|74}T{qed%v+qRD=&ieO;W3ps?4F4d zWX+3s{3Uhe1+z|!U4LqPddi(>(#fk2v3NpokKFrNyR} zR6Tj|V)U;?D52MrR$MTTY}uv_I^SD2ZZx0+zfDXHB2$>XcDhI%r^^mU{Au z@f-D^n@0cfwxTBQL*Us}yZrsn)#CJV#G_=bn;U3p3*NVWnb!dNSmJZL#eZGHvyw^B zl;w+qhMcL9UgzP7eM*{qOKz@6X`aC-ePf)*VX!e z48Ri(lb8wbDx=(?SLM%LX8!NLZ8_urmp|3E(1>r+dvA}KVGa1-%&BvxMCsaW{$K4q B(5e6c literal 0 HcmV?d00001 diff --git a/docs/source/userguide/preprocess.rst b/docs/source/userguide/preprocess.rst index d22e71d..864a7cd 100644 --- a/docs/source/userguide/preprocess.rst +++ b/docs/source/userguide/preprocess.rst @@ -15,10 +15,10 @@ Each of these quantities has to be set for each channel (mostly depending on the .. code-block:: python - data_cube.polly_config_dict['first_range_gate_indx'] - >>> [252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252] - data_cube.polly_config_dict['first_range_gate_height'] - >>> [3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 6.75, 6.75, 6.75, 6.75, 6.75] + >>> data_cube.polly_config_dict['first_range_gate_indx'] + [252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252, 252] + >>> data_cube.polly_config_dict['first_range_gate_height'] + [3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 3.75, 6.75, 6.75, 6.75, 6.75, 6.75] To check the correct setting of the bin height, the ``preproSignal`` and the ``sigBGCor`` can be compared close to the trigger @@ -65,7 +65,7 @@ Where the ``preproSignal`` (left) still contains the pretrigger and the laser pu that contains a 0 count rate. .. image:: img/first_range_gate_indx.png - :width: 100% + :width: 100% The ``first_range_gate_height`` affects the range corrected signal as well as the mapping from time to range. In the original labview @@ -110,6 +110,6 @@ Later, when the range corrected signal is calculated, the the count rates are co .. image:: img/first_range_gate_height_correction.png - :width: 100% + :width: 100% diff --git a/ppcpy/calibration/select.py b/ppcpy/calibration/select.py index 183e5ca..409ebce 100644 --- a/ppcpy/calibration/select.py +++ b/ppcpy/calibration/select.py @@ -1,6 +1,9 @@ import numpy as np +import matplotlib.pyplot as plt +import matplotlib + def single_best(d, name_val, name_min): """select the best calibration constant @@ -40,4 +43,67 @@ def single_best(d, name_val, name_min): min = np.array([e[name_min] for e in l if e[name_val] >= 0]) best[k] = val[np.argmin(min)] - return best \ No newline at end of file + return best + + +def plot_cals(d, param, used=None): + """plot the calibration constants + + Parameters + ---------- + d : dict + the dict as in data_cube + param : str + the parameter to extract + used : dict + the LCused or etaused (will produce a dashed line in the plot) + + + Examples + -------- + + >>> plot_cals(data_cube.pol_cali, 'eta', used=data_cube.etaused) + >>> plot_cals(data_cube.LC, 'LC', used=data_cube.LCused) + + """ + + channels = set() + for k, v in d.items(): + channels.update(v.keys()) + print('all channels', channels) + + for c in channels: + guess_yscale = [] + + fig, ax = plt.subplots(figsize=(8,4)) + for k, v in d.items(): + if not c in v: + continue + + if 'db' in k: + marker = 'o' + fillstyle = 'none' + else: + marker = '.' + fillstyle = 'full' + + time_mean = np.array([np.mean([e['time_start'], e['time_end']]) for e in v[c]]).astype('datetime64[s]') + eta = [e[param] for e in v[c]] + if used and c in used.keys(): + ax.axhline(used[c], color='dimgrey', ls=':') + + ax.plot(time_mean, eta, marker, label=k, fillstyle=fillstyle) + guess_yscale.append(np.min(eta)*0.95) + guess_yscale.append(np.max(eta)*1.05) + guess_yscale.append(np.mean(eta)*1.2) + guess_yscale.append(np.mean(eta)*0.8) + + ax.set_ylim(np.max([0,np.min(guess_yscale)]), np.max(guess_yscale)) + ax.set_ylabel(param) + ax.legend() + ax.set_title(c) + ax.xaxis.set_major_locator(matplotlib.dates.HourLocator(interval=6)) + ax.xaxis.set_minor_locator(matplotlib.dates.HourLocator(interval=1)) + ax.xaxis.set_major_formatter(matplotlib.dates.DateFormatter('%H:%M\n%d.%m.')) + + return fig, ax