From f4d0477c0754118b868b4cbff0646ceb45214981 Mon Sep 17 00:00:00 2001 From: choglass Date: Sat, 1 Jun 2024 04:16:15 +0200 Subject: [PATCH] Debugging JEDJAE ZURVIR XAWJUB --- cell2mol/charge_assignment.py | 25 ++- cell2mol/classes.py | 99 ++++++++--- cell2mol/new_c2m_driver.py | 162 +++++++++--------- cell2mol/new_cell_reconstruction.py | 43 ++++- cell2mol/new_charge_assignment.py | 114 ++++++++++++- cell2mol/test/check_Cell_object.ipynb | 231 +++++++++++++++++++++++++- 6 files changed, 542 insertions(+), 132 deletions(-) diff --git a/cell2mol/charge_assignment.py b/cell2mol/charge_assignment.py index fffc3c58..01d185f7 100644 --- a/cell2mol/charge_assignment.py +++ b/cell2mol/charge_assignment.py @@ -52,11 +52,11 @@ def get_possible_charge_state(spec: object, debug: int=0): if spec.NO_type == "Linear": possible_cs = charge_states[2] ## When Nitrosyl, we sistematically get the correct charge_distribution in [2] (charge = 1)and [0] (charge = 0)for Linear and Bent respectively elif spec.NO_type == "Bent": possible_cs = charge_states[0] else: - spec.get_connected_atoms() - if [a.label for a in spec.connected_atoms] == ['S', 'S'] : - possible_cs = [charge_states[0], charge_states[3]] # charge 0 and 2 - else : - possible_cs = select_charge_distr(charge_states, debug=debug) ## For ligands other than nitrosyl + # spec.get_connected_atoms() + # if [a.label for a in spec.connected_atoms] == ['S', 'S'] : + # possible_cs = [charge_states[0], charge_states[3]] # charge 0 and 2 + # else : + possible_cs = select_charge_distr(charge_states, debug=debug) ## For ligands other than nitrosyl else: possible_cs = select_charge_distr(charge_states, debug=debug) ## For organic molecules ### Return possible charge states @@ -359,15 +359,20 @@ def get_protonation_states_specie(specie: object, debug: int=0) -> list: needs_nonlocal = True non_local_groups += 1 if debug >= 2: print(f" GET_PROTONATION_STATES: will be sent to nonlocal due to {a.label} atom") - elif a.connec > 1: - block[idx] = 1 + elif a.connec >= 1: + # block[idx] = 1 + elemlist[idx] = "H" + addedlist[idx] = 1 + # Sulfur and Selenium elif a.label == "S" or a.label == "Se": if a.connec == 1: elemlist[idx] = "H" addedlist[idx] = 1 elif a.connec > 1: - block[idx] = 1 + # block[idx] = 1 + elemlist[idx] = "H" + addedlist[idx] = 1 # Hydrides elif a.label == "H": if a.connec == 0: @@ -390,7 +395,9 @@ def get_protonation_states_specie(specie: object, debug: int=0) -> list: addedlist[idx] = 1 else: # nitrogen with at least 3 adjacencies doesnt need H - if a.connec >= 3: block[idx] = 1 + if a.connec >= 3: + elemlist[idx] = "H" + addedlist[idx] = 1 else: # Checks for adjacent Atoms list_of_adj_atoms = [] diff --git a/cell2mol/classes.py b/cell2mol/classes.py index 8595b885..c9d79d57 100644 --- a/cell2mol/classes.py +++ b/cell2mol/classes.py @@ -1232,7 +1232,7 @@ def get_unique_species(self, debug: int=0): if debug >= 0: print(f"Getting unique species in cell") self.unique_species = [] self.unique_indices = [] - + self.species_list = [] typelist_mols = [] # temporary variable typelist_ligs = [] # temporary variable typelist_mets = [] # temporary variable @@ -1256,7 +1256,7 @@ def get_unique_species(self, debug: int=0): if debug >= 2: print(f"New molecule found with: formula={mol.formula} and added in position {kdx}") self.unique_indices.append(kdx) mol.unique_index = kdx - + self.species_list.append(mol) else: if not hasattr(mol,"ligands"): mol.split_complex(debug=debug) for jdx, lig in enumerate(mol.ligands): # ligands @@ -1273,7 +1273,7 @@ def get_unique_species(self, debug: int=0): if debug >= 2: print(f"New ligand found with: formula {lig.formula} added in position {kdx}") self.unique_indices.append(kdx) lig.unique_index = kdx - + self.species_list.append(lig) for jdx, met in enumerate(mol.metals): # metals found = False for ldx, typ in enumerate(typelist_mets): @@ -1288,7 +1288,7 @@ def get_unique_species(self, debug: int=0): if debug >= 2: print(f"New Metal Center found with: labels {met.label} and added in position {kdx}") self.unique_indices.append(kdx) met.unique_index = kdx - + self.species_list.append(met) return self.unique_species ####################################################### @@ -1540,25 +1540,26 @@ def get_selected_cs(self, debug: int=0): tmp = specie.get_possible_cs(debug=debug) if tmp is None: self.error_empty_poscharges = True - return # Stopping. Empty list of possible charges received. + return None # Stopping. Empty list of possible charges received. if specie.subtype != "metal": selected_cs.append(list([cs.corr_total_charge for cs in specie.possible_cs])) else : selected_cs.append(specie.possible_cs) self.error_empty_poscharges = False - return selected_cs + return selected_cs ####################################################### def assign_charges (self, debug: int=0): - if not hasattr(self,"unique_species"): self.get_unique_species(debug=debug) - if debug >= 1: print(f"{len(self.unique_species)} Species (Metal or Ligand or Molecules) to Characterize") + # if not hasattr(self,"unique_species"): self.get_unique_species(debug=debug) + # if debug >= 1: print(f"{len(self.unique_species)} Species (Metal or Ligand or Molecules) to Characterize") - selected_cs = self.get_selected_cs(debug=debug) - if debug >= 1: print(f"Selected Charges: {selected_cs}") + # selected_cs = self.get_selected_cs(debug=debug) + # if debug >= 1: print(f"Selected Charges: {selected_cs=}") final_charge_distribution, final_charges = balance_charge(self.unique_indices, self.unique_species, debug=debug) - + print("final_charge_distribution", final_charge_distribution) + print("final_charges", final_charges) if len(final_charge_distribution) > 1: if debug >= 1: print("More than one Possible Distribution Found:", final_charge_distribution) self.error_multiple_distrib = True @@ -1602,14 +1603,13 @@ def assign_charges (self, debug: int=0): print(specie.formula, specie.charge_state, specie.totcharge, specie.smiles) for idx, ref in enumerate(self.refmoleclist): - print(f"Molecule {idx}: {ref.formula}") + print(f"Refenrence Molecule {idx}: {ref.formula}") if ref.iscomplex: for jdx, lig in enumerate(ref.ligands): for specie in self.unique_species: if specie.subtype == "ligand": issame = compare_species(lig, specie) if issame: - set_charge_state_simple(specie, lig, debug=debug) set_charge_state (specie, lig, mode=1, debug=debug) print(lig.formula, specie.formula, issame) for met in ref.metals: @@ -1625,18 +1625,76 @@ def assign_charges (self, debug: int=0): if specie.subtype == "molecule": issame = compare_species(ref, specie) if issame: - set_charge_state_simple(specie, ref, debug=debug) set_charge_state (specie, ref, mode=1, debug=debug) print(ref.formula, specie.formula, issame) + for idx, mol in enumerate(self.moleclist): + print(f"Unit cell Molecule {idx}: {mol.formula}") + if not mol.iscomplex: + for ref in self.refmoleclist: + if not ref.iscomplex : + issame = compare_molecules(ref, mol, debug=debug) + if issame: + set_charge_state (ref, mol, mode=2, debug=debug) + print(mol.formula, ref.formula, issame) + else: + for ref in self.refmoleclist: + if ref.iscomplex: + for lig in mol.ligands: + for ref_lig in ref.ligands: + issame = compare_molecules(ref_lig, lig, debug=debug) + if issame: + set_charge_state (ref_lig, lig, mode=2, debug=debug) + print(lig.formula, ref_lig.formula, issame) + for met in mol.metals: + for ref_met in ref.metals: + if ref_met.get_parent_index("reference") == met.get_parent_index("reference"): + met.set_charge(ref_met.charge) + prepare_mol(mol) + + def assign_charges_for_refcell(self, debug: int=0): + for idx, ref in enumerate(self.refmoleclist): + print(f"Refenrence Molecule {idx}: {ref.formula}") + if ref.iscomplex: + for jdx, lig in enumerate(ref.ligands): + for specie in self.unique_species: + if specie.subtype == "ligand": + issame = compare_species(lig, specie) + if issame: + set_charge_state (specie, lig, mode=1, debug=debug) + print(lig.formula, specie.formula, issame) + for met in ref.metals: + for specie in self.unique_species: + if specie.subtype == "metal": + issame = compare_metals(met, specie) + if issame: + met.set_charge(specie.charge) + print(met.formula, specie.formula, issame) + prepare_mol(ref) + else: + for specie in self.unique_species: + if specie.subtype == "molecule": + issame = compare_species(ref, specie) + if issame: + set_charge_state (specie, ref, mode=1, debug=debug) + print(ref.formula, specie.formula, issame) + ####################################################### def check_charge_neutrality(self, debug: int=0): if self.subtype == "reference": moleclist = self.refmoleclist else: moleclist = self.moleclist - totcharge_list = [mol.totcharge for mol in moleclist] - if debug >= 1: print(f"Total Charge of the Cell ({self.subtype}): {sum(totcharge_list)} {totcharge_list=}") - if sum(totcharge_list) == 0: self.is_neutral = True - else: self.is_neutral = False + totcharge_list =[] + for mol in moleclist: + if not hasattr(mol,"totcharge"): + if debug >= 1: print(f"CELL.CHECK_CHARGE_NEUTRALITY: Charges not assigned yet") + self.is_neutral = None + else: + totcharge_list.append(mol.totcharge) + + if len(totcharge_list) !=0 : + if debug >= 1: print(f"Total Charge of the Cell ({self.subtype}): {sum(totcharge_list)} {totcharge_list=}") + if sum(totcharge_list) == 0: self.is_neutral = True + else: self.is_neutral = False ####################################################### def assign_charges_old (self, debug: int=0) -> object: @@ -1821,17 +1879,18 @@ def assess_errors(self, mode): if self.has_isolated_H: case = 1 elif self.has_missing_H: case = 2 else : case = 0 - elif mode == "reference": + elif mode == "unit_cell": print("-------------------------------") print("Errors in Reference Molecules") print("-------------------------------") if self.has_isolated_H: case = 1 elif self.has_missing_H: case = 2 + elif self.error_reconstruction: case = 4 elif self.error_empty_poscharges : case = 5 elif self.error_multiple_distrib : case = 6 elif self.error_empty_distrib : case = 7 else : case = 0 - elif mode == "unit_cell": + elif mode == "neutrality": print("-------------------------------") print("Errors in Unit Cell") print("-------------------------------") diff --git a/cell2mol/new_c2m_driver.py b/cell2mol/new_c2m_driver.py index 9e17f912..0899848c 100644 --- a/cell2mol/new_c2m_driver.py +++ b/cell2mol/new_c2m_driver.py @@ -8,7 +8,7 @@ from cell2mol.cell_operations import frac2cart_fromparam from ase.io import read from cell2mol.charge_assignment import balance_charge -from cell2mol.new_charge_assignment import print_output, compare_molecules, prepare_mol +from cell2mol.new_charge_assignment import print_output, compare_molecules, prepare_mol, get_unique_indices, assign_charge_state_for_unique_species from ase import Atoms from cell2mol.new_cell_reconstruction import * import copy @@ -94,94 +94,72 @@ if not refcell.has_isolated_H: refcell.check_missing_H(debug=debug) refcell.assess_errors(mode="hydrogens") - - refcell.assign_charges(debug=debug) - refcell.assess_errors(mode="reference") - - if refcell.error_case == 0: - refcell.check_charge_neutrality(debug=debug) - refcell.assign_spin(debug=debug) - refcell.create_bonds(debug=debug) - - refcell.save(ref_cell_fname) - - - #################################### - ### RECONSTRUCTS THE CELL OBJECT ### - #################################### if refcell.error_case == 0: + # Define new cell object for the unit cell newcell = cell(name, cell_labels, cell_pos, cell_fracs, cell_vector, cell_param) newcell.get_subtype("unit_cell") - newcell.get_reference_molecules(ref_labels, ref_fracs, debug=debug) + + # Get reference molecules + newcell.get_reference_molecules(refcell.labels, refcell.frac_coord, debug=debug) if not newcell.has_isolated_H: newcell.check_missing_H(debug=debug) newcell.assess_errors(mode="hydrogens") - # newcell.refmoleclist = copy.deepcopy(refcell.refmoleclist) + + # if newcell.error_case == 0: # Omit this condition + #since newcell.error_case with mode="hydrogens" is same as refcell.error_case with mode="hydrogens" - ref_molecule = Atoms(symbols=refcell.labels, scaled_positions=refcell.frac_coord, cell=cell_vector, pbc=True) - all_molecules, reconstructed_molecules = reconstuct(ref_molecule, newcell, refcell, cell_pos, cell_fracs, cell_vector, sym_ops, debug=0) - - newcell.moleclist = [] + # Reconstruction of the unit cell + reference = Atoms(symbols=refcell.labels, scaled_positions=refcell.frac_coord, cell=cell_vector, pbc=True) + all_molecules, reconstructed_molecules = reconstuct(reference, newcell, cell_pos, cell_fracs, cell_vector, sym_ops, debug=0) all_molecules.extend(reconstructed_molecules) - for mol in all_molecules: - newmolec = molecule(mol.labels, mol.coord, mol.frac_coord) - newmolec.origin = "cell.reconstruct" - newmolec.set_atoms(create_adjacencies=True, debug=debug) - newmolec.add_parent(newcell, mol.cell_indices) - newmolec.add_parent(refcell, mol.ref_indices) - for atom, idx in zip(newmolec.atoms, mol.cell_indices): - atom.add_parent(newcell, index=idx) - for atom, idx in zip(newmolec.atoms, mol.ref_indices): - atom.add_parent(refcell, index=idx) - if newmolec.iscomplex: newmolec.split_complex() - newcell.moleclist.append(newmolec) - - for mol in newcell.moleclist: - if mol.iscomplex: - mol.get_hapticity(debug=debug) - for lig in mol.ligands: - lig.get_denticity(debug=debug) - for met in mol.metals: - met.get_coordination_geometry(debug=debug) - - for mol in newcell.moleclist: - if not mol.iscomplex: - for ref in refcell.refmoleclist: - if not ref.iscomplex: - compare_molecules(ref, mol, debug=debug) + if not newcell.error_reconstruction : + # Get moleclist for the unit cell + newcell = get_moleclist(newcell, refcell, all_molecules, debug=debug) + refcell.get_unique_species(debug=debug) + print("refcell.unique_species", [specie.formula for specie in refcell.unique_species], refcell.unique_indices) + selected_cs = refcell.get_selected_cs(debug=debug) + if selected_cs is None: + newcell.error_empty_poscharges = True else: - for ref in refcell.refmoleclist: - if ref.iscomplex: - for lig in mol.ligands: - for ref_lig in ref.ligands: - compare_molecules(ref_lig, lig, debug=debug) - for met in mol.metals: - for ref_met in ref.metals: - if ref_met.get_parent_index("reference") == met.get_parent_index("reference"): - met.set_charge(ref_met.charge) - prepare_mol(mol) - - newcell.check_charge_neutrality(debug=debug) - newcell.assess_errors(mode="unit_cell") - newcell.assign_spin(debug=debug) - newcell.create_bonds(debug=debug) + newcell.error_empty_poscharges = False + print("selected_cs for reference cell") + for specie, select in zip(refcell.unique_species, selected_cs): + print(specie.possible_cs, select) + for specie, idx in zip(refcell.species_list, refcell.unique_indices): + print(specie.formula, specie.unique_index, idx) + + newcell = get_unique_indices(newcell, refcell.species_list, debug=debug) + print("newcell.unique_indices", newcell.unique_indices) + for specie, idx in zip(newcell.species_list, newcell.unique_indices): + print(specie.formula, specie.unique_index, idx) + + newcell.unique_species = copy.deepcopy(refcell.unique_species) + + final_charge_distribution, final_charges = balance_charge(newcell.unique_indices, refcell.unique_species, debug=debug) + # print("final_charge_distribution", final_charge_distribution) + # print("final_charges", final_charges) + newcell.assign_charges(debug=debug) + newcell.assess_errors(mode="unit_cell") + newcell.check_charge_neutrality(debug=debug) + if newcell.error_case == 0 and newcell.is_neutral: + newcell.assign_spin(debug=debug) + newcell.create_bonds(debug=debug) + + refcell.unique_species = assign_charge_state_for_unique_species(refcell.unique_species, final_charges[0], debug=debug) + refcell.assign_charges_for_refcell(debug=debug) + refcell.assign_spin(debug=debug) + refcell.create_bonds(debug=debug) + # Save cell object newcell.save(cell_fname) + + # Save reference cell object + refcell.save(ref_cell_fname) + output.close() sys.stdout = stdout - surmmary = open(surmmary_fname, "w") - sys.stdout = surmmary - print("*** Reference molecules ***") - print(refcell) - print_output(refcell.refmoleclist) - print("***Unit cell molecules ***") - print(newcell) - print_output(newcell.moleclist) - surmmary.close() - sys.stdout = stdout - # Error handling case = refcell.error_case error_fname = os.path.join(current_dir, f"refcell_error_{case}.out") @@ -190,12 +168,32 @@ handle_error(case) error.close() sys.stdout = stdout + + # Summary + surmmary = open(surmmary_fname, "w") + sys.stdout = surmmary + print("*** Reference molecules ***") + print(refcell) + print_output(refcell.refmoleclist) + surmmary.close() + sys.stdout = stdout - # Error handling - case = newcell.error_case - error_fname = os.path.join(current_dir, f"unitcell_error_{case}.out") - error = open(error_fname, "w") - sys.stdout = error - handle_error(case) - error.close() - sys.stdout = stdout \ No newline at end of file + + if newcell.error_case == 0: + # Error handling + case = newcell.error_case + error_fname = os.path.join(current_dir, f"unitcell_error_{case}.out") + error = open(error_fname, "w") + sys.stdout = error + handle_error(case) + error.close() + sys.stdout = stdout + + # Summary + surmmary = open(surmmary_fname, "a") + sys.stdout = surmmary + print("***Unit cell molecules ***") + print(newcell) + print_output(newcell.moleclist) + surmmary.close() + sys.stdout = stdout \ No newline at end of file diff --git a/cell2mol/new_cell_reconstruction.py b/cell2mol/new_cell_reconstruction.py index 2cce2667..a45e1c8d 100644 --- a/cell2mol/new_cell_reconstruction.py +++ b/cell2mol/new_cell_reconstruction.py @@ -62,7 +62,7 @@ def find_row_index_from_matrix (matrix, query_row): return -1 ###################################################### -def get_fragments (newcell, updated, indices_in_ref, refcell, cov_factor: float=1.3, metal_factor: float=1.0, debug: int=0): +def get_fragments (newcell, updated, indices_in_ref, cov_factor: float=1.3, metal_factor: float=1.0, debug: int=0): updated_labels = extract_from_list(updated, newcell.labels, dimension=1) updated_coord = extract_from_list(updated, newcell.coord, dimension=1) @@ -89,7 +89,6 @@ def get_fragments (newcell, updated, indices_in_ref, refcell, cov_factor: float= # Adds cell as parent of the molecule, with indices newmolec.add_parent(newcell, indices=cell_indices) - newmolec.add_parent(refcell, indices=ref_indices) newmolec.set_fractional_coord(mol_frac_coord) newmolec.set_adjacency_parameters(cov_factor=cov_factor, metal_factor=metal_factor) newmolec.set_atoms(create_adjacencies=True, debug=debug) @@ -267,7 +266,7 @@ def merge_fragments (frags: list, cell_vector: list, cov_factor: float=1.3, meta return None ###################################################### -def fragments_reconstruct (subset_remaining_fragments, target_ref, cell_vector, newcell, refcell, debug: int=0): +def fragments_reconstruct (subset_remaining_fragments, target_ref, cell_vector, debug: int=0): list_of_found_molecules = [] final_remaining = subset_remaining_fragments.copy() @@ -332,7 +331,7 @@ def fragments_reconstruct (subset_remaining_fragments, target_ref, cell_vector, return list_of_found_molecules, final_remaining ###################################################### -def reconstuct (reference, newcell, refcell, cell_pos, cell_fracs, cell_vector, sym_ops, debug: int=0): +def reconstuct (reference, newcell, cell_pos, cell_fracs, cell_vector, sym_ops, debug: int=0): new_structures = apply_symmetry_operations_reference(reference, cell_vector, sym_ops) @@ -359,7 +358,7 @@ def reconstuct (reference, newcell, refcell, cell_pos, cell_fracs, cell_vector, if len(updated) > 0 : # #### make blocks and get fragments #### - fragments = get_fragments (newcell, updated, indices_in_ref, refcell, debug=0) + fragments = get_fragments (newcell, updated, indices_in_ref, debug=0) molecules, remaining_fragments = classify_fragments(fragments, newcell, debug=0) all_molecules.extend(molecules) @@ -380,13 +379,43 @@ def reconstuct (reference, newcell, refcell, cell_pos, cell_fracs, cell_vector, print(f"Group {i}: {[rem.formula for rem in group]}") target_ref = newcell.refmoleclist[i].get_parent_indices("reference") - list_of_found_molecules, final_remaining = fragments_reconstruct(group, target_ref, cell_vector, newcell, refcell, debug=0) + list_of_found_molecules, final_remaining = fragments_reconstruct(group, target_ref, cell_vector, debug=0) print(f"{list_of_found_molecules=}") reconstructed_molecules.extend(list_of_found_molecules) if len(all_found) == len(cell_pos): print("Reconstructed successfully") + newcell.error_reconstruction = False else: print("Error in reconstruction!!") + newcell.error_reconstruction = True + + return all_molecules, reconstructed_molecules + +###################################################### +def get_moleclist (newcell, refcell, all_molecules, debug): + # Get moleclist for the unit cell + newcell.moleclist = [] + + for mol in all_molecules: + newmolec = molecule(mol.labels, mol.coord, mol.frac_coord) + newmolec.origin = "cell.reconstruct" + newmolec.set_atoms(create_adjacencies=True, debug=debug) + newmolec.add_parent(newcell, mol.cell_indices) + newmolec.add_parent(refcell, mol.ref_indices) + for atom, idx in zip(newmolec.atoms, mol.cell_indices): + atom.add_parent(newcell, index=idx) + for atom, idx in zip(newmolec.atoms, mol.ref_indices): + atom.add_parent(refcell, index=idx) + if newmolec.iscomplex: newmolec.split_complex() + newcell.moleclist.append(newmolec) + + for mol in newcell.moleclist: + if mol.iscomplex: + mol.get_hapticity(debug=debug) + for lig in mol.ligands: + lig.get_denticity(debug=debug) + for met in mol.metals: + met.get_coordination_geometry(debug=debug) - return all_molecules, reconstructed_molecules \ No newline at end of file + return newcell \ No newline at end of file diff --git a/cell2mol/new_charge_assignment.py b/cell2mol/new_charge_assignment.py index 3b645109..9cf5933a 100644 --- a/cell2mol/new_charge_assignment.py +++ b/cell2mol/new_charge_assignment.py @@ -7,15 +7,106 @@ from rdkit import Chem ###################################################### +def assign_charge_state_for_unique_species(unique_species, final_charges_tuple, debug: int=0): + + for specie, final_charge in zip(unique_species, final_charges_tuple): + print(specie.unique_index, specie.formula) + if (specie.subtype == "molecule" and specie.iscomplex == False) or (specie.subtype == "ligand"): + charge_list = [cs.corr_total_charge for cs in specie.possible_cs] + idx = charge_list.index(final_charge) + cs = specie.possible_cs[idx] + specie.charge_state = cs + # print(specie.charge_state.protonation) + specie.set_charges(cs.corr_total_charge, cs.corr_atom_charges, cs.smiles, cs.rdkit_obj) + elif specie.subtype == "metal" : + charge_list = specie.possible_cs + idx = charge_list.index(final_charge) + cs = specie.possible_cs[idx] + specie.set_charge(cs) + for specie in unique_species: + if (specie.subtype == "molecule" and specie.iscomplex == False) or (specie.subtype == "ligand"): + print(specie.formula, specie.charge_state, specie.totcharge, specie.smiles) + return unique_species +###################################################### +def get_unique_indices(newcell, reference_species_list, debug: int=0): + + newcell.unique_indices = [] + newcell.species_list = [] + for mol in newcell.moleclist: + if not mol.iscomplex: + for ref in reference_species_list: + if (ref.subtype == "molecule") and not ref.iscomplex: + issame = compare_molecules(ref, mol, debug=debug) + if issame: + mol.unique_index = ref.unique_index + if debug >= 1: print(f"Matched {mol.formula} {ref.formula} {mol.unique_index} {ref.unique_index}") + newcell.unique_indices.append(mol.unique_index) + newcell.species_list.append(mol) + else: + for ref in reference_species_list: + if ref.subtype == "ligand": + for lig in mol.ligands: + issame = compare_molecules(ref, lig, debug=debug) + if issame: + lig.unique_index = ref.unique_index + if debug >= 1: print(f"Matched {lig.formula} {ref.formula} {lig.unique_index} {ref.unique_index}") + newcell.unique_indices.append(lig.unique_index) + newcell.species_list.append(lig) + if ref.subtype == "metal": + for met in mol.metals: + if ref.get_parent_index("reference") == met.get_parent_index("reference"): + met.unique_index = ref.unique_index + if debug >= 1: print(f"Matched {met.formula} {ref.formula} {met.unique_index} {ref.unique_index}") + newcell.unique_indices.append(met.unique_index) + newcell.species_list.append(met) + return newcell +###################################################### +def get_unique_indices_old (newcell, refcell, debug: int=0): + + newcell.unique_indices = [] + for mol in newcell.moleclist: + if not mol.iscomplex: + for ref in newcell.refmoleclist: + if not ref.iscomplex: + issame = compare_molecules(ref, mol, debug=debug) + if issame: + mol.unique_index = ref.unique_index + if debug >= 1: print(f"Matched {mol.formula} {ref.formula} {ref.unique_index}") + newcell.unique_indices.append(mol.unique_index) + else: + for ref in newcell.refmoleclist: + if ref.iscomplex: + for lig in mol.ligands: + for ref_lig in ref.ligands: + issame = compare_molecules(ref_lig, lig, debug=debug) + if issame: + lig.unique_index = ref_lig.unique_index + if debug >= 1: print(f"Matched {lig.formula} {ref_lig.formula} {ref_lig.unique_index}") + newcell.unique_indices.append(lig.unique_index) + + for met in mol.metals: + for ref_met in ref.metals: + if ref_met.get_parent_index("reference") == met.get_parent_index("reference"): + met.unique_index = ref_met.unique_index + if debug >= 1: print(f"Matched {met.formula} {ref_met.formula} {ref_met.unique_index}") + newcell.unique_indices.append(met.unique_index) + return newcell ###################################################### def compare_molecules(ref, mol, debug: int=0): if (ref.natoms == mol.natoms) & (ref.formula == mol.formula): if (sorted(ref.get_parent_indices("reference")) == sorted(mol.get_parent_indices("reference"))): print("Matched", mol.formula, ref.formula, ref.get_parent_indices("reference"), mol.get_parent_indices("reference")) + issame = True + # mol.unique_index = ref.unique_index # set_charge_state_simple(ref, mol, debug=debug) - set_charge_state(ref, mol, mode=2, debug=debug) + # set_charge_state(ref, mol, mode=2, debug=debug) + else: + issame = False + else : + issame = False + return issame ###################################################### def set_charge_state_simple (reference, target, debug: int=0): @@ -161,12 +252,19 @@ def prepare_mol (mol): ###################################################### def print_output(moleclist): + for idx, mol in enumerate(moleclist): if mol.iscomplex: - print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula} {mol.is_haptic=} {mol.totcharge=} {mol.spin=}") #\n {mol.adjnum=}\n {mol.madjnum=} \n {mol.smiles=}") + if hasattr(mol, "totcharge") and hasattr(mol, "spin"): + print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula} {mol.is_haptic=} {mol.totcharge=} {mol.spin=}") #\n {mol.adjnum=}\n {mol.madjnum=} \n {mol.smiles=}") + else: + print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula} {mol.is_haptic=}") # {mol.totcharge=} {mol.spin=}") #\n {mol.adjnum=}\n {mol.madjnum=} \n {mol.smiles=}") # print(mol.adjnum) for lig in mol.ligands: - print(f"|- {lig.subtype}({lig.type}) {lig.formula} {lig.is_haptic=} {lig.denticity=} {lig.totcharge=} {lig.smiles=}") + if hasattr(lig, "totcharge") and hasattr(lig, "smiles"): + print(f"|- {lig.subtype}({lig.type}) {lig.formula} {lig.is_haptic=} {lig.denticity=} {lig.totcharge=} {lig.smiles=}") + else : + print(f"|- {lig.subtype}({lig.type}) {lig.formula} {lig.is_haptic=} {lig.denticity=}")# {lig.totcharge=} {lig.smiles=}") # print(f"|- {lig.connected_idx}") # print(lig.groups) for group in lig.groups: @@ -175,7 +273,10 @@ def print_output(moleclist): print(f"|--- {met.label} {met.mconnec=}") print("") for metal in mol.metals: - print(f"|# {metal.subtype}({metal.type}) {metal.label} {metal.coord_nr=} {metal.coord_geometry} {metal.charge=} {metal.spin=} {metal.get_coord_sphere_formula()} {metal.mconnec=} {metal.connec=}") + if hasattr(metal, "charge") and hasattr(metal, "spin"): + print(f"|# {metal.subtype}({metal.type}) {metal.label} {metal.coord_nr=} {metal.coord_geometry} {metal.get_coord_sphere_formula()} {metal.mconnec=} {metal.connec=} {metal.charge=} {metal.spin=}") + else : + print(f"|# {metal.subtype}({metal.type}) {metal.label} {metal.coord_nr=} {metal.coord_geometry} {metal.get_coord_sphere_formula()} {metal.mconnec=} {metal.connec=}") #{metal.charge=} {metal.spin=} # print(f"|# {metal.get_coord_sphere_formula()}") # print(f"|# {metal.coord_sphere_formula}") # print(f"|# {metal.mconnec=} {metal.connec=}") @@ -183,5 +284,8 @@ def print_output(moleclist): # for bond in metal.bonds: # print(f"|--- {bond}") else: - print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula} {mol.totcharge=} {mol.spin=}\n {mol.smiles}") + if hasattr(mol, "totcharge") and hasattr(mol, "spin"): + print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula} {mol.totcharge=} {mol.spin=}\n {mol.smiles}") + else : + print(f"{idx}: {mol.subtype}({mol.type}) {mol.formula}") #{mol.totcharge=} {mol.spin=}\n {mol.smiles}") print("") \ No newline at end of file diff --git a/cell2mol/test/check_Cell_object.ipynb b/cell2mol/test/check_Cell_object.ipynb index aaf2dbea..ddfa5e1f 100644 --- a/cell2mol/test/check_Cell_object.ipynb +++ b/cell2mol/test/check_Cell_object.ipynb @@ -2,16 +2,36 @@ "cells": [ { "cell_type": "code", - "execution_count": 2, + "execution_count": 30, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "text/plain": [ + "'2023.09.6'" + ] + }, + "execution_count": 30, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ - "import numpy as np" + "import numpy as np\n", + "from rdkit import Chem\n", + "from rdkit.Chem import Draw\n", + "from IPython.display import display\n", + "from rdkit.Chem.Draw import IPythonConsole\n", + "from rdkit.Chem import rdDetermineBonds\n", + "\n", + "IPythonConsole.ipython_3d = False\n", + "import rdkit\n", + "rdkit.__version__" ] }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 14, "metadata": {}, "outputs": [ { @@ -28,13 +48,206 @@ }, { "cell_type": "code", - "execution_count": 28, + "execution_count": 42, "metadata": {}, "outputs": [], "source": [ - "name = \"XAWJUB\"\n", + "name = \"ZURVIR\"\n", "refcell = np.load(f\"{name}/Ref_Cell_{name}.cell\", allow_pickle=True)\n", - "cell = np.load(f\"{name}/Cell_{name}.cell\", allow_pickle=True)" + "# cell = np.load(f\"{name}/Cell_{name}.cell\", allow_pickle=True)" + ] + }, + { + "cell_type": "code", + "execution_count": 45, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "------------- Cell2mol CELL Object ----------------\n", + " Version = 2.0\n", + " Type = cell\n", + " Sub-Type = reference\n", + " Name (Refcode) = ZURVIR\n", + " Num Atoms = 167\n", + " Cell Parameters a:c = [20.504, 12.347, 22.524]\n", + " Cell Parameters al:ga = [90.0, 94.56, 90.0]\n", + "---------------------------------------------------\n", + " # of Ref Molecules: = 7\n", + " With Formulae: \n", + " 0: H8-C4-O \n", + " 1: H8-C4-O \n", + " 2: H8-C4-O \n", + " 3: H8-C4-O \n", + " 4: H64-C44-N4-Cr \n", + " 5: Li \n", + " 6: Li " + ] + }, + "execution_count": 45, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "refcell" + ] + }, + { + "cell_type": "code", + "execution_count": 43, + "metadata": {}, + "outputs": [], + "source": [ + "def generate_xyzblock(specie):\n", + " xyzblock = \"\"\n", + " xyzblock += str(specie.natoms) + \"\\n\" # Number of atoms\n", + " xyzblock += \"\\n\"\n", + " for l, pos in zip(specie.labels, specie.coord):\n", + " xyzblock += f\"{l}\\t{pos[0]}\\t{pos[1]}\\t{pos[2]}\\n\"\n", + " return xyzblock" + ] + }, + { + "cell_type": "code", + "execution_count": 66, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "112" + ] + }, + "execution_count": 66, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "len(['N', 'N', 'N', 'N', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'C', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H'])" + ] + }, + { + "cell_type": "code", + "execution_count": 65, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "112\n", + "112 0\n" + ] + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcIAAACWCAIAAADCEh9HAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dd1xT1/s48CcTEgJJ2MiSLUumIFZUULRibRVb/VitW5x1oVVrVVx1Yh1V6yragluLu7gVAQFZsmWEHXYIgSwyfn+c9n75oZ2oVDjvv+653CQneb18PPfcc56HpFarAcMwDPu3yN3dAQzDsPcbDqMYhmFdgsMohmFYl+AwimEY1iU4jGIYhnUJDqMYhmFdgsMo9m9cvHgRHcTGxopEokuXLqEmcYBhvQcOo9i/0TGMtrS0EE3i4L+vvr6+vr4eANra2srKyoqLi2UyGQDw+XyBQNDdvcPeJ9Tu7gD2Xmpvbz948CAAZGRkAEBlZSVq1tTUdHPP/rbHjx+r1erPPvuspKQkJiZGKBQuXbrU3Nz8xo0bffv2DQoK6u4OYu8NHEaxf4NKpY4dOxYAMjMzAUBfXx8179+/3809+yfKy8tfvHhRVFSEmvn5+QKBoKqqqm/fvt3aL+w9g8Mo9m+QSCQUa3R0dABAU1MTNUkkUkFBgYODg0KhSE9PZ7PZurq6+vr6NTU1YrHY2tq6W3vdGY/HS09Pr6ysRM3s7Ozq6ury8vIPPvigezuGvV8o4eHh3d0H7P1jaGiI4qaurq6FhYWpqSlqNjc35+TkDBo0qLW1dfv27SKRiEajWVhYJCQkZGZmenp6dm+3O8rNzbW0tPz888+5XG5+fr5MJps1a9aQIUPq6uo4HI6NjU13dxB7b+BHTNi/MXToUHQwYMAAJpNJND09PQsKCmJjY4m7++Tk5NjY2NTU1O7pKAAxUDh69CifzyeasbGxVlZWAKCrq+vu7u7j46OlpQUADg4Opqam3dNX7P2Eb+qxHi4vLw8d8Hg8qVRKNIVCobe3NwCYmpp2jJtDhgx5953E3ms4jGJvhkgk+vLLL8VicXl5+YIFC27dugUACoXiyZMnK1asIJFIFRUV3dIxpVIZFhYGAHFxcfPmzSspKUHN4uLibukP1vPgm3rszYiMjPz000/37t27a9cuEomkVCqDg4M9PT2pVGp+fn5OTo6bm5tYLN6xY8fly5dfvHgBACkpKTdu3HhTHYiNjUUjzaysrHv37h08eFChUAAAn89funTpJ5980t7evn79ei6Xq62tHRERYWNjs2vXLolE8qY6gPVaOIxib0ZJSYmTk5OZmRm6Kb579+7Tp09jYmK0tLRaW1sVCoW3t7dSqSwqKuLz+UKhEACam5vf4DrToqIitJy+tra2pKQkOTlZpVIBgEgkqqur2759e1BQUHh4uJaWVm5uLnrJy5cvUajFsK7AN/XYm2Fra5uVlUUsaRo9evTXX38tlUrnzZsHADdu3KipqWlvb0d/PXToUExMTHl5+ahRo95gHw4fPnz16tXy8nK0eH716tVkMpnJZNbW1vbr12/dunUsFuvEiRNTpkwJCwszMDB48uTJG/x0rNfCYRR7M2bOnLl48eIHDx7I5fKFCxcymUwAIJFILBYLAFxdXbW1tcPCwpYvXw4AixYt8vf3v3v3bllZ2Rvsw8KFC4cMGXLv3r2SkhIA2LlzJ51O37ZtG4fD4fP5bDYbAPT09HR0dCIiIgBgzpw5b/DTsV4Lh1HszdDS0oqMjCSarq6uAKChoXHo0KGEhITKykp7e/vu6pu2tratre2aNWu0tbUDAgJevWDr1q0LFy7U1dV9933DegASLmmHvW2JiYnp6enPnz8/cuSIRCKh0Wg0Go1Op8tkMqVSicatXVdfX89gMFgslkgkkslkIpGob9++JBKppqaGzWYzGIy6ujqpVGphYdHY2KinpwcATU1NXC73+vXrL168mDlzJl4uiv07OIxi70hhYaGdnV139+I1Vq5cCQDjxo0bPHhwd/cFey/hMIphkJ2dbWFhgfIDYNg/hcMohmFYl+B1oxiGYV2CwyiGYViX4DCKYRjWJTiMYp2JRCK001wmk7W0tDQ0NKAJdKFQiKoVYRjWEQ6jWGfnz5+Pi4sDgJSUlNOnT69Zs0YkEgHA8ePHUUoRDMM6wruYsNcQi8VCobCtrQ01W1pa1Gq1VCoFgPLy8p07d1IolIEDB4rF4g8//NDMzOzOnTva2tp+fn7d2msM6x54NIq9xi+//LJ//36iWvKRI0f279//9OlTAPjqq6/WrVt34MCB4cOH83g8dPtfW1vb2NjYnT3+txQKRXZ2Npq4qKqqQicrKiqkUimqdQoAu3fvzszMvHPnDgA0NDR03POKYYDDKPZaU6ZM2bBhw4wZM1Bz9erVGzZsGDFiBACIxeI+ffoAgJGREQDs2LEjLCzs7Nmz3dfZf08ikYSEhMTHx3/77bdVVVVr165F51euXKlQKFDtaAB4/vx5fX09yqIikUhycnK6rcfYfxK+qcc669OnD9pyzuFwzM3N5XI5lUrNysoikUg6OjpisfjBgwdSqRRtnRwzZoybm1tCQkJ39/rfuHHjxvjx42fOnAkAUqlUIpHw+XwAUCqVAJCcnIzy5Le0tADAxYsX8/PzRSIR3uyEdYLDKNZZcHAwOnBxcXFxcUHHeXl5lpaWDg4OH3300cqVK9lsdkFBAQDk5uaiBHTIlStXDA0N35fN6QKBwNTU9ODBg19++SUApKamHjt2LCMjA80Cs1gslE9v0qRJAPDZZ5/NnTu3oqJi//793dtt7L8Gh1Hs77p+/XppaSmfz//oo4+qqqoaGxsZDEZGRoafn9+ECRPIZHJZWVl9fT2KQe+FAQMGnDhxoqGhYf78+Wq1mslkTps2LTc3l06nA0BtbW1rayvKl4phfwLPjWIAAKWlpejWtb6+ns/nZ2ZmovOFhYVEtSJ/f/8pU6aMGTNGQ0ODw+GEhITMnj3b19cXAJhMpqam5vnz56lUKo/H665v8U95eHj4+vqmpaV5e3uHhYWhdV0lJSV8Pn/Dhg1kMnn+/PkAsHv37oEDB06YMAEAjI2N0Z0+hhHwaBQDALh48WJAQIC3t3diYqJIJLp27dr58+cBIDIycvbs2TY2NgDA4XBMTExQvSMAMDQ0NDEx6XhH/9VXX7W1tRUWFnbLV/h3pk2bdvPmTfRl16xZAwDW1tbnz5+Xy+UPHz708fEBAAsLC+J6Go1mYmLSXb3F/ptwGMV+k5ubq1ari4qK0CP4lJQUAECPXACAWBNqaWk5fvx4mUyGHkMNHz68422vlpaWu7v7u+76W9DY2KipqZmRkUHkeMawP4LDKPabyspKOp3O5/ONjIzUajUq497c3Iz+am5ujg7YbHbHESgaqHZUUFBw7dq1VatWKZXKBQsWDBkyxNTUNCAgIDk5OSMjIzQ09J18m39gy5Yt6GDhwoVGRkaoaWJicvr0aSMjIzabnZeXV1NTExAQ0NzcfPPmTTab3b9/fwsLi9TUVIVCgaY1sN4Mh1HsNyNHjvT29mYymSKRiEQi/e9//wOAf7H7U6FQoOlUtVrd2toqk8lQQdD29vZ39vSpoqJCJBI5OTm1trZmZmbS6XQbGxtdXd2ioiISidQp9BNFotDNO9EkDmpqagoKCgICAsRicUpKipGRkZmZmYWFRXl5uUwmw2EUw2EUAwBwd3c3NDQEgL59+0qlUlSgGAB8fX07jj3/ptjYWKFQSGQEP3bsWGxsbHV19TvbLZqfn19WVubk5NTU1BQTE6OjozN+/HhdXd20tDQymfzqCPovXb58ubCwsK2tTVNTEwD27t1rYGBQVFQ0efLkt9B97D1DCQ8P7+4+YN3PxsYGhUsjIyNTU1NPT090vl+/fv+o5Byfz798+TKVSt27d6+enl5SUpK1tbWxsfHGjRtNTEz4fP67GbuVlJTcuXNHJBKlp6e3tLRoaGiUlJSUlZUlJiZaWlo6OTn9o3crLS01NDRct26dp6dnfHw8i8WaNm3a9OnTNTU1lUolsbQW67XwgifsjeHz+XPmzDEzM6utra2qqjp58iSTyaTT6TExMQBApVI1NDTeWWdMTEwcHR1tbW1R09ra2tHREdf+xN4GHEZ7uJ9++ok4aG5uvnbtWqfzb9D169fnzZs3efLkmJgYDodTVVU1ceJEIyMjDoezadMmfX39efPmnTlz5vbt21988QV6ybJly7qe06SiogId1NXVyWQy1LSwsDAwMHBycnry5Mnt27cvXbqkp6eXlpZGfO7ff39PT89PPvkEAPT19ZcsWTJ16lRU4nTIkCHDhw/vYuexHgCH0R7u5s2bxIFIJHry5Emn829Qe3s7mjokqNVqNEPa3t6uUqkAQKlUKpVKuVxOvKTrRRVRhWQAOHToEI/HCwsLc3BwGDhw4MqVK5OSknx8fPbu3fvdd98xGAw6ne7m5gYdFnL9HWw2G60VpdPp1tbW5ubmWlpaAKCnp2dgYNDFzmM9AH7E1MNVVFRs3boVAEpLSwHg2bNnqImenr9Zw4cP37Rpk7+/f2Vlpa6urpmZ2ahRowDg5MmTALBjxw6ZTNbU1LRw4UKxWIz2Ar2pnCbV1dUA0NraCgDt7e1UKlVXV1cmk5WWlv7vf/8bOHAgANTV1WVnZ//www8AQGzNwrCuw6PRHq5Pnz6LFi1atGiRmZkZAHh4eKAmjUbrdOX9+/cXL168dOnS0tLSHTt2oJPEwd/Rr1+/xYsXb9y48fLlyyQS6cMPP0Tnx40bBwCenp6Ghobo0TaTyYyIiIiIiBg0aFDXvyMAnD179uzZs1lZWQDQ2NiImjU1Nebm5kVFRcRlgwcPRp/LYDDeyOd29N1336ENtZcvX87Jydm0aRM6v2XLFlzGvGfDo9H3XllZGYlEsrCwaG5urqyslMlkjo6OTCazsLBQW1ubQqFwuVwAoFKpAKChoYGanQgEggMHDly+fFkqlcpksvT0dHSemEz8mz744IMPPvgAHY8fPx4dfP755xs2bCgqKmKz2VQqVVNTk8PhoD/p6OiQyW/g//KOGe2Ibe/JyckTJkyYNWtWXl5eY2Pj+vXrO72qvLw8MjJy48aNXe8AALx8+RKN8SsrK01MTHJzc9F5tD2MRCK9kU/B/oNwGH3vJSUlkclktBr8+vXrTU1Ny5cvZzKZDx8+tLGxIeYNV65cqa+vT2RiJs4jeXl5fn5+VCqVxWKxWKzq6moUiXJzc8vLy588eTJ16lSVSrVr1y4PDw99fX0vL6/c3NzCwkL07OUvhYWFFRUVnT9//pNPPmEymSgDNABs376967+Al5cXOnBwcNDW1vb29ibOa2hoREdHi8ViTU1NMplM5Lg7c+YMAFy8ePG1/6m8Kj09XSaTDRgw4OXLl46OjgDA4/FMTEw6zQVv2LBBU1MzPT19wIABDQ0N6Dck/k/CeiocRnuCtLQ0BoOBZj8B4MGDB3p6ellZWTY2Nt7e3k+ePJFKpQEBAVlZWWhBaH5+vrOzc8d30NPTq6mpIZp9+vRBqTYrKipaW1vRxlC1Wv3ixQu0Sh8ABAIBSgj/d7DZbC8vLyLevVlfffUVOvj88887NokDYukrMfKlUCg5OTlCobCwsPCLL77482AaFhZmamqqr69fU1Nz9uxZlMfk2LFjc+fOtba27njl5s2b9fT0ULDW19dHvyGat8V6MBxGewI6nc5gMIhVmZqamuipNAAsXrzYzc3NwMDg9u3b0dHRKAScPn26UwhwcHAQiUQRERFqtZqY0yTExsbevXt3wIABqHn8+PE7d+68y11Jb4Ozs/PmzZtTU1O5XG5sbCyFQhkxYsTLly8fP37c3Nw8ffp0Q0PDixcvmpqaFhcXo4AIAGfOnMnOzgaAhoaGbu0+9h+Cw2hP4OLiEhgYqK+vf/36dQAYNGiQmZlZUVGRQqGorq4+dOgQuiwqKioxMRH+YHx04sSJoqIitVptb2+/d+9edPK7774TCoWmpqbW1tabN29G9Tbmzp07atSo+Pj41NTUd/QN3xo0QBYKhRQKBQAkEkl9fX1DQwOa5WxqamIymR1TWKnV6k65rwjbtm1Dc76zZs2i0+nEz37w4MGXL1/W1dUNGTJEIBDcu3dPS0urf//+ZmZmqamparWamIXA3lM4jL73nJyc0L2qkZGRv7+/WCxGxYI8PT3RsLTjxVVVVYcPH0b/2mfPnv3VV19dvXr1q6++UqlU8+fPP3bsGLqM2O1jampaX19fXFxMp9Nfvnz5Tr/Yu3X69OmEhIT6+vp+/foBwObNm1ks1osXL9auXSsQCOrq6vT09Nra2shkMvq/hPg1GhsbtbW16XS6rq4un89/+PChtbW1t7d3SUmJvr4+ulIikRQXFw8ZMqStre3Zs2doNZiZmRmPx1OpVDiMvu9wGH3vEXu6jYyMUKpQBP3jFIlEFRUVBgYGbW1tJBLp008/vXDhAsrN0draWlVVhZ5uo2xMr31/Jyene/fuZWVlOTs77927l8lkoof+Xl5enSZY32vTp0+fMGFCZmYm2piwYcMGU1PTo0ePAsDhw4d37twpl8snT548dOhQdL2Xl5e2tnZbW9uyZctWrlzp5uZWW1s7e/bs9evXp6en9+nT5+jRoyjr8/79+0NDQ69cuVJcXCwSidB/bPv27TMwMCgoKEDzudh7DYfRHu748eMHDx4Ui8VTpkxBk575+flyuTwiIqK0tHT//v2lpaUSieRPFjaicRYKH8bGxsR5TU3NTs+p32uFhYVisRgA6uvrxWJxfHz8xIkTAeDFixdDhw4l5kaJha6ffvopAOzZs+fjjz9Gm7Ju3rw5Z84cPz8/Pz8/9F9UbGwsAMhkMgAYN27cnDlzqqqqDhw4AABffvmlq6vrlStXuuGrYm8aDqM9nLGx8bZt29AxGp/269fvwoULAODm5jZ8+PDs7OyIiAilUknsc++FQkJCVq5c2dzc7Orqamtra2dnd/LkyYkTJ86ePXvOnDkymaxTmisej3f37t3Q0FB7e/uIiIgPPvhAoVC0t7fz+fxz586hVK0A0HGtKJlMplKpaAYWACgUCpVKJZPJaI8s9l7Du5h6HSIimJqa8vl89CyFRCKVlJQsXrw4NDS0vLx8+fLl6BrioGejUqkkEun27dsxMTFZWVlUKlUikVy5cuXatWvl5eWvXi+RSNBjuo8++sjExMTa2looFAYGBl65cqWioqKwsFChUJiamo4cOXLkyJEaGhqenp5jxoxJTU01MDCYP3/+lClT0DKJIUOGDBs27B1/WeyNw6PRXufUqVPo4NatW8TJ3NxcHx+fAwcOyGQyGo1GPMr/ozWParX6yZMnQqFwxIgRBQUF/fv3p1AohYWFBgYGxA6l946+vr6+vr62tjYAUCgUNNH8R9tGUT5TACCRSGQy+fjx43379qVQKLGxsUqlct68eWPGjEFXjhs3jsvlPnny5OHDh/v27UPZoRBiES72XsNhFAMAyMvLQ8tC0eLTvLw8tAOnsrLytdevWrXKxsbG1tb25s2bv/766/79+1ksVkxMDCov+i57/gYNGDCgT58+aGM+nU5Hu1oNDQ2VSiWfzzcxMUGznBoaGkqlMiAgYNu2bSqVasqUKQAwd+7cDz/8MCEhISUlZenSpQAQEhKC3nby5MmNjY3p6ekcDkcmk73LpKvYu4HDKAYAYGZmhp6HII6OjhEREbdv30Zh9Pbt2wMGDCgtLfX29lYqlQ8ePCgoKNizZw+6+Ndff3369CmDwSgpKQkICOieL/BPqKTSRzdvisXiEZ99lpGRgfI/OTo6oolLHx8fPT09Yn7T29t76tSpvr6+AoHA19dXrVZ/9tlnxcXFaAXu30Sj0caOHRsTE6NQKHAY7XlwGO3JpFLprl27NmzYAABff/11SEhIeXl5SEhIfX19VFRUx3lPPz+/n3/+eeXKlW1tbevWrUMnidv/EydOREdHOzk5eXt7t7e3R0VFEY9KkJaWFrlc/s4q1nWJXH5m5MipcXEA8PLWre9UKh8fn5CQkGPHjgUGBhoZGf3444/ff/+9paUlujwtLe3QoUMoP9alS5dUKpVSqbS0tPT39wcAMpl87NgxOp2Olu76+PgQ9Vc60tHReXt7YbFuh8NoT6ZUKnk8HjouLi5uaWmpr68HALlcTmSMJxw+fFgqlVIoFBqNFh0dDQCVlZWhoaFRUVE5OTnBwcEZGRlRUVHokZSGhkZOTo6NjQ2qwBwcHMxisQoKCtBbVVdXa2tro0nGd6mhoUFXV5dMJre0tNDpdKFQiOY36+rqiFlIxdmzU+PigEaD9nb7X39VBwWlp6fr6+ujXwYAiAPiPfv06UM0z5w5k5qa2tDQQGyl7fg1qVQqWlSL9Sr4SX0Pl5KSEhYWFhYWJhQKAeDSpUthYWGbN29+7cWampooDymKBTo6OmjPokKhyMvLYzKZzs7OqB7c0aNHf/nll3Xr1hUXF48cORLt3/f29jYyMpJIJNu2bUM7Jt+x8PBwVJLk9OnTKSkpS5YsQeeXL19OpKn+bYERkwkAEokkMzMzNTVVLpe3tbVt3749LCys02YtV1fXjrmlp06dumvXrl6ygAH7m/D/nD3cgAED0NLxSZMmAcCnn346b968qqoqYj35n9DR0fHw8BCJREZGRtOmTePxeBkZGa6uro2NjRwOJzQ0NCUlhcguCgBoYnT79u0DBw5ET7GbmppQhjpNTU25XC6XyzvuT38bWlpaaDQaSm6vUqnQYLljqn+6v3+zhgZHKAQAGDvWnUy2t7cfOnSooaHh2rVrbW1tg4ODRSIRMcZctWrVxo0bz58/j/YgoH0KTCYT7Zc9ffo0n89fs2bNjRs3GAzG8OHDc3NzExISamtr582bp6+vf/bsWVtbWyKrC9Yj4TCK/aHp06cDgFgs1tPTCwoK4vF4e/bs4XK5aA+oQCBISkoilvUQrK2to6KiUNResGBBcHCwgYFBcHBwWlpaQkLCihUr3mqfT5w4wWQyk5KS/Pz8hELhwYMHAYCYbQAAsLVVp6efOHDg9u3b69es+Tg/v6CggEQiGRgYcLnc5ubm/Pz8oqIiDw8PdDmbzd63b1+nT7GxsbGxseHxeDo6Omg3bUtLi0KhAACpVNrY2FhbW4uaAoEAbY7CejAcRnsyLS0ttPUQAE6cOKGpqYn2MpqYmGzZsuUvXx4cHIwOamtrUTwSCARRUVFZWVnt7e2ZmZloKHru3LlPPvmEWF/p7u6emZkZEhJChA+ZTCYWi9FqobdtxYoVBgYGqLdcLhdlvM/Pz+94DdfRcc6RI3NQb728AECtVg8YMGDjxo2JiYlLlizZvXs3yuu8YMGCI0eO/NFnVVRUtLe3E1nuT506FRcXV1dXh7IcbNq0iclkZmZm9qTMA9hr4TDawxE3p+gATX2SyWRU2/LOnTv6+vqenp65ublFRUWlpaVz5sxhMplXrlxxcXGxt7dHr/Xw8Fi1apVKpXJwcNi5c+fQoUMHDhyI0m4CwKNHj0aNGtVxmfrTp087prW/cOFCWlpaVVUVkUXlLXFwcECztObm5hwOx9XVFZ13cXEhk8ltbW3Lly9nMBjt7e1z585NSEhYtGiRRCL5+uuvR44c2adPnxUrVjx8+DApKampqQm9kDh4rSFDhgAAsQ1pxowZ48aNS0tLu3v3LgBs3LjR2Nj48OHDb+3rYv8VOIz2apWVlWiyTyAQlJaWZmdno2lEHo+HlvggNBqNxWJVVFSwWCwOh0OhUAYMGNDe3v7o0aMNGzY8f/6809sOGzZs8+bNiYmJMTExxsbG48ePNzQ01NTUTEhI+OWXX4gaTW/cl19+iQ5QET0ijK5duxYATp48+dFHH3388cft7e1ZWVmtra1RUVG3b99+/vy5WCwuLS1dsWJFQEBATExMeXk52n2AZnj/HN6JhOEw2is0NDSg3JcNDQ1cLrelpQWVzWhtbe2Uyv6bb76h0+mpqaloXSQAcLnc+fPnAwCDwejXr5+vry+a9dPV1d2wYYOtrW2nDaNMJtPS0jIlJeXo0aMffvghk8lUKBR3796dMmWKqanpuXPn3l4YRcRicXFxsZWVFYVCkUgkurq6KpWqtraWqByFhuQ3btyoqqoaO3bsoEGDRo8eTcxyyOVyT0/Pjs/lOvqj4nTEcn0PDw83NzcAQEtr58+fj4vZ9Xh4wVOvsGjRInSwcOHCuro6ovZvVFTU3LlzIyIiFi9ejM5s3bo1IiJi7NixxGtFIlFlZWVxcbG2tjaNRrtz505aWtpPP/2UmZl58uTJtWvXNjU1RUVFAcDPP/8cFha2f7/5jBmzb9++bW6+eOLE3ZMmbfLx8YmPTzt+/HFTU9PbLjVcUlIyefLkxMSkJUs2pqWlnTp1uqwMiouF33wTbmVllZaWLxD8duXQoUODgoLQ5MPYsWMbGhq8vLxcXFzi4uKIzQVo4VdUVJSbm1vWkiVtY8bE29t/NmLErFmz6uvrz507t2DBgvj4eAAgk8loBT6JRKJQKMQ7kMlkHEZ7PBxGewWVSiUSiUQikVKpBID29nbURH8tKSmJiIi4dOlSRUXF0aNH0QrTx48fZ2dnV1VVTZ8+XalUnjlzprS01MrKytTUlEwm79mzZ+nSpQcOHNizZ8+ECRPkcnlKSkpiYmJtbe21a/D990ChUJKTFXl5gJb/6+lZaWoGBQUFvTam2NpCcjIAwNSp8OgRoIc6fD4sW/aPv+mJEyfCw8PHjZujvjDQfcKkmVu3xHxx8vjxOoEAZs+evX//T59/vmLRokVKpZJOp6PsTQqFoq2tzdbWdtKkSSkpKf7+/jt37kTvhrYhPHnyZHJrq/XFi1q3bg0uKlqQkyOXy9H6U4VCgX5SrDfDN/W9QnV1NbpLRdVDU1NTUdPAwGDw4MFTp049cOAAi8UqKys7duyYXC6fOXPmqVOnmpub7927t3Tp0qCgoODg4IqKirq6OjabzWaz29raUN5iOp1ubW194sQJtCfS0tLSzw+ysmDKlPHHjm2srt4cGVns6upQUsIdNMjW1va3/kil0ocPH44ePRo1BwyAPXvgzJk/7H9KCri7A40G2dlgbQ18PtjYAAC8eAH9+zX6brAAAB66SURBVP9/VwoEAkNDQ41rFw8oF2nV1msBTH22/HPGSjIZWlrY06dH02gwdSpkZGTQaLSZM2fOmzcvOTl5/fr1NTU1DAbj2bNnpaWlZDJZrVYvXLgQlVchkUhGVVVaMlkzlcpRKLSUyrS0tJEjR7a2tmppafXmPK0YgkejvYKZmVl4eHh4eDgqHzJw4EDUZLFYWlpaSqXS3Nycy+W6u7uTSKSNGzdu27YN1dIQi8UdNzsKBIKCgoKCggIGg9FxMebHH3+MkpXo6ekBQHg4nD3r1L//1vz8xz4+ypUr9efP/y3d0bx58wAgMjIyIyODeDmdDjNmwO8l4ODyZQgLg98nHgAA1q4FtHZzxw6oqIDTp387/3tC6v/j4+Nz69YtUqtIW1ovZnIAgEZlfvDBl9XVcPkyJCSoZbI2dKVUKg0KCtq1a5efnx+dTg8ICEhKShIKhZ988omhoeGZM2dQAgEbG5uQkBBGnz4AwFEoAKCFQvH09IyOjl61apWvr++Z38P/yZMnGxsbiZEscYD1eHg02iugHZzoQENDgxgWovMdN/kAwKZNm4ill8OGDYuOjvb19UU7zfv164dqB124cOGHH364ceNGXV3dwoUL0a26r6/v1atX29qgsvJZTU0CmbwiMNBBTw8sLMDd3evpUygvBz4/8MGD7LKystbWVqlUSpQhCQ6GWbMATVxOmAALFgCfD0QgMjUFHg9KS39r8vnwyuqA30ybNm3Xrl0/JGUuJ1GY4mYAqNTtT6HQRo50tbT8pbRUeP16cb9+H1ZWCq5cuUKj0RQKRUlJiY+PT2pqqq6u7vXr19GOI1Qf1MzMzN7e/v79+wn6+hbm5kyp9GVFRbSDQx9t7efPn/N4PG1t7RcvXqCPzs3NHTVqVFpaGmoSB1iPh8Nor7Bx40Z0EB4eDgDEA6WNGzeq1Wp/f/8FCxY4OzsrFAq5XF5TU2NgYKBUKquqqiZNmnThwoUpU6b07dt36dKlKKccAAQHBxcWFvr4+CQlJV2+fHnZsmX79u1btmzZ4cOHJ09efu0ao6HhibFxMZX6SUJCalBQmIMDPT4+urp6oIuLTWqq/erVq2NiYtAaTwAwMEB9gyVLgMEAHR0AACoV9PQgLw+eP4eSErh0CVat+u3rNDVBYeHrvymFQjE1Nb2nU7zcLXjFKKe7955Rv7jibMm6cqXg+vU0Q0NvS0vDy5czm5s1WltJ33zzjUQiiYqKUqlUU6dODQgIePjwIRplL1u2rKamBq2uzcvLi0tKQv9V3Pr224mWlnZ2do8ePfL09Jw0adL48ePR6qjHjx8vX768qKgINUuJqI/1dDiM9nYKhSInJ+fAgQM1NTWbNm2aNm3aw4cPXV1dx44du2fPHg8Pj759+44cORLtaPr0009bW1ujv/lmnqbmlatXryQksN3dMzMzFy5cWFRUdODAgSVLltTUZGRnN7i7u0dG7m9sbIyNjVapVBwOSKVlDIYLnw8MBp3LpaMyxQjKXGphATExAAC+vgAABgawfj0AABrtGRmBj89vN/7OzjB5MgD8dn1xMZibA4rJzc3N165du3TpEgC0t7c/4E09v5QFABcvNikU6ry82MrKShMTEzp9Co1G2rJli1KplEgkpqamurq6AIBKH6NeEfX7KBRKxydjVVVVVCq1trYWJdOj0+loohlFT1tb2z9aLIX1VDiMYiAUCmtqagBAoVCw2eza2tqcnByBQEChUE6ePKmtrV1QUEBsDFXn5g64dg14vBAAn5KSoTU16FF1XFxccnJyWVkZl8vlcDirV6+mUqkoT92aNWsoFEpycrKFxRipFPr1+wd9Ky2FESPgwgUAgAULoLEROBwg8uuPGAGVlVBdDVFRgAbcFRUVDg4O6K80Gq2xsRFFt7S0tOHDh4eHhzMYjF27LpiY9KmutkfL9dG3Q9mU3dzc0KrPjszNzZ8+fTpo0CA0szFq1Cg3NzcajYYmQ+RyObEsF+udcBjFoLGxMSkpCQBQ3uXs7OyrV6+yWCyRSPT999/X1dWh3U1o1TqFx/Pk8Vq0tHTa2kyVytGjR//444+odDOFQomIiJg4ceKzZ8867vvcsWOHpqbmt99+6+oKr8Sov2BgAPfvw7FjAAAUym+x8uOPf/vr3LkAANnZQAQxU1PTkpIS4uV6enpobFhTU1NbW7tz504Oh2NoqPvy5SVLSzKHw9m3bx/6Xn9i586d+/btO3/+/LBhw7y8vFDEtLa2ViqVarXa3t4+JiZmzpw5X331FfFxAPB3cmhhPQN+Uo+BtbX1rFmzZs2apa2tffz48fb2diaTyeVyUZ5NQ0PDGzdu1NbWoovVNBoA6LS1AYCCTM7Ly3Nxcbl9+/awYcNQejoAMDAwWLx48Z49e65evfrqxx0/fjwvL++POoNCuVqtVqvVN2/eXLPmyxcvttLpoh9//BFd8MMPP3S8vrQUmpqgtfW3pq6urq+v79y5c+fPn5+SktIxiXJNTc3u3bubmpoGDRooFFauXLkSABISEpYuXfrnkZTJZH799dcHDx6cMGHCqFGjUIo8d3d3Ly+vc+fOEXXqjYyMqFQqMSfQcTct1rPh0WhvRyaTiX2fgwYNunHjxurVq0eMGIGGoqgWUysRpQAUgYGx/fuPotFUKtVmlcrS0tLX19fKymrlypW3bt0KCwszMTHJycnZvHmzQCCwsrIaPHgwehy/ePFiBoORkZEhEAgaGhpe7Ulubu769ettbW2pVKqdnZ1UKn327FloaOj169clEklcXNysWbMA4OHDh2hzKtK3Lxgbw+8P/AEAli1bhvZKkUgktH4eAJYvX37r1q1r167179/fw8Ojqanp3Llz586dy8rKem1tD4lEEhYWhjYmoTwmCxYskEql69ev3717N3GZo6NjVVXVa78O1nvgMNrbUSgUIipNmTLlxYsXI0aMAABtbW1TU9OLFy8CQMeE8Dps9uCEBNDSIgN88fLlt99+O2PGDBKJZGNj4+vri+5kJ06cqKurix7aoJWkAKCjowMAv/zyi5aWFo/HI2I3YceOHUeOHEGZPk6dOpWTkzNx4kQymWxoaMhisZ48eTJq1CgLC4vq6urKykozM7O6urr6+npnZ+eOMRR5da8Um81uaGggkUgNDQ08Hs/Y2Bh19dViKsiPP/44YsSIkJAQmUyWm5uLdnapVKpOEdPd3d3d3f3VpKtYr4LDKPZ/9PX1KysrUe53VLAeVcsgYo1EIqHT6WgZEADY29s7ODgsXbpUrVYvWrSIWMCEdpTL5fLQ0NDjx493vGXetGlTdXU1ugsGgJycnD179rDZbBcXl7S0NJSTNCYmpqCg4NmzZyUlJWw2e8CAAT///LO2tvaePXtSUlIyMzOzs7PNzMyKiopSUlL+ZjZPOzs7oVDo6OhIIpFeu69/586dK1asILpaUFCAgiN69HTr1q36+nqUkwXDOsFzo73X5cuXiQOhUHj//n0A2Lt378SJExcsWJCcnEyUvnBxcWEwGGKxeMaMGZ1KFa1du/bw4cM//PBD//79if08Z8+eBYDIyEhHR8dXY1afPn369u2LjtevX79v3759+/ZNmTJFQ0MD5VcWCoUymWz16tVqtdrDwyM+Pj4mJsbd3d3V1RXd1yckJFy4cKHjzfVfunHjBpvNzs/P9/f3z8vL+/bbbwHg2LFjKDXJ1atXO0VJS0vLjt80ODg4IiKiYxJVDCPgMNp7XUDLiAAuXLjQ3NyM6tQ7Ozvr6OgcOXLE398fPYQBgFmzZunp6Z08eXIyWq75/3ttEqP29va6urqSkhKiNOlroSVWAMBgMD744IPFixd7eHh8++23TCbT0NDw3LlzOjo6Li4uNBqtX79+33//PQCMHTuWzWYnJyf/Vpzub8jIyIiJidm/f39oaKiNjc3Tp0/Rptj79+/b2Ni0t7efP3++rq6u4/bW2bNnR0ZGrlq1auHChfD7LEFlZSX60Hv37oWHhxNxtqSkBA9UezN8U9978fn8Q4cOAQBKGJqeno6af5SyiM1mP378GAAcHR07nlcoFC0tLbq6ulKplEajUSiU9vZ2tVq9fv36vLw8Ozu7P+mDWq0mllKxWKxbt27Z2NicPn3awMAA5ZZHe/A3bNhQXFxsa2ublZVlampKpVJ5PJ5KpUJ1ShoaGkQikZWVVUVFxa5du0gk0tChQ5ubm8eMGWNsbHznzp07d+7873//I5FIzs7OAoEgOTl5zZo18PujMxqNdubMmQcPHhBbZgGAw+GcPXtWLpejmQpUmunUqVMbNmxQKBSRkZE6Ojrr168/f/48AOzevXvbtm1oLhjrhXAY7b24XC7KsXTnzh0AsLa2Rk0UK181bdq0zz77rFMd9idPnhw4cMDV1ZVOp1Op1ICAAG9v719//VUoFI4bNy4qKqqtrc3V1bVv374cDsfLywtlNvn495Wfy5Ytmzp1qpWVFVGwpKPFixf3798/NDRUV1c3PT395s2bDg4OXl5et2/fXrJkSXJyMlpp/9NPP1EolKVLl65atSoiIsLU1LS6uvrAgQNoAVZ1dbVSqUTPiBBvb2+UpHnq1KnEycDAwFc7QMz2EkpKSlJSUsRicUNDg5aWFpqFQPVHsV4Lh9HeS1NT09raGh0AAJvNRs0/STPcseASsm/fvujoaHS+02Tlli1bPvroI39/fx6P9+jRI/RYqampqePy+OHDhwcGBra1tWVmZk6fPp1CocTHx4eEhBQXF+fn57e2tiYmJoaGhvbv35/BYNDp9NGjR2toaNy8eVMikaAgeO7cObTECgDQtk4A6NOnDwBs27ZNR0cnNzd31qxZJ06cYDAY5eXlY8eOJZPJHZ96JScnx8XFhYWFlZWVyWSysrKyoKAgdN7BwQHNORBSUlJQDdGAgID79+/HxcUBAJ/P/+c/P9Zz4DDaey1YsIA40NfXnzhxYqfzf4dSqewYW/ft22dkZFRSUjJhwoT8/HyULM7KyurRo0cnTpy4d+8en8/38fHp+A4kEonFYu3fvx/VGT1w4MCOHTvOnTunr68vEonu3bunq6tbUVExfPhwLS2t1atXo32lxBojPz+/hoaGoqIiAFAoFCqVCi32BIB169ZZWVmdOnVKW1v77NmziYmJXl5eRkZGy35PB/3NN98AwI4dO9D2+aysLIFAcOPGDRRGr127Nn369E5hdNKkSTY2NnQ63cDAID4+fu7cuYCTOfV6OIz2XkRJS3Tg/ftOdeL838HhcMrLyy0sLFBz2bJl3t7e169fFwqFGhoara2tLBYL/WnOnDmjR49OTExEG087aWpq+vXXXwFAIBAAAI/HGzZsWFFRkYODQ0RExIMHDwoLCwFg586dGhoa6Dk7YmlpaWlpiZbQT5s2LTQ01MPDIyoqaujQoQAwefJkPz+/mJiY0aNHoyKgd+/eJeZ2XV1di4uL6+rqysvLIyMjCwoK0PIp9PDt1a1W7u7uKPuqtbW1jo4OsfTVx8cHrYvCeiccRrG/kJGRYWtry2KxSkpKNDQ0SkpKUPhITk52d3ffunXr5s2bNTU17ezsjIyM0JJSVD108eLFw4YNc3JyGjhw4KuzAc+ePRMKhaNGjUJNVMIIfp9SQM2Oa1T/jkmTJgUFBfH5/IcPH65atYrNZqtUquDg4Pz8/LKyMj6f39bW1uklNjY248aNKysrc3JyEvxepwmlbuo0DgUAYsCOfgEiiUnHbFVYL4TDKPYXLly4MHv2bBaLFRcXp6+v/9NPP6EgEhkZuXXrVnNz8+PHjxMXC4XC0tJSf39/hUJx6dKlefPm5ebmTp06ta6uDi0VcnNzs7e3F4vFaWlpCoWCCKNcLhfdSp84cQIA+vbtGxQUNHToUDQBOmTIkMGDB3/xxRdo0BcWFvZHu+DR7qmqqio0Yi0oKGCxWM+ePZs0aZKGhgYRKDv54osvfHx86uvr0QW+vr4AcP369Tf1G2I9Gw6j2F97/PhxQUHBixcvAgMD5XL5rVu3AADVg+skLy8vMTFx+fLlIpHoxx9/5HA4AKCjoxMfH19bW2tvb89kMplM5tGjRxkMRscpxZCQ36qMjB8/nsPhBAQEAACdTkfn0fIA4rn5X95Bm5ubd0z6OXjw4AMHDrS0tKD6UZ0MHToUpRFxdHSUSqXEIHTUqFE4/R32d+Awiv01Op2uoaGBYhmJREJRjKgh3El5eXlSUpJIJLK2tnZwcFAqlUlJSQUFBSikImg1KLp3Rog67+iAGKUS5/8RIvzp6+tTKBS0olNHR+e1i1iJzVpoTT6xu/TVXf8Y9lo4jGJ/zc/Pz8bGBmV7otFow4cPBwCUYf5VTU1NZWVlaBaSRCKhJUR1dXUdwyjy2qWab8Thw4fRAdpQgEqnEE0Me7NwGMX+QmBgIBrNubm5MZnMCRMmoPPBwcFMJvPV693d3SdOnCgQCNAT+cDAQD8/P21tbSJjKYb1MK/PdoNh/w6Px+PxeIGBgWKx+OrVq0ZGRlZWVlZWVjk5OSKRiKiIh2E9CQ6jGIZhXYIzPGEYhnUJDqMYhmFdgsMohmFYl+AwimEY1iU4jGIYhnUJDqPYWySRSHbt2hUWFvb06dM7d+5UVlY+evQoPT09NTU1IyOjqakJFYACAOIAfk/yBL+nQyaar90RL5fL6+rqAKCpqempqWm+ufmv/fs/S0iIiIhobGxsbm5etGjRW/2OGIbDKPYWrV69euDAgbt27SKTyampqfX19ZGRkYWFhYWFhbdu3aqoqDh27Bi6kjgAAKLg89KlSyUSCdHsWJ6ewOPxUI2m2mHDfGtq+lVWfpiVpX/3bkVFxenTp9Vq9R/VRMGwNwXvYsLeorKyMlRSadCgQY8fPz5y5EhGRsbLly/9/Pxqa2vHjBkjFApRgdJOCeRRpWWUFEqlUqHmH61xbm1t5fP5JkVFNJVKwuEwmpvFDQ0PHj/W1dXtmGkfw94SHEaxf0+tVt++fTs4OBgAbty44enpWVdX5+7uLpVKExISAgMD1Wq1Wq0mqpLMnDmzra3N1taWQqGganQ0Gs3ExAR+L2SCNDc3b9q0CQAyMzMBICcnBzUbGxtf242srKzo6OgJNBpHImE0NwPAMx6Pw+F8/PHHM2fOXLt2LQDU1ta2t7ebmZm1tbVVVFSo1WoLCwstLa3KykoajWZkZPQ2fyesh8NhFPv3lEpldHQ0CqM//fQTh8NJT093d3dvbW29fPlyYGDg2LFjv/7666CgIJSpRFNTk0qlOjs7q1QqKpWal5cnl8sHDRoEAFwul3hbDoeDsoVWVVUBgKurK2qirHcA0NbWVlxc3NLS4uXlJRAIrKysli9ffo/DSQwL4+jo5DEYAZs3l166BABkMhm98/Hjx7Oyss6fP19dXR0VFaVWq2fMmGFnZ/fw4UMul/vRRx+9458O60lwGMW6pKSkJDIyEgBEIhEAJCUlsVgsVLgYAObOnZuVlVVUVOTt7e3k5GRgYPDZZ595enqSSCRDQ8PExESUmw4AJk+eTLwnqvcJAHZ2dhQKhWiiAx6PR3Zz66dSpavV2dev33z4MDc3t7GxcdGOHaHr1hl7e/+6fXuYpyeLxcrNzWWz2du3bxeLxdXV1UKh8OnTpyhPFQCkpqbW1ta+fPkSJWnGsH8Nh1GsS/T09FC+zqtXrwKAhYXFgAEDmpubUeVhAHB1dXV1dSWuJ8Z9hw4dYrPZRHGRcePGEdds3rwZHaCSc0QTHTB9fY3a2kCl8gV4GhlZSaVSqdTa2lqpVIqK4qG8U/b29kTRZpRpn8fj1dTUNDU1oZNNTU0aGhoo+mNYV+AwinUJm812cXGB31PTm5iYuLi4NDQ0/OULt27dKpfL7969+08/kS6RgEoFJibA54vb2rIrKpqamq5cuaJUKlEuvlcfzQ8ePFggENjb23/66aeFhYWnT58GgKCgIDs7O2LgjGH/Gg6j2L9HoVCmT5+OjmfNmmVjY6OjowMALBaLqP72J+h0OlEnGVEqlSQSiUwmq1QqlUqlVqtRzSWFQkEmk1Hl5DY7O256OvD5CjLZYvToSUJhRkbGggULfv7555CQELVaHR0d/epnWVtbo2dZXC530KBBarUazZk6Ojq+Wm4Pw/4RnCgP+w85deqUkZERUYc5Pz9/+/btXC73u+++8/PzQ+lKla2td2fMaGluThKJ5kdFZWdni0SikJAQFEnPnz8fFRWFq9Fh7xIejWLvGQqL9eGlSwCAhrtEeaWff/65pqamvb195MiR3dc7rDfCYRT7bzlx4sS9e/f4fL6Pjw8ArF+/XkNDIy0tzc/P7y9f29jYKBaLExIS5s+f/0cVmDHsjcObQbG3aPny5cSBXC5fs2aNWq1etGiRh4fH9OnT+Xz+F198gS6YMmUKOpgzZ05ERMSXX36Jmlu2bImIiPj444//zsc5OzuHhoYePXoUx1DsXcKjUewtqq6uJg5UKhWfz79//76lpWVDQ8OpU6fUarVCoUAXoAMOh8NisQCAyWTq6emZmJigMs56enqvLZ/3WugxF4a9M/gRE/YWBQYGenh4AMDNmzczMjK8vLzMzMwYDIZYLL5z545KpXJyckIP6zMyMjomecKw9wgejWJvkYGBQUREBACgvUPe3t5z5869f/9+bm4uusDDwwNdQGz0xLD3Dp4bxd6pwYMHy2Sy9PT02bNnNzY2ElvpdXV1AaC8vHz27Nnd2kEM+8fwTT32FikUCiqVShx0ar56fUREBI1GW7JkybvuKIZ1AR6NYm8RESvRQadmJ5WVlXK5PDExsVPuUQz7j8OjUey/paioyNbWtrt7gWH/AA6jGIZhXYJv6jEMw7oEh1HsPfPTTz8tWrRo2bJldXV1nVKRYli3wOtGsfdJfn5+XFzc8ePHBQIBk8nMyclB54kDDHv3cBjF/tPEYjGFQtHQ0JDL5e3t7fHx8cOGDQMAGo1GJpNzcnLCwsLglcKiGPYu4TCK/aedOXMm8MkT69xciUr1y7hxl589Q5VCfvjhh6FDhzo7O+NNUFi3w2EU+0/j7t9vVlYGIhEbwEetPs3hpKSk7Nix4+TJk05OTp0ufv78+Z49e86dOycQCC5dusTlcj09Pa2trRMTE9vb24cMGdItXwHr8fAjJuw/TYPPp4tEYiYTANTt7U1NTfb29pWVlUKh0MDA4Pvvv0eXoYPo6GipVAoAUqk0JyenqKhIKBQCQE1NDZFrCsPeODwaxd41iUSycuVKCoWiUqlCQ0OfPn06ePDgwsLChISE9evXR0dH6+vre3p62tnZJSQk5MlkHwEwxWIAEHC5/IICiUSioaFhbGwMAAYGBug9DQwMGhoaampqUlJS5s+fDwCampoAsHfvXkNDw6Kioo4FnDHszcJhFHvXIiMjAwMDJ0yYIJPJ8vLyhEKhQCCoqKhoaGiQyWQFBQWtra1oI1N9ff2PpqYTJ0yw1NE5ePo0bfJkpwsXdu3axWAw/Pz8tLW1O76tvr7+vHnzzMzMdu/ezefzd+7cCQArVqzw8PD45ZdfZDJZ93xbrBfAYRR71woKCkaPHg0AGhoaAHDr1q2MjIyEhAQdHZ3w8HCUuH7fvn2GhobFxcUcLtdy2zYAiEtN3RUczCsrQ4U8LSws+vXr9+qbW1tbv9Mvg2EAlPDw8O7uA9bTxMbGmpqa0mi058+fq1Sq+/fvOzo6AsDVq1cdHBzy8vIoFAoab9bU1JBIpJCQECsrKzKZvG3btri4OB0dnc8//3zGjBkMBsPY2HjQoEEA4ObmZmZm5ujoyOFwUFNPT6/T5xoaGqKCyQwGw9PT09vb29TUlEqlmpmZ2dnZoaiNYW8cfsSEvXnXrl1ra2sDgISEhLKysnPnzqHz58+fV6lUs2fPjoyMXLVq1cKFCwGARCKhv6KD+vr6R48eEW+FissDgL29vYaGhqWlJdF89XOZTCZKYEqlUo2NjfX19dEMqba2Nq4sgr09+KYeeytOnz6tpaUVHx/v5eXV0tLyww8/AEBRUREAsNnsc+fOyeVyOp0OAO7u7lKp1MPDQ0dHRyaTubq6UqlUV1dXAAgKClKpVN37RTDsL+Ewir0V/v7+XC63qqoKABgMxvDhwwHg5s2bxAUohpaWlpJIJEtLS4lEkpWV9fz5c21t7YyMjMLCQl1dXWIoimH/ZTiMYm+FtbW1vr4+WpBEo9Hs7OwAoNOzdQBITk6mUCiWlpbl5eU3b95saGgICwsbNGjQ06dPbW1tcRjF3gv4ERP25nE4HCsrKyqVymKxzMzMzM3N0QN0LpdrbW1NTIYCQE5OTlpamkQiyczMlMvlYrFYLBYrFIr4+HhHR0f82B17L+DRKPbmoW3vAODs7AwARkZGqBkYGPjqxQwGQ1tbm8lkikQiANDS0tLW1kaPhjDsvYCf1GPdzNnZediwYd7e3qjp5+c3bNiw1z6Ix7D/JjwaxbqTs7MzmUwGAGNjY39/f7FYjFYmeXt7o5LLGPbfh2sxYRiGdQm+qccwDOsSHEYxDMO6BIdRDMOwLsFhFMMwrEtwGMUwDOsSHEYxDMO65P8BDVs58gcTphMAAAVmelRYdHJka2l0UEtMIHJka2l0IDIwMjMuMDkuNgAAeJzt030w1VkYB3CXSF7ulSRFws1beWfphd/5HkRL6W5pspQNm8gokZd2rVBSd1TMpljWa3Oz28pbWWrd3y/tblaUVTaTkpJK5CWVatFe087sX/vn/rU9M2c+z5w5c+Z55jxnWFreLScLgdw/EfW3ybyZcqunE4V/U+m9Mz7438h7L++D/09VtaaZ56Uoy+XeL7mo+IB0mI4UMspxQTjsZg3x8npm3y0O3Wsk+L7ImCG6Q2jrugmeWRnjc6oAC03z4arSId1UvhiScymyfTO2eDtLbrAOiD0e0lCt3kcMNFwQ9XQmGc50RazFW3KwtIrcjYpDpPUtYl/tSKaoBNfqxphxm2ym37gOEUREeriL0hjze5iVrgZj43xpy5teMCKKkOBUttdMgYo9j0B0clg62q5Ec/wlsjN6rKRLh/rOuI5Knm1DfJQWtexsna6NsS99h7SZQxDtzyPONX1odX+FWj9Nsi/6AjbkqtP1PoUuVollmB2hTL/rfym10E3GgRV9SHaVSIfqYnHUuhnXqyLZBScssbCSj8GMRDYxToeaiKdI1rYmJsjAiH62X0LOOpYzAWed6MJADyyTGHDbtcchcmtCyy59rl/MgGm+jwuv3WT3lBMb3mMsPWhKTk35E2cjfcI/30h85FdBeYcBEfrr4HbpYzJfoMXe81GEKCSShNqYwtDWDmf1UqA2+JDM2BGPsKUn0TVogGcCCQRcDt6EqdHC3lLC13iHqHXq1PkbFfinzKEB6w2on7SfKFYb0oHJKtmbLudURtIRc7cC1z46wkV8ngFH0ztQ6zNHw7oelNq+xsbMnVBybwHnIkeDJC7kxZIeqXXxA3hLS5BU3IrxsB2YFRjEpYlUaeLETnhkfculX5Cn/U5i9NRmsV1fOhCNuT9Lte//zuQe/wGjVqHMI0NPdmLyOXIMClGyK4WV1zKhRTaPoent5XxUqE9n/6FLr+o+cmm/dQ7iDZpU7BbC8qv9YMerJhJBDftFjDYV+lkhVZNlkwRCusbbAbvs/cgKqZAq3akl54qDyd1aHYrLJ0je2nxyWtmFtrWNMbIaGjQOOdLx0SHiZaQhjTxGaaoSHwGr7LiEN8NYGhqOsg5D7sDubmQ/r8fwmWXclkUCyJn2QSKw5xyqCAraJ6EyN4uVjtYTfns7njwoYv9UtiUPPK/CpyaVdMTakBtdL1AnPCGb80bSe7uWCW8yIpOq9qi6adyw8skwSVkWDvdqPdK4eAPqL8tDceo8+dpQHVdCJWRvUyXbV6EBYWENkU7y2V/nmMEpg89Yemqysn9BhifMyKs1nyIycAlaOzMxz9EBlUoeUO37GNEJ8aR1IAkPN6fjYUQzWeWyDcH6NaiN3wircyVoVy/Du21WcIqvQE9MAtbxpkibykUE+0oQfUVANb2bmZW/KNB8hXm0+20LiRkbQFsan97b4o6RMAEdjJOjnWJ5bDuiS/X0F1Gj0jASmGpCE6wsKUkaJwe8F1A7KwMauFMHMLeihp9cgu98IZdpUoSpn7KxN0nINUcdwpu3xdhcepyT9YtH4iJYnM/hokbCsda2BbUZOVxuXgaOVjYie9QYWkPXMZZ8E7LZQ8LlCRxaOYtu9t2N37K6ULR0EvxNW7EwtwIhX/FpeHsL83p3trQjdAyJXDSTWn2MncrgUaOCG+TikBPLJndjq8oZLDIfgEXrS9Rn/ojZhpdweqQJQZEl6HRlsdorCy/3iLjpfg8+9cAx52DO+PY8mpngBCXD05z7wBTSvA5DolLAiUKeITB2P4yaKrh9Feo0Z48CIsN8G3hf5hHtvwCq8Hdpn/CdPgAACBV6VFh0TU9MIHJka2l0IDIwMjMuMDkuNgAAeJx9WFuOXccN/Ncq7gbmoEl28/FpS0YUBJaBRMke/J/9I1VsKcNBgMgGfH1Vc06TrCpWz6cX//z9y9/+/Pfrv3/sy6dPr9f6P/9W1etfttb69PuLH16//vaXv357ff7+y68/v/n8xz+/ff/HS0TxE/3PR+wv3//4/ec38vr22o9qlZzXekzUol72HNXUev9JBU6fsNwVxLmEx8ufpW7L33EG3HmkPNOAE+C2v+LZe6UO3AbOHz+Spnye5ua3T504J95x5/UZ32ZUxHm9rWeLr9o4yglbOUrzBq5jtpRAP56O+p5dx84ABoD6pG5dSeA657gAGFl7PjEbWGfLQiOeY6tU0Jq116pRSgEnxO1UfLCMXAZcnuU2WiirgZ4ixl6HmKNUdgtHOAMoAKLU8CV8s5ZYFoBxtmsOoDYwxbZ3LeYemIY/hUnJqEWsi8GbYnd7bGGQgsEIGrFtIDeQKFLXUXu9ybMkXeOV+OFghe/IO5stuzBrzsajjgC51XfNyjkc1pnnVA8HNMY56jnHXGfpnI7jSaVxp7NNcaR6NGJPYkhepGk3+Scx4qlMj3lMzgcEBHXxM/LgcbUPK4+UDxRfDfQ6WT3wEHfHuUFczGIAOZ98IIUyknxFRBKo4POZT9QGHnWxrhv80krydJnJqEat3y1uVt1L0x3QC1hktWWUo7vrNghntSIK0mlFhECz85mnZw5tcxhEph98ia5ror8D6ZfCmaeEaiwBx0HYR0xrjfloNBL+k7mJNBwDXRWWtnw+My+PMGic7g2ClSXd98AEcj6zuiI9qWgJkJsWYNRehE6+250QFFvW7Ag/B+zAgSDhyWKTRqpkFrXmZ1uXBgZCpwN4JQSbCqFbHBWBW+DggbFPY7MLxJegkULm6yzqPRfNbgA5oDfSJw2q1gfGV01R9FxrtMjugByeix5ggAkpFtlhaIEM4J1PoOdnkxN+Eq660coNHQ/gNTf0ZB8SUmBpexMYhq8GMJuZMFUaCkQBy8pWT6yYxmF1uW58FB0EvuqtR6kP/rvb3RaNQ2BBQtUIlEISrVM+KLyl7QCvVlSBgePRgqNjD4FFs/CtjYSFBy2GZqMpFBBcFoQZSGsScVFAA8qB+yLb4NWofQB3MwOF41ToZWjQlKCoxOoYFNqnK1+gA8X/BDfYahUbuDHf7ZeW5C+7vvJg0DQnReHz3dGacJAI43njlsUO2BwEuv6hnKsenBy2xsJLrZJICey4MXKM/yv7IWXQNvwau4IeggODoDJmDlp/bWsBH4OADXGiogO9m8xnYgRfuXfgQXl1RoPjtoDwzhoVHe23V5Dad1tg6MWpI0fklPkxQIMPxWKma4L2aEg9SAo+mQSKf8XXYRVkzYMnGnpzmBZskhjqIhCuYtVWmIklwaFanvTRTnT7K12iFnrPY6KZDDyFHqBzg5tYMV8pRd1yC9pQLaokjQvkm9BsaGIijAPYBzjIIVLh4nsiq1+/Ff66iEw48KIR7+17tNNXA2FGabQ1hRtJ22vK9rkG/I4Ie0lv3pFF8bQTk4ijeNd+aCIb3O0Lcwf9kR7gzvCZYXJuXRImKLutOPJsKA/rk1QZ6nAOCVTmMKXtXZPbEtEDU5qL2jklLFFZyd4rq78LS7H7c9DOvRmCzAgW9dshZi6Ip/DQGr7g0UjQtigPkh7MOOwxlm3NWMgh+XMSDezJ60G26JUAAsyV4dVIw4pUPmkpch+zrmNZ1uhnrH55Z9fd09Ly07tjrblb4s4I3Wh1PQ4a47uuy9doUehVJpLlzvZOpDTkXfYXsxoECetHggqr+PfM7b2OkJNkhtLYPygnKVxCiHNI4nd1r0mk4HioiJT7JMXRjIeANBAdBtLvKc/pGKfchcoy0P3SD/XEDyQLL7YIe/p0LxFp1xh5ZNMI+jWEJaU1LmTVzcy0fWg96ofYTq4O4VgZEKA9co3+PbWv7lCBjBCrIdQYJ4/dwtAyXp3SHToghvBB+IC9yr2EnDtTXOo9o4Eu2va+sZAY96CdGbjSfnoHZMW1dTcK9jUWx4ykeYcDQ0cR1MbGjlmdkIPlvQM5nGRwk7Z09EkWtwAc0eeiTG8rRMi1voSgzXvx3mXRd4N3YNwnYqc5nXKhS/waDTORWXU2EHkVUZUEgx0rFz8iGXb8AFa/GkEjec3CoxPWTRvcfPh4d7W7wf4UWYOCKfaUPujbts0rlVwkt0nctCVNyvZO9uodqi3btZHbbs7ciFHCOa0da77eelfhkmLc40AezJnpCFeImE2q3Za1kPTQWGxfyBB848ihnj0Pes1tCW4Nh0gkVhyEysXVZs1z+rUX+CDoRqQHfJ6rG+lw6qziVpQbVzCmOSR8VIw4UbyxDOD1NqdF0FYwQq6VQJ5KuvI7sHqWJwLJkrECcQNmgLjklM+8ol5vg6mc9p+DsAE9kssBsk6k/CDS6b2EYaMzfXXBdfXD4sdmuoaJ6zBXCMiLS2qbMW4sbh+gdllXMOnm8eH9riOIJL1xQHdXT+cQpi4sAIyBlMdZJ0MQP+/7ed2j3MAkxCq2GGY/BczgymcqKBK8mWNPYlXjAwy7PtynVzSZDvNDMMrpdfn/FRL2XvsH/ruwNZn6DGxyyrX4ggn94XJY1RjUG/im4D+DF/xMPzxVVh+A8uSqRgazxBWMqxp7KeaNXqShPACqemOeRpyn6GHk9mH+on0ARkxosH87Ikiz99rscwX/9u3Lh18k8f9//hIKnz/9B5KIk4KAONJXAAAFtXpUWHRTTUlMRVMgcmRraXQgMjAyMy4wOS42AAB4nO2WQa4jNwxEr5JlAvg3SIoiKXg5m1nlAkFuMofPK3kwBwiQXWDANsSWRBariv3t+faffP78+fnr+/e////611+//fg9n+nT/fqyJ73s5Cue3cvmnY/ttRSoXVMvf/Lstd/xTGS4Arb3Lt+Eek4OobPT/WXPXnZ4Zj2W5m/Xeg7PrelhdbZVslzjvni8fdXu1378cOKbXLrMD5E4vuYQ6Z21iIyvrNHtq6qtXvWcXuNczubOG7E87a9+nEzzTRZhO15f/phPrdfwaKdKTM9zD8vqs51IRu23Epldt3g7J+N1nr0XCRRPnrhbLFd4bkLRnaXQisrXLyT7OTPV7yYwnOEPz991b8pgvc5Wkbu9qihy2D/veSbPES7W3VNOhKy8iOwovzXmWJA6XbK1grO81jq3+siOEsrrpJPW8lm3XySet8PtEQfE+N43MLVDleyYs9SvmZvZOtzm9vjiedpnNqPlNVDk5a40XAADEpfDHTcI0U9T6nB37IlQIHvTYRjSHSod7pxbSNfeuuRZcCmJhM/o7trJ9fyZLfbQ+/KmrU6WziXc2hOt9ajy4fA1JjaOnTrvLwE7a28CNnHUmRMdp6m8PCA3SNCYLRQXCb5F5EkdDSCj9Y40cZ4UUx1yiJ5a75VJP7o3hKKzY2vU8bZ2dRDQD+c74ilRxA+aAkmoqBLotMdF1gq67QrSACaOd50fACseojkU8yWmxgQBFBTimmSrC0CpRHUzZAB65GNqfgf/EAedWuRDV0yd6wSKy58FeKNGDIBTmQ26ZwN4nKKlJVRfX0tor0IZqkwBKvdSSicWDERMbTsRE2w6IXaALblpZ8wCPPS5tU7/qrnCopZfkdWlkxMAXXRCdtIYZB/xf0HYRdGAgC1QvJ91dYwpLHGEm6lMrV50uLUnTVlsqLugEw6CHOYK9qeY+tkHH+NUMFofV0APagZpHzEwELnkRyHZdG9m1bwPhZwuJUxgCzQS5ajz4AOXzDNYgUqn4mkp5phQ1J85spE0ihAmkX4vTzi5LxOwOrmr9yV9Jm6rdYgq7WV8JIky8R+kCiHloWhiLgMRRUiS4/nTQsNuq90iP2Kl1+wYjPN6m6eJJI4+TovlsNQvDXt2is74UIgLqQkgecfcDuKVS7wdt0+ArFJ0jFVxrhthBPeko/TQPVsUoCunPi1cKbEmlnbo+Z6Rs1HIrkujIpcUGxbeJPe0uLOlAPKIuR+WCAFaf23B0AxdxE8vO4omAeZNovx9PdvyzrLpjH6pMjC41hJ2JHw4DIxx3f2i7qMTGCCuO3BADE/9HL/PRVnLvehrfDwH49hXmVGhfMkluOPaEf+PMsbzxKNiytkB4MVo1WmIDthSNk8TRYateck0cYkaQNTaw7RAxQubXh+3IC1xZAOh1vkdudU+Gkwcv0DyajxTwzQYJj9JeETjxCRkeqNE5f6ZGjpYoX2mcmsqMSX8MxRw+UT4iKDEfsbb0nSnmpS+FhbR2mDSO+mTrabuQfYhM4pj4oQhJ2miA6ZyDtZ6jQKrRGnifrJVrsl+XXBUCBsq5b5yU70/yPe975yChpp4lv7hMSxOmYFZYsv5MOXXHUeMFJoEcrW6IbKV/BpTgxBMCiAyuwZsIHXdjrnlwpbWXVJC/bjr1a45jN8c8RiMTSNnM20lFV5FRjlViighEELSGr0caKJ3m8Dnp/U8IO8tofin95vXML896jtxjjEnr89R1H0/4x2JuiUIKzkg7mv1scmARwL8SFI0CJZd5ZMyOnhLyCuuW9LqUOu4Rm8svJmJML0K/CgKDaq4AL87LFrjcWOhqBRct+ivYRFwPH+1Gsba2GeKLEaZ/E+kELMxIbIH6Ii+kxPuh6YI3Nh3vKy5FoTIm3WdNNqQWIjGDjK7r6Pgc87nddRnXzHy8lHnjx//AJBznjIkaSh6AAAAAElFTkSuQmCC\n", + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcIAAACWCAIAAADCEh9HAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3daVxTxxYA8JMEEiCA7IuAIhDFpWrBHWpdqIqCrZVFtCrWCoobiki1yloVrFjXVtxxreBW1IdI1aIVWgFRka0QUVE2ZSeBEJL7PkyblxcWkVXk/D/0N95M5s5Fe5h778wZGkVRgBBCqLXoXd0BhBDq3jCMIoRQm2AYRQihNsEwihBCbYJhFCGE2gTDKEIItQmGUdR6mZmZT548AYDXr1/HxcV1dXcQ6hoYRlHrZWRkSMLo77//3tXdaaX09HRSSEtL69qeoG4Kwyhqk5iYmL17954+fbqrO9J6AQEBMgWE3olcV3cAdW8WFhb29vY5OTn37t3r6r60UnV19ePHjwGgpqamq/uCuiUcjaI20dbWNjY27t27t8xxSVSNj4/v9E69m8rKypSUlJSUFB6P19V9Qd0ShlHUepaWlqNHjwYAAwMDW1tbLy+vkpIS8tHu3btlCu+t3r17L1y4cOHChVpaWl3dF9Qt4U09ar0+ffqQQq9evQ4dOnT48OFHjx5dv35dTk6Ox+PFxMQAQFFRUZf28R8JCQmVlZVTp0599erVtWvX3NzcAMDf39/U1HTixImkjqSA0DvB0ShqH4GBgQYGBjdv3ly3bl1X96URpaWlJKDX1NQ8e/YMAHbt2hUQEPDNN9/Y2tqSOsuWLevCHqLuC0ejqH3o6+tHRkZOnDhx9+7dH330EZvNnjp1KgAcPXq0q7v2j7Nnzz569Ki8vFxXV/f06dNeXl40Gu3AgQN9+/bt6q6h7g1Ho6jdjB07dt++fQCwbNmympqay5cvP3r0aMqUKZIKvr6+MoV2V19fT/oAAD/++KP0Ry4uLqGhoRs2bMjJyXF1dRWLxTt27Fi0aFEH9QT1HDgaRe3pm2++OXz4cFZW1tWrV69evUoOrl+/3tjY2NraurCwkBzJysrqoA6IRKL79++T8p9//pmYmFhRUUHKZWVlioqKCQkJly9frq+v9/X1Xbt2bQd1A/UoNMx+j9pRcXExh8OprKycNWuWUCjMzc3Nzc3l8/kAMHXq1Lq6uo8//hgArl27lpmZ2REdEAgEFhYW06ZNA4DHjx/X1NQ0nNCqqKj4xRdfnDlzpiM6gHogHI2i9uTv719ZWfnFF19cvHhRcrC4uDg3N1dOTm779u2hoaEA8PLly47rg6WlJTmLs7OzoaGhoqIiACgrK8vLywNAcnLy06dP5eTwXz5qN/iPCbWbv//++/DhwwwGY8uWLdLHdXR0dHR0uqRLJJ5Ke/HiBYfDOXPmzKZNm/r3798lvUIfGgqhdvLFF18AwNKlS5uqUFpaKlNod2KxODc3l5S5XG6jdRYvXgwAixYt6qA+oJ4Gn42i9vHnn3+OGzeOzWZnZ2fr6el1dXea8/z58/79+4tEooyMDA6H09XdQd0eTnhC7WPdunUURXl7e7/nMRQA+vbtO2/ePJFIFBIS0tV9QR8CHI2idhAZGenk5KSjo5OTk6OiotLV3Xk7Lpdrbm5Oo9GysrL69evX1d1B3RuORlFbCYXC7777DgCCgoK6RQwFAFNT03nz5gmFwm3btnV1X1C3h6NR1FZ79+5dtWqVubl5ampqN5pIhANS1F5wNIrapKqqikxv2r59ezeKoQBgamrq4uIiFAqDg4O7ui+oe8MwippUXFwsU2ho27ZtRUVF48ePt7e376x+tRtfX185Obljx46RnE8ItQ6GUdSkFStWkMLKlSsbrfDq1avdu3fTaLRuOqAzMzObM2eOUCjEV/aoLTCMoiYJhcKKioqKigqxWPzw4cOHDx9WV1dLV9i8eTOfz3d2dh47dmxXdbKNyID0yJEjz58/7+q+oO4KXzGhJg0ZMsTJyQkA7ty5w2Qyo6OjAUBPT4/D4XA4HGVl5X379snJyaWnp5uamnZ1Z1tv3rx5Z86c8fDw2L9/f8NPd+zYQRJRnz17dvz48QYGBp3eQfS+wzCKmuTk5BQREQEAzs7O2tracXFxOTk5tbW1kgrq6uqDBw++e/du1/WxHWRkZAwZMkReXj4nJ8fQ0JAcpCiqvLy8trZ2zZo1Z8+epdFofn5+Li4u5ubmXdtb9B7Cm3rUpCFDhkgK+/btS01N5fP5e/fuvXXrVlhYmLW1dXl5eUJCwu3bt7u2n200cODAL7/8kmTY43A4+vr6KioqdDpdQ0Nj3LhxALBu3TovL6/Y2Niu7il6T+FoFL0b6SHqgAEDgoKCdHR0UlJSGu6x3F1QFDVt2rRbt27V19dLH1dTUzMyMho0aBCORlHzcDSK3s2zZ8927NixY8eOwsJCf3//KVOmFBcXz507VyYGdSOBgYE3btyg0WgAoKend/36dZIfoKys7PHjx13dO9QNYBhF/1NfX5+enk7Kjx49arSOvr7+vHnz5s2bp6mpSafTT548aWBgEBcXt3nz5k7sabu5fPlyYGAgnU738fHp1auXUCg8fPgwk8mUVNi+fTuJsB4eHrjYCTUKwyj6n6qqqj179pDy1q1bG63DYrH09fX19fVJMnkdHZ3IyEh5efmQkJDLly93Xl/bQ0ZGxsKFC8Vi8fbt28eMGRMQEJCXlyezyr5Pnz6koKury2KxuqKb6H2HYRT9n9evXycmJiYmJjZ1k7569WqZwtixY7///nuKor7++uvc3NxO6miblZaW2tvbV1ZWfvXVV15eXgCgrq5OdhxB6J1gGEX/p7KyksvlcrlckUjUaAUrKytSIG+xCW9v71mzZpWVlTk7OwsEAnKQJJkHgJUrV9bU1HRkr99ZfX29o6Mjl8u1sLAICwsDACsrK1tbWwAwNDR0d3dv9Fs//fQTKTQ6wxT1WBhG0f8xNTWdM2fOnDlz3ukGlkajHT161MTEJDExkUxWBwDJkqfq6uqumhASExMjUyA8PT1v3bqlp6f366+/KikpAYCampq2tjYAKCgoSG7kZcTFxckUEAIMo0gai8UaP348KX/22Wfv9F01NbWLFy8qKiru27fv1KlTAFBeXu7l5eXl5ZWYmNj+fW2Zo0ePyhQA4Pjx4/v371dQULh06ZJkvn1LvHz5klxRB+0Ojbqp7pTZDHU0JSWliRMnamhoaGpqZmdnv+vXhw0btmPHjuXLly9btszS0lJNTY1szFlaWtoBnW2RkpKS06dPA0BhYSE5Eh8fv3TpUgDYt2/fmDFj3qk1Q0NDckVkjSxCBI5G0f+h0WhlZWUyKUhazsPDY/78+dXV1Y6OjpJpQ4qKimTOUOdTVFQcOHDgwIED2Ww2ALx48YIsWPLy8pI8ukWojTCM9gjXr18nOUNTUlKePHnSTE2SelkoFLb6XPv37+dwOGlpaVevXp09e/bGjRvHjBnz5MmTysrK58+fb9q0iVSbP39+q0/RqOrq6rKyMgAQCoWSsaeSkpKFhYWFhYWKikpNTY2Dg0NRUdFnn33Wusx+q1atunDhQnJyckBAQHt2HXVzeFPfIyQkJBgbG+vo6GRmZioqKkoWyzdEwmhbliSpqKhMmjQpNze3vLz84sWL0h+tW7eOTv/nN3ddXV2rT9GohISEZ8+eLVmypLCwcPfu3Tt27ODz+dOnT09JSXn9+nWfPn1GjRr15MkTDodz7ty51iXq5/F4Dg4ONjY2uL4eScMw2lMcOHBAS0srNTXVxcWlmWpkUn1bwmhVVVVERER9ff2xY8cUFRUzMjIyMzOzsrKysrJ69+599uxZMkmzqQlVbVFRUZGfn19UVCQUCpWUlGRmWbFYLEVFxUuXLqmrq7eu/devXwOAjo5OO/QVfUAwjPYUc+fO5XA4Fy5caL5a20ejBw4cKCsrGz9+vKurq/RxsVj89OnTwsJCkmre2dm51adoyv3792k0WkVFhby8PEVRioqKmpqaWlpaWlpaGhoa165d4/F4fD6/1e2TByNkXlRTsrOzk5KSXFxcamtrd+/e7ePj0+rToe4Cn41+OCR5PxtNAKqqqqqurk7etDSDvAuqr69v3U23QCDYtWsXAGzYsEHmIzqdzmKxJDswq6mptaL95n322WdeXl5LliwBgNLSUj6fn5eXl5KSEhsbe+7cOTKp/siRI61un4xGmw+jVVVVL168AACRSJSVldXqc6FuBMPoh2Pfvn0yBQkXF5fS0tKEhIThw4ePHj1a5lORSJScnLx7924nJyc9PT06nU6j0dhs9ogRI1avXn3ixIm0tDSQWmUfGhoqWaokIzw8PD8/f9iwYVOnTm34qZGRkeQVE1k71I60tbX79u0LAIqKioMHD264rHPJkiU0Gu306dOtnofQkjAKANeuXfPy8tq4cWPrzoK6Hbyp/3C8fv362rVrAFBRUZGWlkaj0czNzckrHXNz8ylTpsTGxsbExAwcOBAA6urqEhMT4+Li4uLi4uPjZSKLWCym0+nJycnJycnkyLBhwwYMGEDKT548afTJpkgk2rFjBwBs2LCh82c4DR8+nBS0tLQWLVrUsIK5ufmYMWMSEhIuXLiwcOHCVpyihc9GZ8yY4ePjw+PxmtoKEH1gMIx+OBgMBhmC0en0wMDAiIgIFRWVoUOHWltbW1lZVVRUAEBeXl5ISMgff/xx586dyspKyXdNTEysrKysra1tbGxMTU0ZDEZZWdmjR4/u3bv3xx9//Pnnn0ZGRhUVFeTtUFOrks6fP5+dnW1iYjJ79uxOueJ3tnjx4oSEhCNHjjQVRr29vX/44QfpwosXL7KysjIzMzMyMn7//XcAuHr1au/evUeMGCGZdYB6Ogp9KJycnCSF1atXGxkZNfP3TqPRBg8e7OHhce7cuYKCAkkj5G6dyWRKtywWi8vKyiTtu7q68ni8hh2wsLAAgLCwsA67xLaqrq4mD2fT09MbrSD9M/zqq6+UlZVlfm6S0Kmpqeno6BgWFpafny/dgkAg+Pvvv8PCwn755Zc3b97ItC/54bzPPyX0rvDX6YdDMgacPXv2rl27Xrx4kZ+fHxUV5efnZ2Njw2AwWCyWoaGhm5tbeHj4ypUrY2Nj9+/fX1NTI73uk9ytMxgM6ZZpNNpb3whFR0c/ePBAT09vwYIF7X1l7YbNZpN1nMePH2+0QmZmJlk1n5eXR1FUdXW1jo7Op59+6u7uvnPnTjU1NbFYvGjRIhMTk5KSksjISHd3dyMjozFjxgQGBiYmJorFYiaTWVxc7O7uvnfvXk1NTZn2b968KVNAH4KujuOokwwdOhQAHj58SP7o6en58uVLiqIOHDhw8+ZNSTXykJTNZjdsQTLsKiwsFIvFMp+SnCbbt2/vkN63n/j4eADQ1dWtq6tr+Kn0aPTVq1clJSXSn5KsV7W1tRRFcbncsLAwOzs7BQUFyf9NWlpajo6O3377LQCMHDkyNTX1119//fHHH1euXDljxoyAgIBx48atX79+/fr1w4cP74SLRZ0Dw2hP8fHHHwPAgwcPyB89PT3d3NzWrl1rY2Nz8+bN+vp6cpw8QlVVVX2nxv/8808A6NWrV3l5eTv3uwMMHjwYAC5dutTwo1mzZpEfkZOT0+3bt6U/4vF4AKCkpCTzFR6Pd+XKFQ8PDxMTk+aHLFOmTHFychIKhUKh0NHRscOuD3U2vKnvKchDPbFYLDni6+sbGhrq4OCQl5dnZGS0e/fu+vr6Rm/q32rLli0AsHLlyl69erVrrzvE119/DU1MIN28efOVK1cA4NixYwcPHpT+iCzYl1kBtWHDBiUlJTs7Oz09vbS0tMzMzJ07d44cOZLBYNDp9AEDBtja2i5fvjw0NPTSpUtkgxY5OTk5ObmuytWCOgK+qe8pZMKogoICOcJkMmNiYgoKCjw9PX/++WdyQ/pO76AzMjKuXbumpKTUXeb3LFy4cOPGjdHR0eT3h/RHTCYzIyPj4sWL0r9viPLycmiwaoDL5ZJCbm6uSCQaMGDAgAED+vXrN2vWLHt7+6ioKJlGJBv/ddMdAFGjcDTaU8iE0W3btunr6wPAokWLTp8+fe7cOQ6Hk5WVRWZcCgSC//znP48fPy4qKnpry1u3bhWLxYsXL+4ui801NTXt7e1FItGJEycafqqsrKyrq6urqytzvGEYra+vJ5PAZFJTNzO9VJIUppnsMKjbwdFoT9Hwpl6CRqM5OTl9+eWXR48e3bx5M8k3OmPGDEkFdXV1fX393r17N/wvjUY7d+6cvLw8mVLaXSxevPj8+fNHjhzZsGGDzNDbyMjIyspKJBL9/PPP0selb+rj4+O3bdumrq7eq1cvksiZPCggWrjYCX0wMIz2FKNGjWKz2aqqqk1VkJOTc3Nz4/F4a9eu7du374ABA169elVYWFhSUlJWVlZWVibZwl4anU4Xi8V2dnZkIWZ3MWXKlL59+5bl5t67evWTadPg3wzTurq6lpaW5eXlvXr1kklxT0ajAoFg8uTJt27dAgA1NbUJEyY0bBzDaE+DYbRH+O6773bu3Emn0w8dOqSkpNSvX7+mapJMmgEBAZJ1Ph4eHkuXLhUKhYcOHVJRUWGz2QUFBfn5+QUFBSTO0un0mzdvXrlyxd7evpOup83odPrSOXM8Q0IUPv+8ok8f+bt3s8vKTp065eLicuTIETs7u0WLFllZWe3Zs4fH41VUVFRVVZF1seTno6amtnLlytWrV0sy/O/Zs0eS9gXDaI/T1VMFUGeYM2cOmdLk7++flpbWVDU+n6+oqEin0wsLCyUHV6xYsX///pMnTy5YsODOnTvS9Q8cOAAA5Bkrg8EIDg7uuEtoOYFAUFZWRlGUWCwuKipqqppYQYECqAOgABxb9j+LmpqagoKCn58fab8pU6ZMAYDo6Oj2vzb0XsLRaE/h7e1No9ESEhIcHZsMGrGxsTU1NaNHj5Z5wWJubq6urt5w9xEyJSg0NLS4uHjt2rXffvstl8v96aefWpdbvr2kp6dHR0dv2LChtrbW29s7PDy8YR0ul2tYV8cC4AP0AmAwGAry8goKCsbGxk+fPq2srJwwYYKJiQmbzVZSUiIJBgsLC7ds2WJubu7v7998B0aOHAkAzS/GRR+Uro7jqDO0cDTq5uYGAIGBgdIHV6xYQcZ08+bNI6PRS5culZaWkin3WlpaZFVPZGQkSYwyZcqUiooKiqIkS4BKS0sbrnrqOCkpKb6+vmVlZQUFBQsWLGhYITo6Wl1d/RwABUABFNLpYi732bNn3t7eFEWRRFYNF93X1dUpKSnRaLTi4uJmzr569WpS8PT0bKcLQu87nPCE/uf69esAIP2OHgBmz56toqLy5s2b58+fm5qaAsC9e/eqqqpIwtDFixeTJZIODg63bt3S0dG5ceOGtbX1ixcvli1bRlrw9vauqqrqzAtJSEjYu3cvedV+5syZhw8fkuMURYWEhNjZ2ZWVlZ2ZMePl5cuTe/UaJhbHFxTQaDQyJZ7M8Wo44UleXn7UqFHUv8tJm5Kfny9TQB++ro7jqDNUVlaSAo/Hk6z7lFFSUjJ+/HgjI6OGI8cHDx4AgKqqqq+vr6+v77hx41JTU8nQ7O+//5auyeVyST5TfX19GxsbHo/H4/FcXV3J+LRzpKSkbN26laIoPp/v4uJCJidYWlqGhYXNmjULAGg0mo+PT2pqqoaGxsSJEwFAMmitra2l0Wjy8vKNDp/JnPl169Y1c/YhQ4asXbt27dq1VlZWHXF16D2EYRRRFEVt27aNFL7//vuGn5KF9gwGg7yzXr16NQkoU6dObVi5tLSUTANiMBhOTk6bNm0aPnx4Z4bRnJycM2fOUBQlEAh8fX09PT2ll6iyWKxjx469fv2ajKynT59Op9MVFRWTk5Pv3LlDdkDR0tJq9JfNjRs3AGD06NHNnF2yWF6S5QR98DCMIor6/8xGjVYga3K4XG5UVNSgQYMUFRVJdv1du3aRTFESsbGxw4cP7927N/z7+n7x4sWdGUYbOn/+vPSeIkwmk8yiV1dXHzp0qGTSkuRTMpr28fGRGWtXVVWRFfFVVVVNnQvDaA+EYbQHKS8vJ9nhyL229Efm5ubkVnTChAkNv1hfXz9o0CAAaJjGGADodLq1tfWhQy9jYlLIXB8AMDIymjhxIlkgNHLkyK7K/CQWi4ODg0k3nJ2dY2NjHR0dGQyGhoaGzFXIyclxOBzyhNTAwIAcpNFo48ePP378eHV1NWlt8eLF3t7ezbxlEgqFMgX0wcMw2oMEBgampqZSFHXu3LnIyEjpj5oajT558sTHx4fMDCUGDRrk5+d3/vx5kpNYS0uLzWazWKry8m4ADABgMpVXrfI7fFiQlESdPXtWXn4um63PZDIHDRo0f/784ODgqKgoLpfbCddbUVHx2WdLAVgMBmP9+p8ljzu9vTcCgIkJ5+OPP5bMmZekblq6dKlYLL57966rq6vk14aqqqqHh4eTk9O6WbMoR8d4Q8MfJk/uhEtA3QKG0R4kMDDwr7/+qqysDA8PP3369I0bN7hcLhk0yYTR9PR0Pz+//v37S6InScnh5uYmFouLi4szMjJOnTpFFuqoqqrKy5P7YhbAGoA3dDo1dy7l4EAJhZSFxW0APZmhn6KihqYmZWlJzZ9P7dpFxcZShYVUYCC1bBlFUdSPP1L377f+Mrlc7rlz59atWzdkyBCACzk0Tb6+/j3juTev161ZQ40eTdHpKwHAyOgz0hlTU1MnJydra2tJ8rpPPvkkIiJCKBTy+fyIiAgbGxsajebg4LBy+vTHurpkmlQdnU6Vlrb9LwV9AGgURTW8TUMfpKCgoIKCAh0dnbS0tBEjRpCcePLy8n379jU0NBw8eLCJiQmXy713796jR4/IV3r37u3s7DxnzpxTp07t3btXTk5OJBLJ/Jths9l8Pt/BwcHPb9vTp6aRkfDsGZibg7095ORAWNi1rKxvdu/eMGzYsPT09MePH6enpwsEg//66yeZ7q1eDeXlsHQp/PknWFnByJGy/X/1CoRCMDaGoiKg00Gy2DI/H5KSICcnJTZ2Y2JiYklJCTmupqaWVKFgShWSP/40LWr5dXsAoNG+oKhfXVx++PJL45MnT/76668AMG/ePH9//yNHjhw6dKi0tBQA9PT03N3dly9frq2tnZWVJRKJfnB2PvbkCY/BkBeJmABQUAB6sr8hUE/U1XEcdR7pm/pt27ZNmDDB0NBQJn8weZWkpqY2f/78qKgooVB49+5dKysrACDzQwFAU1NzwIAB48aNmzlzpq2tLQD0799f+kRiMbVkCVVbSy1aRCkpXQbQe/LkiUxnXr2ibtygdu6kFi+mRo2itLQof3/q8WPqyy+pHTsaH40eOkSNHk0JBNQvv1AHDlBBQdTnn1MGBmR0SI0cWSS5hBkzZpAtoURAowAqWVoUwMPAX/39qd9+o4YPtwCA+/fvUxTl4OBQWlpaWlrq4OBABuY1NTXh4eEfffQRaY3JZDo6OsbHx1MUtdrWVkSjkfMJGAxKakloQkJCXFwcRVF5eXmnT59u37849J7DxaA9SN++fcmjQB0dHV1dXTIara2tffr0aU5ODnkL//vvvxsbG6enpysqKj558mTu3LmRkZEAoK2tPWvWrPPnz5eWlhoaGl65coUsduTz+SoqKs+ePaurq5O88pZE5hUrio8dm6CiokImk0rr3Rt694bPPvvfkaAgUFCAxYvBxwfGj2/8EuztITQUTEygoAACAv452KsXjBgB1tYqPj7nR44c2adPHwB49OjRiRMn0hj9PhI9VRG84dPYQ2eZDhsCBQUFublcADAwMDh9+nRlZeX+/fsB4O+//yatKSgoLFiwYMGCBX/88ceePXsuXboUGRkZGRlpaWmpoaHx88SJy83M7t2796ex8TImU1BWVlFRUV9f//z5c5KEsKqqStIU6im6Oo6j90VWVlZNTY2xsTEA7Ny5083NjbzgVlFR8fHxIRP4MzMzORwOABgaGj5+/JiiqD179pCv7Nu3jwzHEhKoefOou3cpkYi6cOECwBgbG7uWdGDvXio3l6IoasECKiGBOnpUtsKhQ1RcHLVkCRUSQp06RXl6UqdOUVlZVKMLTV+9egUAxiprqLAw8YGwEKekkhIqNZW6f/8+jUZjsVi1tbVfffWV5KHw3LlzG323zuVyvby8JNmahw4dqq6u3nDGQnBw8IwZM9auXevq6urn59eyHzn6QGAYRf9YsGABn88nSzxJmFBUVPT29pbZHbOkpGTixIn6+voeHh5UUFDc2LFrDAwAYNmyZeHh4SEhYWpqJwAi5s6NiIiIIKnzfH1936knIhE1diwFQO3a9c+R0lJq82YqLIyKi6MKCigTE+rXX9/SCJ/Pp9FodDp92bJljo6O48ZNYjKHMRhGnp7fAYCGhoafn99XX30VFBRE6m/dupU89m1UTU0NSbkinf1eQUFBXV29T58+pqamP/zwA5nzT97OvdP1ou4Owyj6x4IFC5KTkx8+fEjW/EyaNIns1d6QQCDIzs7+ZcgQis2mAMQAUwBMTU0DJLfZUtTV1Y2NjcPDw5tahNqos2cpOp2i0aiTJ6nUVMrMjAKgXF0pMh0+MpL6d6PoRpD1/gMGDGAwGJLnuRL6+vMAwNrauri4+Kuvvmp5l1asWAEAn3zySUlJiWRxrcSlS5cwjPZY+GwU/U9iYqK8vDyHw0lKSnr+/HnD9BwEk8k0MzPTzsmB2to38vJaQqERwGtV1ejoaBqNbmg49eOPlUn4EgqF9+/ff/bs2cKFC7du3bpx48a5c+e2JI3enDmQlwfr18O2bcW5uZo1NYwRIyAoCAwNAQAcHBr/1osXL/bu3Xvo0CGyepXD4Xz++ecmJiYaGhpaWloAmh4evaqq5gGAiYnJO21imp6eHhYWxmAw9u/f33DqPgBMnjyZoigAMDU1Xb58ectbRh+Cro7j6H1Bbuopivr+++/J487Dhw9LVxCLxaWlpS9fvszJyUlKSipisah/c83N/3cNJZvNfvj/A8W6urrw8HDJFFRjY+Ndu3aR3HrNE4vFU6asodFmA4inTy/7/1VXsh48eDB//nx5eXlyFisrq/SUN0QAABI4SURBVIiIiIbj3zt37pK9o21sbPLz85OSklr4w5k6dSoArFixooX1UY+CYRT94+DBg2Sp6JUrV3788UcAUFFRGTRoUL9+/dTV1RUUFGR+AX8CkAiQBBBGVi8BMJlMW1vbRp8wikSiiIgIc3Nz8t0+ffrs2rWLRO1GlZeX29nZAQCNRgOw1tfXz8zMbLTZqKgoGxsb0iyDwbCzs0tISGjmMklGEhJwmUymvr6+paWlo6Ojj4/Prl27IiIikpKSXr16Jf2VixcvAoCGhsabN29a+tNEPQmGUSQrLS3t448/hn93B5FGtgg1MTGxtLQcM2aMjY2No6Pj3LlzycJQMtBzcXFparBJoh5pHAB0dXWDg4NlVvfX19dfvnyZTFpis9kzZ85UUlKSl5dXVFS0s7MLDw8n9ckDUElcVlVVXbVq1fPnz5u5rpKSktmzZ5P6Ojo6ze+VpKysPGjQoGnTpn399dfk4cb+/fvb8YeMPiQYRtH/iESinTt3koGnmZnZ2bNnnzx5wuVyS0pKmhk5UhQ1f/58AFi7di154Dh27NhmkneIxeKoqKgRI0aQgKWtre3n51dRUfHy5cvg4GBDQ0MAUFJSkg5qDAZDskxATU3N0tJSssVpv379goODm98ciaKoW7dukZbJLC6yDKG2tpbL5d69ezciIiI4OHjVqlWOjo6WlpbSr+PJj4JOp2/duvWd3pKhngPDaE93+/ZtMqXp8ePHkt1A3dzcmskF11BISAgArFq1KjU1lQwkTU1NG70NlxCLxZcvX5YEUyaTKbOeioxGZbLYSWOxWDNmzEhMTGy+b0Kh0M/Pj4yUR48eTRYgvDXjVHl5eWpq6tWrVw8cOODw7yutUaNGNVyOhRCG0Z5u06ZNWVlZFEWdPHny4MGDJiYmly9fftdG/vOf/wDAxIkTKYrKz8+3tLQEAA0NDZmdRGVkZ2f7+fmR5afSwdHOzi4sLIxMtxKLxVwu9+LFi/7+/p9//rmcnJyCggKHwyFDS4IknWo0a9TTp0/Hjh0LAHJycj4+PoWFhQCgoqLyrhcYHR1Nfj2wWKygoKBnz56R4+RHh3o4DKM93aZNm3x9ffft27dw4UKyiL4VjeTl5QGAlpYW+WNVVRV5QcRisRouMBcIBJK0SSQOmpubT5s2jQQ46b2dG57FxsZm9OjRJOVdUlLSqlWrJFGYTqdbWVnt2rXr9evXpH54eDhZR6Cnp7d06VI7OztlZWUWiyUvL+/m5vau40oej+fj40On05cvX/7tt9+Sg5ibGVEYRtGmTZtu3rzJ5XJ37twZFRXV6nbIbEpJEKyvryfTJ2k0mmQ6emZmpo+Pj+TdjoKCgqOjY2xsLAmLkydPBoAdO3bItCwUCr/77jtSXr9+vcz9uEAgiIqKcnZ2luS3Z7FY9vb2klX8DZ8VSI5PnTo1KiqqmcVLDd29ezcjI8PNzS0jIyMjIwPDKKJw+j0CAENDQxMTk+bfXDdPLBbr6uqWlpbOnTt34MCBysrKbDbbyMho9uzZFy9eDAgIiI2NpdFo9+7dI/UtLCyWLFkyd+5cyZsiAAgNDY2OjiabPOfm5rJYLLITiVgsfv78Oanz7NkzmWnzTCbT3t7e3t6+pqbm6tWrJ06cuH79+pUrVyTrlxgMRv/+/TkcTkxMTG1t7cqVKz08PH7++ecjR47ExMTExMQYGBh88803bm5u2dnZn376KQDExsZ+Jp00RYq1tfXTp0/z8vLu3r0LABTmmUSA0+97vOTkZLJRUk5OzosXL97169JzmGRer0sjN9eqqqpubm7NTHrfuHFjdnY2RVHh4eFXr14lBwUCwcCBA8keJzY2Nm/tUn5+PtnyZM6cOQkJCWT2Fcmbt2TJEkkO/NLS0tDQUDMzM9JDW1vbJUuWkI+aH2NyuVy8qUfScDTa01lYWJACmZfecvX19WfPnt2yZUtWVhYAGBkZiUQiPp+vq6u7ZMkSABAIBBUVFbW1tfHx8Tk5OTNmzDh37pxk046m3L17NycnJzU1lWwvSlhaWoaGhgKAs7PzWzumr6+vq6ubnp6upKTEZrNZLNaFCxemT58+ZswYd3d3yT2+urr62rVrPT09Y2Ji9u3b5+rqGhYWFhERAQC1tbXNtK+iojJ8+HBS/uSTT97aH/TBwzCKWqS6ujo7O5uMOm/fvp2Xl7dlyxaSWNPY2NjT03Pp0qU1NTU2NjbJycmRkZG///673r+Z4ffs2bN69ep+/fq9NYYCAJ1OZzAYZHZnq/H5fACQk5OjKAoAKIoiqZ4aPZ2tra2trW15efmZM2f69u0L/64jaIq2trYkmpN8JaiHwzCKWqS4uDgiIoKE0R9//PG3336rqanp16/ft99++/XXX5NsIywWKzo6etKkSQKB4NGjRySMPnz4kDwAraqqasmJrKyszMzMCgoKJEeYTCZJ3wcAx44da0kjJIzKy8uHhoZqaWllZ2eTBQLN09bWHj16NPy7VBShFsIwilrq6dOnZH4oi8Xy9/c3MDCYM2eOzMBNW1v7t99+E4lEa9asIek8tm3b5uTkBACVlZVvPcWECRPIG/+PPvpI+kmrpNzM41dpJIwymUwvL6+hQ4eeP3+eetu7IEVFRclS0Xnz5rXkLAgRGEZRS0mn71y/fn1T1cgK9Jqamr/++gsAysrKWj4albwfl6y7bx0SRjU0NMgiKCUlpbeGURaLReI+AMycObMtZ0c9DYZR1FJ9+/YlUzslt9jNEAqFZJYS2awJWjYabS8kjK5YsYKsjp8+fXqnnRr1QBhGUYuoq6tPmjSJlB2aSpssRVVVldzLX7hw4Z2ejbYLEkZb+AQAoTZq0/tQ1HOoq6tL7rhJfGzelClTJAUSRjttNCoUCoVCoZycXDNpTRBqR7S3PjNCqI0qKirU1NRUVFQ6J5KS06mqqpKtRBDqaDgaRR1ORUWFRqPxeDyyk3tHwzt61MkwjKIOR6fT2Wy2WCzm8XidcDoMo6iTYRhFnaEzX9ZjGEWdDMMo6gyd+bIewyjqZBhGUWfA0Sj6gGEYRZ1BZjRaWFiYlpYGACT/U6Nf4fP5AoEAAOrq6lryUJVk4AcA1tmzsQDb09Ph2bP26DtCb4FhFHUGMhrNzs7Oz88vKytLT0+/ffs2AFRWVp49e7bRr5w5c+bOnTsAcP/+/aNHj771FOvWrQMAOH7c+MQJG4CRxcViV9f2uwKEmoRhFHWgmpqaK1euuLu7x8bGKikpLVu2zMDAQENDIz4+/tGjRxcuXLh27VozXxcKhQKBQCgUtuRcYrE4Ly/vt4MHewsEhTQaADz4448FCxZcuXJFJBK1z/Ug1BgMo6gdxMTE1NfXA0BiYmJxcfGbN2+OHz/+xRdfaGpqzpw58+DBg3w+n8/ny8nJ6enpqaurKykpqaqq6uvry2wLKuPUqVMBAQEkOZ6fn19SUlIzlV++fGlhYfFrQgIA6FEUABSIRCdPnpw5c2a/fv18fHweP34MAJJGmm8NoXfQpbn30Qdi4cKFfD6foqht27atWbNGkj2PTqePHTs2JCQkIyODpNkPDg6mKOrOnTt79+6lKKqoqGjFihWNtnno0KEbN25QFHX37t01a9aQBgcOHBgcHFxcXCyp9vz58507d5I8oQBAA/BXUKh0caFWrvz7r782bdrUr18/yb92GxsbR0dH8kXc/wO1F0xNgtrHyZMn5eXlk5OTJ0yYQKPRrKysHB0dHRwcDAwMSIUffvhh8uTJW7ZscXV1NTY2Jm/SlZWV7e3tG21QV1eX5GdSVVXlcDhr1qw5depURkbGt99+6+fnN3PmzCAVlf4XLypVVUWJRH8BMBiMqVOnWltbs1gsuWXLQFGRAxA0alRQUFBycvKJEyd++eWXwYMHZ2Zm/vHHH9C5GafQB66r4zj6ECxcuDApKSk1NdXT0/P+/ftkj7yGbG1tAWD58uWtO4tAIDh//vyMGTMYDMZ6AD4ABUAB3AcwMzMLCAj46aefmvl6XV1dWVnZpEmTLl68ePHixcmTJ7euGwjJwGejqH0MGjRoyJAhurq6cnJy0tsmSwsNDZWTkwsLC0tPT2/FKZhM5uzZs69evZqXl2dvY6MIUESnA4Bar14hISFN7UcvIS8vr6ampqmpOWvWrFmzZmlqaraiDwg1xPD39+/qPqBuT1lZ2czMjE6nKyoq9unTp6mp79ra2q9evUpMTHz58qWLi0szDdbV1ZEHrJKCNBUVlT6qqlREhDJFAYBg8uRR339PUdSbN29GjhzZfFd5PB7Z15PP50s2+ESoLTBRHupUxcXFHA6nsrLyt99+I7n0G+Xs7Hzu3DnpQkPn3d1jDx78ePRom+PHzczNc3JyCgsLra2tO6rrCDUBb+pRp9LR0fH29gaAdevWtTFvXrqBwUGA/ClTzMzNAcDMzAxjKOoSGEZRZ/Py8urTp8/Dhw9Pnz4tfZzH412+fPmbb76JjIxMTU318vLy8vIqKipqqh2yQpTNZnd4jxFqFk54Qh3u9OnT5ubmlpaW8fHx+fn5Dg4OgYGBrq6uGzZsmD17dnV19fXr169evRodHV1dXQ0AJSUlH330UWhoKAA4Ozs31SyprKys3GkXglCjMIyiDldWVkayLvF4vLy8vJCQkPLycl1d3VevXnE4nIKCAvKAnszVt7e3//zzz48cOUK+a2ho2FSzOBpF7wkMo6gz7N+///Llyy9evBg6dKivry85yGQy8/PzlZSUJk2aZG9vb2dn17t3b/IRGYpKFxrCMIreExhGUWdYvnz5J598Ehsbm5aW5uPjo6qqyufzt2zZoq6uXlBQwGKxWtEmhlH0nsAwijqVsrJycHAwAFRXV2/ZskUgELQuhgI+G0XvDQyjqMN5eHiQJUY2NjaSecpsNptOp9fU1IhEooYT7FsCR6PoPYETnlCHo9PpJIzSaDQ6/Z9/cjQajc1mUxRFBpWtgGEUvScwjKIuQ1Lit3qfOwyj6D2BYRR1mTaGUXw2it4TGEZRl2njrss4GkXvCQyjqMu0ZTQqFAqFQqGcnByTyWzvfiH0bjCMoi7TljBK7uhxM3r0PsAwirrMO4XR0tJSSeHRo0fjxo2Tl5fn8XgjRow4ePBgSzayR6iDYBhFXeadwmjInDkwYgSMGBE+YYKlpWVmZqZQKGQymcnJye7u7kZGRqtXr05LS+vgLiPUCAyjqMu0PIzGbt++Mj4ekpMhOXl1aqqiSKSqqhodHf369evDhw+PGjWqrKxsz549Q4YMsba2PnHiBFkoBQAhISEdew0IYRhFXajlYfR5Wpohj1eroCCi0egA7kuWPH36dNq0aWw2e/HixX/99Vd6erqPj4+Ghsa9e/d27tyZkpJCvvjgwYOOvQaEMIyiLtQwjAYFBZHC999/L12zSlm5hkZTqK1lUNRrFmvH3r0yG9KR/etfvHhx+PDhgICAly9fkqzPGRkZHX8dqKfDNfWoywgEAgCIj48PCwuj0Wh9+/ZNTU0lH0kKRIm6+s+zZ68dMSIlJSWcx9vVRDYTMjgFgDNnzrw16zNC7QXDKOpsJSUlFy5cOHHiRHx8vIaGRlJSUlJSEgDY29uXlpZ6eXkBQMMdmO23bgUOJ+3Uqcm9enVBpxFqGu4MijpKUVHR06dPx44dW19ff/369bFjx164cOGXX36Ji4sjm9mRdZzV1dUWFhYjRowYNmzY77//HhERAQ02BM3Ly9PW1lZQUHjz5g2DwVBXV2/+1C9fviRp8yUFhDoOPhtFHSU/Pz8uLg4A6urqTp48qa+v7+7ufvv2bSaTOWvWrHPnzhUVFd24cYNOpz9+/Njd3d3Dw6OppoyMjBQUFABAS0vrrTEUpLYewRiKOgGORlFHSUlJCQwMnDp1qlAofPDgQVpaGpPJXLBggbOzcy+pG/OVK1fu27fP0dExIiKiqKhIV1cXACQFhN5/+GwUdSAzM7PJkyfX1tY+ePAgISGh0fTM27Zt09HRGT16NADo6urGxcV9+umnGENRN4I39agDaWpqcjgcU1NTAGgqxb2ysvLmzZslW4H+9NNPndc/hNoDjkZRR+nTp4+trS0AsFgsV1fX5iu/fv36ypUrAFBRUdEJfUOoHWEYRR1FU1OTTJJnMBiffvpp85UZDAaZjS/ZZQSh7gL/yaL3goaGxoQJEyZMmECCKULdCIZR9F5wdHSUKSDUXeCEJ4QQahMcjSKEUJtgGEUIoTbBMIoQQm2CYRQhhNoEwyhCCLUJhlGEEGqT/wKyPX4WgffylQAABsB6VFh0cmRraXRQS0wgcmRraXQgMjAyMy4wOS42AAB4nO2Ta1CU5xXHzwK7LAt7Ybnf98JyZ7kGvMD7nAcQElAkipWAEoGAIEVRREhDEfCCYVBhKgarEZDMapoqN6WiZd832jaUQKSo1ZGIREVDMFxCDGoFu04+9Eu/9GNncmb+85s5c+bMuU7rz46CwaTwHys0aKdBFTxThfC1479C8Av+Z/B+wf81zE0A+CYg4IOpAEwBhKZgJgQzHojMwFwEFuZgwQOxBUjEIDECqQRkUrA0AksZyC3Byhis5GBtBTbWYGsDtsZgZwJ2tmAvBgd7cHQAJzE424GLHbi6gJsQFG6gVIBKCGoVuKtBIwcPDXh6gJczeHuBjxx8fcDPFfx9QesEAVoI5EOQCILNIUQKoTJ4wxrCbCDcHpbYw1IHWOYAyx0hwhEiHYFxAuIE6AzUGaJcINoFYlxhhSvEukGcG7ypgLcUEK+EBCWsVMIqFSSqYLUaktTwtjuscYe17pCsgXUa+JUHrPeAFE94xxNSPSHNCzZ4wUZvSPeGd71hkw9k+ECmL2T5wnt+kO0HOX6w2R9y/SHPH7ZoIV8Lvw6AggDYGgDbBGDNN0yYDyZ8gSmfZyYyt+AbSaQyS76xlbWNrUxoxrMQS4ws5VbGtoa5mdrFG8J58LOgsCStGj1nmhjhrgz8MDoAa5b2MLtvcTi6Uod/aNYwxGkKh0ZuIM/rDJP4yUl08TyBUaKb+vVn3VF3vtLg92JbNrPkOhuKxUeyejvF40Qhi8TC70zJdF0UFvu9IPtbO8jdwl1YEHCLhHSGkUWqw68uzjHzgQ3MhOYi5pMkMsZd1hd530OzagvUaE7oB54/QCaJYlZmFfvAy5jWxB3EpI+m9bPDAtqYqjPEOLO6EXuabHIN23lBvSWF1tT/9uDr2piQ1le4z3QKk/YcJxFd4zgY8xN2p8jJ7h2XcO0xMV2T2BSpLTuDlvlC+unEU72fUwXuXTaOFVE6/dTFYjwU0I/XOgpYx6P+6NIuwSe1ZWzZLnvqUbNI6rP7mAyFir67R0fOhZ1l0s6FU5f0WFyiU3CbbecxKboPB7a5chM1DDL93+ClZ9GGPGdJIO8x+u73JJ8sppIIlSuRXLhCEo1WoHCLgqhT7fFO62PiILVm7yXyMSmrgOQEeqIyKBjPOVeixZOHxGRLCeb6foQjTxT4vVSHUq4Rn+da0KYHrUQie4WFq8U04vciTK20omlrFDRFP0H4nUo6udBh2OlSTjRTjUV32/CrNw5y+e/VYpjn12gx7o29q8ewNegZrqvbioKYAeQigWboIsmPPmP6gJb7mKA/heUtgzifuwXN0jO4fUnmtOzlVoyt/5irvmREJ8JrcKy7nh35IJTIbP6it/3mH8yxI3/EWW0O80gZx75c+AEbFU14alsla2TtQZsDH6M8IT7ikNqVWv7TiX7p9Chy+NZ5rFkrpzXRWaykMwWDeZ1EJ+1if1NkS9UpWqySs2y5VE1XJoTitpAUskyvpoKvu8n5lkxyt9ue4tWj5PiqE+S0MJIODc0xhhp6ZQfC6PzsFIlXyfQFhymtEkgwbUUwV/p8Gn1z8vDMTSW3d/soNvzQg9OfLeE2ukkRPMdRJw3hQjsInhxeQJFNPauf7SGS4WH89n4z+y9hELkf9yUmdlWRm8WB5PrIj3hRfdRw51fIgzvdTF6fiiyYh2DHDU3v8m+nSeWSPIzpdCZX3Ndiz1Uj5C9eIL9TivGLHB15v6+dHW+Tobqpi+gXJOzfrLwwvFbC+MfJWcNfkOmXXuSnle9gQboPDt6uQ7uwUGwXxKL5+Fu4o7SEDE6W48MN1fgwv5+siMzGTNcu7C5Zh9rzp3BYfAZfZWsxvKQNx4pKcTVvkQyJLmNmsg53fCGl8oR+ZvlfjekJYzs6+mKAFM1N4tA+Cb23MQZncqX0yS6gt2uMMPugE3V2daOq1lySXuVBS7X+lJTPk70JjjRYq6DpW+0RvbVU+fbnmOyg5uo8mnHxzw34frma6y88gM9ftOCG1iOcoV98VNOMfhcaucKZPFwVNIDdtY3cseO1eKj9CjbMatB66hrOVdxAw+1h6dWXeGC5Gd2QvB3/Xj+Czb4LKFm/CV2OtWHWbyU0b3iAeba9QX8zZw7LuB1MVedhdrGWR1Unr5PLU+EsWzGKm0SfoZv3JPoNPsWeuj+hpfJzPD3ThxkFp/B2FItvxtfj051J3Ot+938Xi4cjMjnNHTtaVxqOAuVpLmZyEffFf4g60UkuKet7TC/eg6q+Nm53m5g27jTGgtzkXt4Hx4ntvwEmYo8wsia8sAAACqx6VFh0TU9MIHJka2l0IDIwMjMuMDkuNgAAeJx9WV2OHrkNfJ9T9AWmIZKSSD2u7UVmEewYSJzcIe+5P1JF2S1+WCD74x3LbHaRLBap3reLf/3j29//89/r+cu+vb1dV/s//6y1rn9ba+3tz4s/XF9+/9sfn9fXH799+XXy9fu/Pn/88xJR/DvxFP5+tf7tx/c/f53I9Xn1W3UtGVe7TdR8XXYP1dB1nlTY6e0WfTntpvj0a95Np7V57Ax245Y1Iwx2Ars+L797b6HFrsNu3nNImNKfRufpvYaP4cduXF9xGr7cx/Xe7i6zrQ4ow61FCW2mYRtmTWk4x4yJ+O6+ho1i6DDUO7RrCxq2McYUGHqsXj1GGq7RpSER97C2VJCa1ltbJZQFO6FdD6T8tvBoBrsYbVpJobQ0nCFizLWLTYTKbAHCKIYCQ4TqswnfrEssFgx99KlRDDUNQ6zPjMXmdFRj3guVkhKLWAaDN3nP9FhDIQWFESSiW7HssESQ2oba9S53k5jqV+BhZ4THctemS1+oNWszfQ2BZdfZV42cxWGcMcbK4oDIwLHuMWxqDZ3VmfC01Hd1uikgrVvdeyWGxLY0zST/IobfK2J6hcn6gICgLp6RG+5WH4zcQ14o3tJwrhErC+4y5wRuEBe1KIasT9xohWUkeXP3oKGCz6N61DQcOsUybvBLV5CnzUxKNGr5bplmK3Np2h39AhbZ6lLC0Z5xGxqnZUcstE52hAt6tvocWXP0NotBy5gDh8i6BvJbLOemcMRYwm5cAo6DsLeYrlbqo56WUKCITksDDGRVGFqb1WdsHqHQQPeOhpUmmXdHBaL6XBmRjlCkBJadEmDsPXetfLddIXTssmSHzzHADgBCC1cWm6SlSsRir83RLUMDA9GnxXC3EGTKhWoxVARqAeCOsldhs22IQ9BI0eZtNPZ7NIpdMWSB3kmfMHS13hC+lRRFznWVFNku0ITmIgcoYKAVF9lhSIEUw10fR85HJyfmCKhqRyo7+rgYbnFDTvogIQWS1jsN3XBUDCOZCVGloKApIFmR3ePNq3DY2lw3uqKCQFdn9qOsF/3tqW6NwiGQIGHXCDqFJGpjzULhLikHeLUiChQcrgXQMYfAohp417SEhDslhmKjIWwgqCwIUywtScRBgR5QFnw2sg1ajdiLYU9mIHCgQi5dnaKEjgqMjkKhPjLyBjqw+W/nBGvZxQZu1HfPTUvyl1lvMVBoipMi8Ppuz56YIBHK884pixnQWQhk/SWc3T1ADllj4EttBS3FMeNKyVH+D+ZDlqG3odeYFdQQAAZBpdQctP5IaQEfnQYdzYmIBvrdpPpECT44d6BBsfuMAsdpgcYbrUQ0NN++nNTe0wJFX6w69oiobT4Mpk6nGMxUTdAeCVk3NoVZmQSKf+DYbTlZc8OjITeD24JVEqO7aAhVsZVSGIEhwaJajJglncj2B1ViNeSeMJFMLjwLOUDmCjcxYj7YitplB9TRtYiSNF4gXzWNNA1UhOsA5gGADFoqVLxXy5Wv7wp9bbQMKHCjEPc+e0nnbGkIMQqjrCnUSFJeQ/qsY2DuEmEu6d53pLF5UolJxBL81HQa2A329IW4g/7YHqDO0JkictMyJFRQekqxx+joPIxPUqV0x2SRQGUWU1LeNTgtsXqgSnVQT1YJQ1RaMPfK6PfAUsz+KLSbMxmCnREsyrejmTkg7gWnq+jC9LQEbRfbg6QHMwZzjGG76lrIIs17BBKYldeB3SJHAghQR8ZcaWkYkUpPTbH3cdedGJar5NNbvjx3157V0jVHzo7W6mzxXSNkI7vrnqAxzjKu2UqKXHdnYrPskdqJLQ37LvOLWhWCuKVLUKEt/jn39hxH2JOkLqXef1JOQjiEsM5hE9+ju1UiOcvDjgjZnhTQjCDQGlgdiuXcKMfINU45C5VhIPtLX+Lxn5YMfDFFmNMjc4mVtpWSeySN0L+GZUkpjQ27aufO1GfpdV8/m21EyyUcIwMNaLdsoT9be8sMLZARzWpYaoyVx2zh0lJeHZIZGiCG0BF+wFzlXMKeW7e40I3RQBdNee8YSFz30Dt14Qr7pR1oK46tPVEwrzE46koauzgQdATB3uiYMS03ZGd4x5DFCS5ukpKOPEnjFIAizjooY6YUYsm1vIQgzb3x3mWed4Nj6NsjZtqkUjZkicdImInUqCMNsa9iVSXBIMfKwY+VDDO+GK58NRaN4DULrgPSTRnsdF7evVLdIH+KXYMNs5hT6uDs1q1eqWRbcpr43rYkSZnayVwdU822bR17294zO9YoYZ1a91ZfbzmrcEkxznFYDtSZ2xGuEF6TtHpKVsOmh8Ri+qINwTeWHN3TK9Atbk1waxi0xMYKIOxcXG1axTm3vEAHQTdaTofOc3RjO6x9tnxHFB1XMG5z2PARMdaJxRtLMdzaNikRlBWUkGPFsU8FVfkYrqzlcMdmybUC6wbEAOvSZPvUK+rWNojKSP0ZWDbQj+Syg6zVUn4SaeRcQrGRmby64Lr6MvgxmbZg4jrMEQLy4pKaYowby7QXU9usWxDp5PHg/S5XEAlqYzHtGT2VQ7h1YQCgDKQ8sFaGYP3c7+d1j+0GJmGtYooh9rWBubjSp4Iizps55iRGNX6AYK+X+3TzJNPg/uBc5XSr/F8bCXMv9QP/bZia3PoMbJps18UXVNOfKodRjUK9g28K/nPxgp7pi1dpCYDtyVGNHcwCVzCOaswlrzd6kTQlAET1zn0a6zybHkJuL/UXTQBcMdGD+XVEsM3ua/OsI/j3z28vH5L2p6Uv3z+/nU9L6JhLzhckuVCn86EIYV3zfA+SKy4/n3342+dZdP6F3J5POLxXXOt8qcnf6/kiI/yIIvXTi1z5Qax8YxF+FZHnIf4xT84BIfJ2aeXbiPCzh/TyDSRPjgmh8pcnUEm0+OWJVRJvlIfWPpAnfk3AOH5i5lWKx8+r2R9p96RBCVh5WD4WCA31VIGAYahPmDr3yXFMwHjyiYA3RB7oA1gJWGsist5wpafACRh3sSduS8RRn0rE2CXOCREn8nJLFjqzeh0WOnsCNyKGL3uCMkKm1j2QOZLo63moEzFc2YO4EzF82WElEdPXeYqI4cuewnQi7nxBuTkKnek5IWJuIU8dOiHT2bEhZF7WH8idkOHMnqcGMfPm/2RnEDOc6YOHlxU666eZiBnOnrAGIdPXcUzI9PVkkBeE9HWeImT6Oo4JGc4OvQYh86vH44ffYOnsEGOmMHg9IWQ4O1SZhMwl85wQ8yQ7y8VC6KycEDOclRNihrOTnknMcFZOiBnOTrmcmOGsnBAznB0SeupZ1BNi5lL1ZIMLN52VE2J2KUx1YoazckLMcHa468Tslc1OzJzQDxF496OzckLMcFZOiBnODuUjZXjVk1TiVpogUoylnhAzt+snz0HMcFZOiBnOygkxw9lpnUXMcFZOiBnOTjMtYoazckLMcNbrviZ01uteJnTW6/4ldHZacBEznJUTYuZGeCYJMcPZOeGKJPRWj4ga7k7rcucR+qtHxA2Hp3m5w0h6rGc5Uho7v+4lkj7r2Z6C9nK2JyEiOMOm5XCB13K25yG81rOciPB6hEP2TGzxcpZjEV57nLOMQ9rLWcYh8nKWccDr445rRF0a+Ptf/zcLP7/9D0u7EPSHzQhdAAAFw3pUWHRTTUlMRVMgcmRraXQgMjAyMy4wOS42AAB4nNWWO5J1NQyEt0IIVXdOWS9brhtOQsQGKHbC4vlaF4gI/pRo7li2Hq1W6/z+6x/f9v3z77/+8cu3/6afn38+J//5c3780K1v++0ff9/2sf+PvfsPPvyBSz/9+fN53HPf13qWeeR53ac7dr/3w7H36wtLhlsWJj8n97seu11bpr3uTcdSFTve+dzjJ15fhjtzf9l64tyLIS3veMt9btmrn/Rd42t3E9+24fx1nszVofjhOzGkESUx3O593lyIJqY9JDXndtri3c85tUPnvaJ1/6xjnN/epgqDuCuJ7FGEvnLERaMKS+NBHBDgQYf1fTkp1d4q2+5qOVpnn1c8a7cHfu4Oc86X3xXcX0H53D/uy9/UvnwZ/tu3DRa2NmFvZhf54PTE1XniUfdxq3xuRCV+dkbm+DGT/2gg0v3uMhDqc4/y3JVb+a91j5Bbnas4Lxp4XiB5QYjznYv3rrxdALUQBVBwW5u66B1t3cQJ4VnnUBbH4Cg3wFjF9WtWzXmtCnvVk31Ik/Dt1vy/uuAAYb3uJuytgAr1RNZywcwD9XHTNo5PbNpJ1lF+ue7nxrg5BUxFG/c1Xb+lLtaxTVdKqKZ4uss61F7vjBcULKrgvu0A3S8Z8rjq67gpQhT9yOFi07rmzVoRTmug+RiaDPPyJLp6KJH3KghIHHAgPOw1dZMWni3+Mj+lvDmPIi9IBDw4owsiqj+HkZgKqwHky588VUcNOceDWWhbvWWgwanOemxeqBaDSRjoeIsKFy5fXtR2r3lBaUklTk2KvpIuyVCpaQHOlUxDPoxfTPCqBFBe7KAdtBHgmFCNqrdjiAY2DGtrHDDsVB2UvtbS6C7bd8tA+hA8KL3EP3hKo3W+jwneuHFVxqq6g+/hJunawySauHPD5MoeC2iCIWmSsKpu9YqpLvQInCFDZfAEsrUIAfvjGj+aioYREKkkEMtN57svXRaDhinEOH4Pf8tVB/NrtSRBB+Xx4UWc9s/Ek580iIohFya/hNUUMDIihMdZ2ydzh5DydpUoqgLtRyTBIJmznHk94/8jafTgBMpCD/xOcegJeqnZckcE6XPl1FCxrpu0J5dx3p4+5AJSUhRiBwo0Jqby3o/JmMT8KOWWCGOaItA1eGZqOVSlPzacEPrMiJBYBShvU/yUJMEFjXzXokYc2Dk2CIplZAVeh/NgOsfBZha4MYTd9p4tsnKCohYuFYVNS2+2bwg28rZURwvAieHryjnYwC6fPURKFN7y7NSnLga9mQjQN0q6RxevNsd1Nfr9pfRtXPheA0iuYK70Bo2p0VDfLqBJ30lr3JmANK2B0oRv9sS6IExmoSE8TFZJPVgb9Q4tmFCUvaX7KZ1hXdKsj17ScpMMwd/9Vr86VV5v6RzD4qg2dV/krTVNF/2YkabCy4NiEnTO35YA1j3jnklK/c+Upu6fo7VBOmyRGD1mbjW2CJbSJ9Dt0e2QvHerYFrNNQkbWa0jd0gNx6i0ieigsetzDhxyE7PFAi5srYeFkIia2M8MeSwWj/aSQVaxBr3PPR8BsTdjwyMWFqtp6SsEPGZaPgJ7KC8ERx5gF6WYvGjVccAv1KGgzydmPDU+JZ0xVc7eLf98hfTWpqaWmQ3mrEaSAd31RVPeN6gy7tqjP3FRCwwwxP4Gy+aDJmFTzRInj4Gr6/NVYbOTTYxSr89oUSajpXOocoe1n+u8YvcSmJ7le9g7nz/B5CFCZlpK9pZInRErd1GXc7q/P3j4GpVCZPPzAsLzogF1BMcYLBhrTMIdYuor7Y7ms0YKbNHrDx4UM1visHD0ScCa87dmc38WEclOc1k3sWXQtL6+JFAcgKy+saRtpZGSwZmb/PfjCLavXp8XgdwIWX1GYYClNtvjOqCrSahOCkLqAHII637mQ46oLpXkM6vmRfQki5wfzhWi9SDpvzwx01G//PkXWyqNytIk2gYAAAAASUVORK5CYII=\n", + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcIAAACWCAIAAADCEh9HAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3daVwTRxsA8GcTkhDCfd+nUbwtqKhY64EoClargCfiASiKoqjUVkG8sWI98MCr3hbwKuCLStUiClZAVOQSEBW55SaEEJJ9P0zfNG84VE6R+X/ob9jM7s6m+DC7O/MMQZIkYBiGYa1F6eoGYBiGdW84jGIYhrUJDqMYhmFtgsMohmFYm+AwimEY1iY4jGIYhrUJDqNY66Wnp798+RIASkpKoqOju7o5GNY1cBjFWi8tLU0URv/666+ubk4rpaamokJKSkrXtgTrpnAYxdrk9u3bhw4dunjxYlc3pPX8/PwkChj2WaS6ugFY92ZmZmZnZ5eVlfXo0aOubksr1dTUvHjxAgC4XG5XtwXrlnBvFGsTNTU1Q0NDbW1tie2iqBobG9vpjfo8VVVVSUlJSUlJHA6nq9uCdUs4jGKtZ25ubmFhAQA6Ojo2NjZeXl6lpaXoowMHDkgUvlja2toLFy5cuHChqqpqV7cF65bwTT3Wevr6+qigoKBw4sSJkydPPn/+/NatW1JSUhwO5/bt2wBQVFTUpW38R1xcXFVV1aRJk/Ly8m7evOnq6goAW7ZsMTExGTduHKojKmDYZ8G9Uax9bN26VUdH5+7du+vWrevqtjShrKwMBXQul/vmzRsA2L9/v5+f39KlS21sbFCd5cuXd2ELse4L90ax9qGlpRUaGjpu3LgDBw4MHDiQxWJNmjQJAE6fPt3VTfvH5cuXnz9/XlFRoaGhcfHiRS8vL4Igjh07ZmBg0NVNw7o33BvF2s3IkSMDAwMBYPny5Vwu98aNG8+fP7e2thZV8PHxkSi0u4aGBtQGAPj111/FP5ozZ05AQMDGjRuzsrKcnZ2FQuHevXsXLVrUQS3Beg7cG8Xa09KlS0+ePJmRkREREREREYE2btiwwdDQcPTo0YWFhWhLRkZGBzVAIBA8efIElR8/fhwfH19ZWYnK5eXlTCYzLi7uxo0bDQ0NPj4+a9eu7aBmYD0KgbPfY+2ouLiYzWZXVVXNmDGDz+fn5OTk5OTU1tYCwKRJk+rr67/55hsAuHnzZnp6ekc0gMfjmZmZTZ48GQBevHjB5XIbD2hlMpnTp0+/dOlSRzQA64FwbxRrT1u2bKmqqpo+ffq1a9dEG4uLi3NycqSkpPbs2RMQEAAA79+/77g2mJubo7M4Ojrq6uoymUwAkJWVpdFoAJCYmPj69WspKfybj7Ub/MuEtZtXr16dPHmSSqXu2LFDfLu6urq6unqXNAnFU3Hv3r1js9mXLl3atGlT7969u6RV2NeGxLB2Mn36dABYtmxZcxXKysokCu1OKBTm5OSgcnZ2dpN1lixZAgCLFi3qoDZgPQ1+Noq1j8ePH48aNYrFYmVmZmpqanZ1c1ry9u3b3r17CwSCtLQ0Npvd1c3Buj084AlrH+vWrSNJcv369V94DAUAAwODefPmCQQCf3//rm4L9jXAvVGsHYSGhjo4OKirq2dlZcnJyXV1cz4uOzvb1NSUIIiMjAwjI6Oubg7WveHeKNZWfD7/559/BoBt27Z1ixgKACYmJvPmzePz+bt27erqtmDdHu6NYm116NChVatWmZqaJicnd6OBRLhDirUX3BvF2qS6uhoNb9qzZ083iqEAYGJiMmfOHD6fv3v37q5uC9a94TCKNau4uFii0NiuXbuKiorGjBljZ2fXWe1qNz4+PlJSUr/99hvK+YRhrYPDKNaslStXooKHh0eTFfLy8g4cOEAQRDft0PXq1Wv27Nl8Ph+/ssfaAodRrFl8Pr+ysrKyslIoFD579uzZs2c1NTXiFTZv3lxbW+vo6Dhy5MiuamQboQ7pqVOn3r5929Vtwbor/IoJa9aAAQMcHBwA4MGDB3Q6PTIyEgA0NTXZbDabzZaVlQ0MDJSSkkpNTTUxMenqxrbevHnzLl265O7ufvjw4caf7t27FyWivnz58pgxY3R0dDq9gdiXDodRrFkODg4hISEA4OjoqKamFh0dnZWVVVdXJ6qgpKTUv3//mJiYrmtjO0hLSxswYACNRsvKytLV1UUbSZKsqKioq6tbs2bN5cuXCYLw9fWdM2eOqalp17YW+wLhm3qsWQMGDBAVAgMDk5OTa2trDx06dO/evaCgoNGjR1dUVMTFxd2/f79r29lGffv2/eGHH1CGPTabraWlJScnR6FQlJWVR40aBQDr1q3z8vKKiorq6pZiXyjcG8U+j3gXtU+fPtu2bVNXV09KSmq8xnJ3QZLk5MmT792719DQIL5dUVFRT0+vX79+uDeKtQz3RrHP8+bNm7179+7du7ewsHDLli3W1tbFxcVz586ViEHdyNatW+/cuUMQBABoamreunUL5QcoLy9/8eJFV7cO6wZwGMX+1dDQkJqaisrPnz9vso6Wlta8efPmzZunoqJCoVDOnz+vo6MTHR29efPmTmxpu7lx48bWrVspFIq3t7eCggKfzz958iSdThdV2LNnD4qw7u7ueLIT1iQcRrF/VVdXHzx4EJV37tzZZB0Gg6GlpaWlpYWSyaurq4eGhtJoNH9//xs3bnReW9tDWlrawoULhULhnj17RowY4efnl5ubKzHLXl9fHxU0NDQYDEZXNBP70uEwiv2fkpKS+Pj4+Pj45m7SV69eLVEYOXLk9u3bSZJcvHhxTk5OJzW0zcrKyuzs7KqqqubPn+/l5QUASkpKaMURDPssOIxi/6eqqio7Ozs7O1sgEDRZwdLSEhXQW2xk/fr1M2bMKC8vd3R05PF4aCNKMg8AHh4eXC63I1v92RoaGuzt7bOzs83MzIKCggDA0tLSxsYGAHR1dd3c3Jrc68iRI6jQ5AhTrMfCYRT7PyYmJrNnz549e/Zn3cASBHH69GljY+P4+Hg0WB0ARFOeampqumpAyO3btyUKiKen57179zQ1Nf/44w8ZGRkAUFRUVFNTAwBpaWnRjbyE6OhoiQKGAQ6jmDgGgzFmzBhUnjhx4mftq6ioeO3aNSaTGRgYeOHCBQCoqKjw8vLy8vKKj49v/7Z+mtOnT0sUAODMmTOHDx+Wlpa+fv26aLz9p3j//j26og5aHRrrprpTZjOso8nIyIwbN05ZWVlFRSUzM/Nzdx88ePDevXtXrFixfPlyc3NzRUVFtDBnWVlZBzT2k5SWll68eBEACgsL0ZbY2Nhly5YBQGBg4IgRIz7raLq6uuiK0BxZDENwbxT7PwRBlJeXS6Qg+XTu7u4LFiyoqamxt7cXDRtiMplozFDnYzKZffv27du3L4vFAoB3796hCUteXl6iR7cY1kY4jPYIt27dQjlDk5KSXr582UJNlHqZz+e3+lyHDx9ms9kpKSkREREzZ8786aefRowY8fLly6qqqrdv327atAlVW7BgQatP0aSampry8nIA4PP5or6njIyMmZmZmZmZnJwcl8udNWtWUVHRxIkTW5fZb9WqVVevXk1MTPTz82vPpmPdHL6p7xHi4uIMDQ3V1dXT09OZTKZosnxjKIy2ZUqSnJzc+PHjc3JyKioqrl27Jv7RunXrKJR//nLX19e3+hRNiouLe/PmjYuLS2Fh4YEDB/bu3VtbWztlypSkpKSSkhJ9ff3hw4e/fPmSzWYHBwe3LlE/h8OZNWuWlZUVnl+PicNhtKc4duyYqqpqcnLynDlzWqiGBtW3JYxWV1eHhIQ0NDT89ttvTCYzLS0tPT09IyMjIyNDW1v78uXLaJBmcwOq2qKysjI/P7+oqIjP58vIyEiMsmIwGEwm8/r160pKSq07fklJCQCoq6u3Q1uxrwgOoz3F3Llz2Wz21atXW67W9t7osWPHysvLx4wZ4+zsLL5dKBS+fv26sLAQpZp3dHRs9Sma8+TJE4IgKisraTQaSZJMJlNFRUVVVVVVVVVZWfnmzZscDqe2trbVx0cPRtC4qOZkZmYmJCTMmTOnrq7uwIED3t7erT4d1l3gZ6NfD1HezyYTgMrLyyspKaE3LS1A74IaGhpad9PN4/H2798PABs3bpT4iEKhMBgM0QrMioqKrTh+yyZOnOjl5eXi4gIAZWVltbW1ubm5SUlJUVFRwcHBaFD9qVOnWn181BttOYxWV1e/e/cOAAQCQUZGRqvPhXUjOIx+PQIDAyUKInPmzCkrK4uLixsyZIiFhYXEpwKBIDEx8cCBAw4ODpqamhQKhSAIFos1dOjQ1atXnzt3LiUlBcRm2QcEBIimKkk4e/Zsfn7+4MGDJ02a1PhTPT090SsmNHeoHampqRkYGAAAk8ns379/42mdLi4uBEFcvHix1eMQPiWMAsDNmze9vLx++umn1p0F63bwTf3Xo6Sk5ObNmwBQWVmZkpJCEISpqSl6pWNqamptbR0VFXX79u2+ffsCQH19fXx8fHR0dHR0dGxsrERkEQqFFAolMTExMTERbRk8eHCfPn1Q+eXLl00+2RQIBHv37gWAjRs3dv4IpyFDhqCCqqrqokWLGlcwNTUdMWJEXFzc1atXFy5c2IpTfOKz0alTp3p7e3M4nOaWAsS+MjiMfj2oVCrqglEolK1bt4aEhMjJyQ0aNGj06NGWlpaVlZUAkJub6+/v//DhwwcPHlRVVYn2NTY2trS0HD16tJWVlYmJCZVKLS8vf/78+aNHjx4+fPj48WM9Pb3Kykr0dqi5WUlXrlzJzMw0NjaeOXNmp1zxZ1uyZElcXNypU6eaC6Pr16//5ZdfxAvv3r3LyMhIT09PS0v766+/ACAiIkJbW3vo0KGiUQdYT0diXwsHBwdRYfXq1Xp6ei38fycIon///u7u7sHBwQUFBaKDoLt1Op0ufmShUFheXi46vrOzM4fDadwAMzMzAAgKCuqwS2yrmpoa9HA2NTW1yQri3+H8+fNlZWUlvjdR6FRRUbG3tw8KCsrPzxc/Ao/He/XqVVBQ0O+///7hwweJ44u+nC/5W8I+F/5z+vUQ9QFnzpy5f//+d+/e5efnh4WF+fr6WllZUalUBoOhq6vr6up69uxZDw+PqKiow4cPc7lc8Xmf6G6dSqWKH5kgiI++EYqMjHz69KmmpqaTk1N7X1m7YbFYaB7nmTNnmqyQnp6OZs3n5uaSJFlTU6Ourv7dd9+5ubnt27dPUVFRKBQuWrTI2Ni4tLQ0NDTUzc1NT09vxIgRW7dujY+PFwqFdDq9uLjYzc3t0KFDKioqEse/e/euRAH7GnR1HMc6yaBBgwDg2bNn6EdPT8/379+TJHns2LG7d++KqqGHpCwWq/ERRN2uwsJCoVAo8SnKabJnz54OaX37iY2NBQANDY36+vrGn4r3RvPy8kpLS8U/RVmv6urqSJLMzs4OCgqytbWVlpYW/WtSVVW1t7f/8ccfAWDYsGHJycl//PHHr7/+6uHhMXXqVD8/v1GjRm3YsGHDhg1DhgzphIvFOgcOoz3FN998AwBPnz5FP3p6erq6uq5du9bKyuru3bsNDQ1oO3qEKi8v/1kHf/z4MQAoKChUVFS0c7s7QP/+/QHg+vXrjT+aMWMG+oocHBzu378v/hGHwwEAGRkZiV04HE54eLi7u7uxsXHLXRZra2sHBwc+n8/n8+3t7Tvs+rDOhm/qewr0UE8oFIq2+Pj4BAQEzJo1Kzc3V09P78CBAw0NDU3e1H/Ujh07AMDDw0NBQaFdW90hFi9eDM0MIN28eXN4eDgA/Pbbb8ePHxf/CE3Yl5gBtXHjRhkZGVtbW01NzZSUlPT09H379g0bNoxKpVIolD59+tjY2KxYsSIgIOD69etogRYpKSkpKamuytWCdQT8pr6nkAij0tLSaAudTr99+3ZBQYGnp+fRo0fRDelnvYNOS0u7efOmjIxMdxnfs3Dhwp9++ikyMhL9/RD/iE6np6WlXbt2TfzvDVJRUQGNZg1kZ2ejQk5OjkAg6NOnT58+fYyMjGbMmGFnZxcWFiZxENHCf910BUCsSbg32lNIhNFdu3ZpaWkBwKJFiy5evBgcHMxmszMyMtCISx6P95///OfFixdFRUUfPfLOnTuFQuGSJUu6y2RzFRUVOzs7gUBw7ty5xp/KyspqaGhoaGhIbG8cRhsaGtAgMInU1C0MLxUlhWkhOwzW7eDeaE/R+KZehCAIBweHH3744fTp05s3b0b5RqdOnSqqoKSkpKWlpa2t3fi/BEEEBwfTaDQ0pLS7WLJkyZUrV06dOrVx40aJrreenp6lpaVAIDh69Kj4dvGb+tjY2F27dikpKSkoKKBEzuhBAfKJk52wrwYOoz3F8OHDWSyWvLx8cxWkpKRcXV05HM7atWsNDAz69OmTl5dXWFhYWlpaXl5eXl4uWsJeHIVCEQqFtra2aCJmd2FtbW1gYFCek/MoIuLbyZPhfxmmNTQ0zM3NKyoqFBQUJFLco94oj8ebMGHCvXv3AEBRUXHs2LGND47DaE+Dw2iP8PPPP+/bt49CoZw4cUJGRsbIyKi5miiTpp+fn2iej7u7+7Jly/h8/okTJ+Tk5FgsVkFBQX5+fkFBAYqzFArl7t274eHhdnZ2nXQ9bUahUJbNnu3p7y/9/feV+vq0mJjM8vILFy7MmTPn1KlTtra2ixYtsrS0PHjwIIfDqaysrK6uRvNi0fejqKjo4eGxevVqUYb/gwcPitK+4DDa43T1UAGsM8yePRsNadqyZUtKSkpz1Wpra5lMJoVCKSwsFG1cuXLl4cOHz58/7+Tk9ODBA/H6x44dAwD0jJVKpe7evbvjLuHT8Xi88vJykiSFQmFRUVFz1YTS0iRAPQAJYP9p/1gUFRWlpaV9fX3R8ZtjbW0NAJGRke1/bdgXCfdGe4r169cTBBEXF2dv32zQiIqK4nK5FhYWEi9YTE1NlZSUGq8+goYEBQQEFBcXr1279scff8zOzj5y5Ejrcsu3l9TU1MjIyI0bN9bV1a1fv/7s2bON62RnZ+vW1zMAagEUAKhUqjSNJi0tbWho+Pr166qqqrFjxxobG7NYLBkZGZRgsLCwcMeOHaamplu2bGm5AcOGDQOAlifjYl+Vro7jWGf4xN6oq6srAGzdulV848qVK1Gfbt68eag3ev369bKyMjTkXlVVFc3qCQ0NRYlRrK2tKysrSZIUTQEqKytrPOup4yQlJfn4+JSXlxcUFDg5OTWuEBkZqaSkFAxAApAAhRSKMDv7zZs369evJ0kSJbJqPOm+vr5eRkaGIIji4uIWzr569WpU8PT0bKcLwr50eMAT9q9bt24BgPg7egCYOXOmnJzchw8f3r59a2JiAgCPHj2qrq5GCUOXLFmCpkjOmjXr3r176urqd+7cGT169Lt375YvX46OsH79+urq6s68kLi4uEOHDqFX7ZcuXXr27BnaTpKkv7+/ra1teXn5palT39+4MUFBYbBQGFtQQBAEGhKPxng1HvBEo9GGDx9O/m86aXPy8/MlCtjXr6vjONYZqqqqUIHD4YjmfUooLS0dM2aMnp5e457j06dPAUBeXt7Hx8fHx2fUqFHJycmoa/bq1SvxmtnZ2SifqZaWlpWVFYfD4XA4zs7OqH/aOZKSknbu3EmSZG1t7Zw5c9DgBHNz86CgoB9++AEACILw9vZOTk5WVlYeN24cAIg6rXV1dQRB0Gi0JrvPaMz8unXrWjj7gAED1q5du3btWktLy464OuwLhMMoRpIkuWvXLlTYvn1740/RRHsqlYreWa9evRoFlEmTJjWuXFZWhoYBUalUBweHTZs2DRkypDPDaFZW1qVLl0iS5PF4Pj4+np6e4lNUGQzGb7/9VlJSgnrWU6ZMoVAoTCYzMTHxwYMHaAUUVVXVJv/Y3LlzBwAsLCxaOLtosrwoywn21cNhFCPJ/89s1GQFNCcnOzs7LCysX79+TCYTZdffv38/yhQlEhUVNWTIEE1NTQCwtrb+8OGDhYXF9u3bW3673aGuXLkivqYInU5Ho+iVlJQGDRokGrQk+hT1pr29vSX62tXV1WhGfHV1dXPnwmG0ByJIkuywBwbYl6WyslJGRoZGo6HVMWVkZEQf9e3bd8qUKQDw9OnT+/fvS+woEAgGDRqUmpoqKyvbeCEjCoUyatSohQt/19cvCQjwRl02PT293r17h4aGzp8/n81mL1q0qK6urqSkpKGhYfr06R17nf+TlZXFYDAoFIqJiQmPx3N0dFy6dOnx48evXbumoKBQVlYmXllKSsrIyCgrKwsAtLW18/LyAIAgiG+//Xbx4sWzZs1CyQNdXFyUlZXXr1/f3LDQhoYG0eqqXTtiAes8XR3Hsc6zdevW5ORkkiSDg4NDQ0PFP2quN/ry5Utvb280MhTp16+fr6/vlStX0JsleXl5FovFYMjTaK4AVACg02VXrfI9eZKXkEDm5eVNnLhMdEPv6+v766+/RkVFiefb7yAuLi7+/v6rVvnGxDxetmyZo+NqoVBYXU1yueT69T8BgJHRWGvrI6Ix86LUTcuWLRMKhTExMc7OzqLs9/Ly8u7u7g4ODutmzCDt7WN1dX+ZMKGjLwHrLvBfy56ltra2urq6rq5OSkoqKirKxMREX1+/cacpLS0tODj48uXLr169QlsUFRUrKipcXV2PHTv24cOHN2/eoJR6VVVVBEEQBEUoPA7AAFhVX/9zYKDK7Nlw6xYcOCBVUMCvrQV5eRAIBP7+/nV1deiAKioqY8e+0dCQHTgQ+vWDAQPg8GEoKIAjR2D/frC0hGHDWnmNr1+/5vF4NBptw4YNM2dCv+/ZRxWFscIR9ZyGo0dpMjIQEcEBgJyc6fX1LhyOu6GhIequPnr0iCTJY8eOpaSkeHh4nDhx4siRIxEREcePH797925xcbFWba1TYiIUFY0EGJqfD+Xl8P9587CeCYfRnuXMmTPq6uopKSlDhw6dN28eANBoNAMDA11d3ZUrVxobG6uqqg4ZMuT58+eovra2tqOj4+zZsy9cuHDo0KHTp0+fOHGC/P8HQSRJkqQAwF5efpe5uQkA8PnAYsHs2RAcrM7jlT948KeRkQKTydy3b19ycnJKSsrLly+rqrjh4bL19f8eZ/VqqK2Fx4+bbXxeHvD5YGgIRUVAoYDorjo/H549S4iL+yM+Pj4hIaG0tPTIkSPoYe6RKLZydRaUwSh4Td6bPWaMnZwcnDr1BgAGDdLT16+sqGC9efOmpqbG2tr69OnTp06dOnHiRExMTExMjKamppub24oVK+zt7TMyMgQCwS+OjgOLijhUKk0goAuF0Mwq01iP06V9YaxTid/U79q1a+zYsbq6uhL5g1H0UVRUXLBgQVhYGJ/Pj4mJsbS0BAB0Fw8AKioqffr0GTVq1LRp02xsbABATc2of3/yf+PZyV27SBcXsq6OXLSItLLKO3DgzPHjxyVu5PPyiu/cIfftI5csIYcPJ1VVyS1byBcvyB9+IPfuJZ88aaL9J06QFhYkj0f+/jt57Bi5bRv5/fektjYJQI4d++/SRhoaGmfPnp06dWp9fb2QoJAAVQxVEkB44w90nH79zAAA4AmDQVIoNBaL9d13302bNo3P55MkyeVyz549O3DgQHQ0Op1ub28fGxtLkuRqGxsBQaCL5FGppNhLs7i4uOjoaJIkc3NzL1682FH/C7EvEu6N9iAGBgboUaC6urqGhgbK0FxXV/f69eusrCz0Fv6vv/4yNDR8+fJlRESElJTU3LlzQ0NDAUBNTW3GjBlXrlwpKyvT1dUNDw9Hkx1ra2vl5OQqK/Pev6/PyqKHhkJoKPzwA+zdCwCwcmWxuTkrKemXkpJkiXitra2mrQ0TJ/67Zds2kJaGJUvA2xvGjGn6EuzsICAAjI2hoAD8/P7ZqKgIRka9R4z4cdiwYUOHDtXX10cX6+npuZKu1JdXKsf7UEuwmCYmQiH4+ha8f58NAHZ2OoqKF7OyLCZPnkihUNBlAoC0tLSTk5OTk9PDhw8PHjx4/fr10NDQ0NBQc3NzZWXlo+PGrejV69GjR48NDZfT6bzy8srKyoaGhrdv36IkhNXV1aInIVhP0dVxHPtSZGRkcLlcQ0NDAJg4cSIAoMAnJyfn7e2NBvCnp6ez2WwA0NXVffHiBUmSBw8eRLsEBgai7lhcHDlvHhkTQwoE5NWrVwFGWFnZcrlcicFDjR06RObkkCRJOjmRcXHk6dOSFU6cIKOjSRcX0t+fvHCB9PQkL1wgMzJIiZHyDx8+PHXqFBq+Osx4JxkUJDwW5O+QUFpKLltGAjwBIBgMRl1d3fz580Wv1ObOnYt6oxKys7O9vLxE2ZoHDRqkpKTUeOHl3bt3T506de3atc7Ozr6+vp//9WPdGA6j2D+cnJxqa2vRFE/R++vFixdLrI5ZWlo6btw4LS0td3d3ctu26JEj1+joAMDy5cvPnj3r7+/fu/dEZeXLc+eGhISEoNR5Pj4+gYGBFArF3t6+uQXixRUXkyNGkADk/v3/bCkrIzdvJoOCyOhosqCANDYm//hDci+BQBAWFjZq1CgAkJGRef/+PUEQFApl+fLl9vb2o0aNp9MHU6l6srI/AwCTyfL19Z0/f/62bdvQ7jt37hQIBE22JyUlpbKyEqVcEc9+Ly0traSkpK+vb2Ji8ssvv6Ax/6mpqTiM9jQ4jGL/cHJySkxMfPbsGZrzo6OjAwBNjiHn8XiZmZm/DxhAslgkgBDAGsDExAQtbCdBSUnJ0NBw2rRpaFg7hUJxcHB4/vx5c83IzibZbHLCBJJCIQmCPH+eTE4me/UiAcgFC0jUow0NJf+3UDRJkmRlZWVAQIAob7SKisqUKVPYbDaVShU9zxVRVZ2NCseOHZs/f/5Hv5aamhodHZ3evXs7OzsDwLfffltaWiqaXCty/fp1HEZ7LDz8HvvHwoULR40aRaPRjh49mpCQYGBgUFRURKVSc3JymhxqXslkKtTVfaDRVPn8pQDXlZUHDx7M5XJ1dXVR9ngKhcLn8588eYKSdBgbG5uYmDx8+JDL5RIEMXXqVB8fn2GNhjU9egQTJwKFAhMm/CcsbEq/fsU5OSpcLnXAAKirg59+gimQpSoAABOCSURBVEWL/q1cUFAQFBR08OBBtMKHgYGBiYlJQkJCVVUVALDZ7O+//97Y2FhZWVlVVRVAxd1dobp6XkFBHABs27YtIyPj/PnzLX8tmzZt2rFjR//+/V+9eiUUCpOSkkRvn8SheU3y8vL19fWVlZU4Z3PP0tVxHPtSoJt6kiS3b9+OHne6u7vn5+eLKgiFwrKysry8vKysrISEhGJpadG7+QUAKioq4eHhPB4vNzc3Li6Oy+Wiverr68+ePdu7d2/0+6anpzdhwgTRBCorK6u4uDiJlly/Lhg4cClBEAMGWBPETAChuflFK6s7AKSl5T9PQhMTExcsWCAa8TpkyJAxY8bQaDT0o6WlZUhISON58Q8exKC1oy0sLPLz8xMSElr+TvLy8qSlpQmCsLCwAICVK1e27TvGvk44jGL/OH78eH19PUmS4eHhv/76KwDIycn169fPyMhISUlJWlpa4g/wtwDxAAkAwYqKNAqFIIiysrLmDi4QCEJCQkxNTdG+urq6kydPFi0MZWlp+eeff4rX/+WXX9ALLoIgAMwBgCAoXl7ZNTWCsLAwKysrtCONRrO0tBR1aalUqq2tbeO4LA5lJEEBl06na2lpmZub29vbe3t779+/PyQkJCEhIS8vT1T/jz/+mDlzJgAoKyt/+PChzV8z9hXCYRSTlJKS8s0338D/VgcRp6ioqKWlZWxsbGZmNmLECCsrq1mzZm3fvh29eEHxq76+HiVybgy9BUIHBwA1NbVx48bJycmhH9Fz2IaGhhs3bqBBSywWa9q0aTIyMlQqlUaj2dranjx5Eo1slZOTGz9+PIqJACAvL79q1aq3b9+2cF2lpaUoIAKAurp6y/fdsrKy/fr1mzx58uLFi1Hu0cOHD3fEt419BXAYxf4lEAj27duHOp69evW6fPnyy5cvs7OzS0tL0f1+c1Cm5wkTJhw8eNDZ2XnkyJEtpIgXCoVhYWFDhw5FAUtVVdXKykpJSennn3/evXu3rq4u/H/aFNTNFA07ZTKZ+vr6op6skZHR7t27P5o+6t69e+jIaDllNA2hrq4uOzs7JiYmJCRk9+7dq1atsre3Nzc3F38dj74KCoWyc+fO5lK1Yj0cDqM93f3799GQphcvXohWA3V1dW0hF1xj/v7+AODh4ZGfn4+iVe/evTMzM1vYRSgU3rhxQxRM6XS6xPh81BuVyGInjsFgTJ06NT4+vuW28fl8X19f0SNRNJaroqKi5b0qKiqSk5MjIiKOHTs2a9YsdMbhw4e/fPny078WrIfAYbSn27RpU0ZGBkmS58+fP378uLGx8Y0bNz73IP/5z38AYNy4cSRJ5ufnm5ubA4CysrLESqISMjMzfX190U26eHC0tbUNCgrKzc0lSVIoFGZnZ1+7dm3Lli3ff/+9lJSUtLQ0m81GwRpBSaeys7Mbn+L169cjR44EACkpKW9v78LCQgCQk5P73AuMjIxEzxkYDMa2bdvevHmDtqOvDuvhcBjt6TZt2oSGxy9cuBBNom/FQXJzc9HtOfqxurra1tYWBZ3GE8x5PF5ISIiVlZWo+2lqajp58mQU4MTXdm58FisrKwsLC7TCR0JCwqpVq0RRmEKhWFpa7t+/v6SkBNU/e/Ysmm6kqam5bNkyW1tbWVlZBoNBo9FcXV0/t1/J4XC8vb0pFMqKFSt+/PFHtBHnZsZIHEaxTZs23b17Nzs7e9++fWFhYa0+jrKyMgCIgmBDQ8OKFSsAgCAI0XD09PR0b29v0bsdaWlpe3v7qKgoFBYnTJgAAHv37pU4Mp/P//nnn1F5w4YNEvfjPB4vLCzM0dFRlN+ewWDY2dmhJaHgf1NaGyMIYtKkSWFhYc1NXmpSTExMWlqaq6trWlpaWloaDqMYiVOTYACgq6trbGzclhHjQqFQQ0OjrKxs7ty5ffv2lZWVZbFYenp6M2fOvHbtmp+fX1RUFEEQjx49QvXNzMxcXFzmzp0relMEAAEBAZGRkWiR55ycHAaDoa2tjQ7+9u1bVOfNmzfiCysBAJ1Ot7Ozs7Oz43K5ERER586du3XrVnh4uGj+EpVK7d27N5vNvn37dl1dnYeHh7u7+9GjR0+dOnX79u3bt2/r6OgsXbrU1dU1MzPzu+++A4CoqKiJ4klTxIwePfr169e5ubkxMTEAQOLZKxjg4fc9XmJiIlpvLisr6927d5+7u/gYJonX6+LQzbW8vLyrq2sLg95/+ukn9GLq7NmzERERaCOPx+vbty9abtPKyuqjTcrPz+/Xrx8AzJ49Oy4uDo2+cnJyAgAXFxfRkp9lZWUBAQG9evVCLbSxsXFxcUEftdzHzM7Oxjf1mDjcG+3pzMxQ8k0QjcH8RA0NDZcvX96xY0dGRgYA6OnpCQSC2tpaDQ0NFxcXAODxeJWVlXV1dbGxsVlZWVOnTg0ODhYlPWlOTExMVlZWcnIyWl4UMTc3DwgIAABHR8ePNkxLS0tDQyM1NVVGRobFYjEYjKtXr06ZMmXEiBFubm6ie3wlJaW1a9d6enrevn07MDDQ2dk5KCgoJCQEAEQp+pskJyc3ZMgQVP72228/2h7sq4fDKPZJampqMjMzUa/z/v37ubm5O3bsQIk1DQ0NPT09ly1bxuVyraysEhMTQ0ND//rrL7Q4KAAcPHhw9erVRkZGH42hAEChUKhUKhrd2WpozT4pKSmSJAGAJEmU6qnJ09nY2NjY2FRUVFy6dAnlN0Gjo5qjpqYmiuYrV65sSzuxrwMOo9gnKS4uDgkJQWH0119//fPPP7lcrpGR0Y8//rh48WI0t53BYERGRo4fP57H4z1//hyF0WfPnqEHoNXV1Z9yIktLy169ehUUFIi20Ol0lL4PAH777bdPOQgKozQaLSAgQFVVNTMzc8GCBR/dS01NDc2dF83Nx7BPgcMo9qlev36NxocyGIwtW7bo6OjMnj1bouOmpqb2559/CgSCNWvWTJo0CQB27drl4OAAACjrUsvGjh2L3vgPHDhQ/EmrqNzC41dxKIzS6XQvL69BgwZduXKF/Ni7ICaTKZoqilapwrBPhMMo9qnE03du2LChuWpoBjqXy/37778BoLy8/NN7o6L346J5962DwqiysjKaBCUjI/PRMMpgMFDcB4Bp06a15exYT4PDKPapDAwM0NBO0S12C/h8PhqlhBZrgk/rjbYXFEZXrlyJZsdPmTKl006N9UA4jGKfRElJafz48agsmmPeAnl5eXQvf/Xq1c96NtouUBj9xCcAGNZGbXofivUcSkpKojtuFB9bZm1tLSqgMNppvVE+n8/n86WkpFpIa4Jh7QgvIoJ1uMrKSkVFRTk5uc6JpOh08vLylZWVnXA6DMO9UazDycnJEQTB4XDQSu4dDd/RY50Mh1Gsw1EoFBaLJRQKORxOJ5wOh1Gsk+EwinWGznxZj8Mo1slwGMU6Q2e+rMdhFOtkOIxinQH3RrGvGA6jWGeQ6I0WFhampKQAAMr/1OQutbW1PB4PAOrr6z/loSrKwA8AjMuXowD2pKbCmzft0XYM+wgcRrHOgHqjmZmZ+fn55eXlqamp9+/fB4CqqqrLly83uculS5cePHgAAE+ePDl9+vRHT7Fu3ToAgDNnDM+dswIYVlwsdHZuvyvAsGbhMIp1IC6XGx4e7ubmFhUVJSMjs3z5ch0dHWVl5djY2OfPn1+9evXmzZst7M7n83k8Hp/P/5RzCYXC3NzcP48f1+bxCgkCAJ4+fOjk5BQeHi4QCNrnejCsKTiMYu3g9u3bDQ0NABAfH19cXPzhw4czZ85Mnz5dRUVl2rRpx48fr62tra2tlZKS0tTUVFJSkpGRkZeX19LSklgWVMKFCxf8/PxQcjxfX9+EhIQWKr9//97MzOyPuDgA0CRJACgQCM6fPz9t2jQjIyNvb+8XL14AgOggLR8Nwz5Dl+bex74SCxcurK2tJUly165da9asEWXPo1AoI0eO9Pf3T0tLQ2n2d+/eTZLkgwcPDh06RJJkUVHRypUrmzzmiRMn7ty5Q5JkTEzMmjVr0AH79u27e/fu4uJiUbW3b9/u27cP5QkFAAJgi7R01Zw5pIfHq7//3rRpk5GRkei33crKyt7eHu2I1//A2gtOTYK1j/Pnz9NotMTExLFjxxIEYWlpaW9vP2vWLB0dHVThl19+mTBhwo4dO5ydnQ0NDdGbdFlZWTs7uyYPqKGhgfIzycvLs9nsNWvWXLhwIS0t7ccff/T19Z02bdo2Obne167JVFeHCQR/A1Cp1EmTJo0ePZrBYEgtXw5MJhtg2/Dh27ZtS0xMPHfu3O+//96/f//09PSHDx9C52acwr5yXR3Hsa/BwoULExISkpOTPT09nzx5gtbIa8zGxgYAVqxY0bqz8Hi8K1euTJ06lUqlbgCoBSABSIAnAL169fLz8zty5EgLu9fX15eXl48fP/7atWvXrl2bMGFC65qBYRLws1GsffTr12/AgAEaGhpSUlLiyyaLCwgIkJKSCgoKSk1NbcUp6HT6zJkzIyIicnNz7aysmABFFAoAKCoo+Pv7N7cevQiNRlNUVFRRUZkxY8aMGTNUVFRa0QYMa4y6ZcuWrm4D1u3Jysr26tWLQqEwmUx9ff3mhr6rqanl5eXFx8e/f/9+zpw5LRywvr4ePWAVFcTJycnpy8uTISGyJAkAvAkThm/fTpLkhw8fhg0b1nJTORwOWteztrZWtMAnhrUFTpSHdari4mI2m11VVfXnn3+iXPpNcnR0DA4OFi80dsXNLer48W8sLKzOnOllapqVlVVYWDh69OiOajqGNQPf1GOdSl1dff369QCwbt26NubNS9XROQ6Qb23dy9QUAHr16oVjKNYlcBjFOpuXl5e+vv6zZ88uXrwovp3D4dy4cWPp0qWhoaHJycleXl5eXl5FRUXNHQfNEGWxWB3eYgxrER7whHW4ixcvmpqampubx8bG5ufnz5o1a+vWrc7Ozhs3bpw5c2ZNTc2tW7ciIiIiIyNramoAoLS0dODAgQEBAQDg6OjY3GFRZVlZ2U67EAxrEg6jWIcrLy9HWZc4HE5ubq6/v39FRYWGhkZeXh6bzS4oKEAP6NFYfTs7u++///7UqVNoX11d3eYOi3uj2BcCh1GsMxw+fPjGjRvv3r0bNGiQj48P2kin0/Pz82VkZMaPH29nZ2dra6utrY0+Ql1R8UJjOIxiXwgcRrHOsGLFim+//TYqKiolJcXb21teXr62tnbHjh1KSkoFBQUMBqMVx8RhFPtC4DCKdSpZWdndu3cDQE1NzY4dO3g8XutiKOBno9gXA4dRrMO5u7ujKUZWVlaiccosFotCoXC5XIFA0HiA/afAvVHsC4EHPGEdjkKhoDBKEASF8s+vHEEQLBaLJEnUqWwFHEaxLwQOo1iXQSnxW73OHQ6j2BcCh1Gsy7QxjOJno9gXAodRrMu0cdVl3BvFvhA4jGJdpi29UT6fz+fzpaSk6HR6e7cLwz4PDqNYl2lLGEV39HgxeuxLgMMo1mU+K4yWlZWJCs+fPx81ahSNRuNwOEOHDj1+/PinLGSPYR0Eh1Gsy3xWGPWfPRuGDoWhQ8+OHWtubp6ens7n8+l0emJiopubm56e3urVq1NSUjq4yRjWBBxGsS7z6WE0as8ej9hYSEyExMTVyclMgUBeXj4yMrKkpOTkyZPDhw8vLy8/ePDggAEDRo8efe7cOTRRCgD8/f079howDIdRrAt9ehh9m5Kiy+HUSUsLCIIC4Obi8vr168mTJ7NYrCVLlvz999+pqane3t7KysqPHj3at29fUlIS2vHp06cdew0YhsMo1oUah9Ft27ahwvbt28VrVsvKcglCuq6OSpIlDMbeQ4ckFqRD69e/e/fu5MmTfn5+79+/R1mf09LSOv46sJ4Oz6nHugyPxwOA2NjYoKAggiAMDAySk5PRR6ICUqqkdHTmzLVDhyYlJZ3lcPY3k80EdU4B4NKlSx/N+oxh7QWHUayzlZaWXr169dy5c7GxscrKygkJCQkJCQBgZ2dXVlbm5eUFAI1XYLbbuRPY7JQLFyYoKHRBozGseXhlUKyjFBUVvX79euTIkQ0NDbdu3Ro5cuTVq1d///336OhotJgdmsdZU1NjZmY2dOjQwYMH//XXXyEhIdBoQdDc3Fw1NTVpaekPHz5QqVQlJaWWT/3+/XuUNl9UwLCOg5+NYh0lPz8/OjoaAOrr68+fP6+lpeXm5nb//n06nT5jxozg4OCioqI7d+5QKJQXL164ubm5u7s3dyg9PT1paWkAUFVV/WgMBbGlR3AMxToB7o1iHSUpKWnr1q2TJk3i8/lPnz5NSUmh0+lOTk6Ojo4KYjfmHh4egYGB9vb2ISEhRUVFGhoaACAqYNiXDz8bxTpQr169JkyYUFdX9/Tp07i4uCbTM+/atUtdXd3CwgIANDQ0oqOjv/vuOxxDsW4E39RjHUhFRYXNZpuYmABAcynuZWVlN2/eLFoK9MiRI53XPgxrD7g3inUUfX19GxsbAGAwGM7Ozi1XLikpCQ8PB4DKyspOaBuGtSMcRrGOoqKiggbJU6nU7777ruXKVCoVjcYXrTKCYd0F/pXFvgjKyspjx44dO3YsCqYY1o3gMIp9Eezt7SUKGNZd4AFPGIZhbYJ7oxiGYW2CwyiGYVib4DCKYRjWJjiMYhiGtQkOoxiGYW2CwyiGYVib/Bfo/YWlwRtSUwAABvh6VFh0cmRraXRQS0wgcmRraXQgMjAyMy4wOS42AAB4nO2Ue1DU1xXHz+6yy7KwD5b3ex8sLK/lGfABv3suIiSgSBQrESULBFy0KIoKaSgCPsOgwlQMVisgmdU2VR4qFS37+0XbhhqIFrU6Eh9R0RANSIhBrWDXJDP9o+1//afT3pkznzvf+c6dc849945Zjt4E65LDP1aJNdZbo4pna9RbKRD8E1XBVvL4NqIfhH/BV4aXAoFIpbJu/vM+1Y/6v+WPtv8O8L5Pm8f7P/83aW8DJr7QxigSmvi2IqMtGMW2YCcGO56JL7Ez2ktMfAd7owPPKHUAmRRkfBNfLjMq5Ca+I9/oqDAqHcFJoOc7KcHZSeXirOe7uqhcrYKbDbi5grsUPNzB0wO8pODtBj5u4OsDfmJQ+YFaBRoxaDXgrwWdEgJ0EBgAem8I0kOwEkKCIdQXwkLA4AXhBogQQqQEouwhWg4xCnjNGWJdIM4dZrjDTA+Y5QGzPSHeExI8gfEC4gXoDdQbEn1gjg8k+cJcX0j2gxQ/eF0Fb6ggVQ1papinhvkaSNfAAi1kaOFNf1joD4v8IVMHi3XwkwBYEgBZgfBWICwNhGw9LNPD8iDICYK3g8D6J+QGQ14I5IfAO6FQEAqFobAiDExhUBQGKw2wygA/DYficFgdDmtE4Cy0NlsINkKRrZBnJ7F3EPJlcoWjUODk7OKqENvxHKQyvqPSSeBq7ZutW6rVzoMfAkrKsrdj4ONmRrwxF9+fE461M3uYTVc5vDnPjL9p0THEaxQvDl1Gnv4Ik/7hQfQJPICJkiuWJUf90Xyi2qrr2dYVLLnExuCGPfm9XdJholIkYMlXtmSsPhE3hD4n29o6yY2SjVgcfpVEd8WSaWrGz05NMJMRjcyI7hSuIhnkNnfGUhp0C+22O6BOd8DS/+wuMhkU8/Nq2Lt6Aa1N2YkZH4xZxgdFtGmp2erxZs1D7jTT5gJ28CJ7y0qcadi1gVe5MdFtL3Gr7ShmbN5P4o8P40DSd9idpSSb1p3GRfukdGF6c4Kh4gg6rhLTX488sYR6VeGWWcNYlWi2jJ7agLvCz+OFzmLWc28Y+nTI8FFdBVux0Z0G1E6ThoI+JleloW9vNpNjsUeZ7GNx1CcnGWeYVdwK10nMmNOH/Wt8uZFaBpnzX+Dpp3Os5xwlEbwHGLItkHw4vZTEa3yJ7ORZks6fi+KVKqJd6o7X2x4QD7kzeytdiBn5xaQwIhDVkVF4zLsaHR7dIzYry9AU8gEOPVLh13IzyrkmfGZyoM1324hM8RJLFkhp/C8luLTaiWYvVNEsywgRdqnpw6lO653O5CSPt2PpjXb87LWd3Kp36jA28HN0GA7C3gW3sS3yKS6uX42ipH7kEoDmmhPIt8G3LeGtdzDNcggrWwdw0rQS7XJyua0Z9rTixWpMbvgVt/00n47E1eLt7gZ26L0YonD5g8X1i78w+/b8FscNhcx9dQr7YuobbFI146E11SzfOYC2RDxAZVpq/C6tL3X8qxf91Ot+wuDVE1i7SElr5+Szsq4sjOJ1EbP8OPuzUleqzTJgjZJlK+VaOi8tBtdEZ5FZFi0Vfd5NTrTmkRvd7hTP7SX75x8gh8UJ9OLFCcaaQ69iRyydHB8lqRqFpXg3pTUiGWbPjeLKn41hSGERHrmi5rasvYmN3/Tg2EczuOV+coTAYTTLo7mYToIHB6dQ4tLAWsZ7iGxwEL+808L+TRxJ7qR8iunHa8iVDRHk0tC3eEq71zrnZ8nd691MUZ+GTNlHY+dlXe/sL8dI9YwiTOryJmf9F2HPOT4Kp0+SX6il+Emhmbzb18EOtytQ23ycWKZk7J+c9BhXJ2PCUpSs9V2QsRd68t28t7A4JxgHrtWjW2wMdoiS0X74DVxXXkYGHlbivWXb8d6q82RuQgHm+R7H7rLFaDhxCAelR/BlgQHjytrxdmk5LuBNk4uSM5iXacZ1n8ipMu08M/uPAnpA4EZvPu8npRMP8eJWGb21PAkfm+T00Uag12r5WLDTi3r7+lFNm4nk1ATQckMYJZWTZEuaJ40yqGjOanfEIANVv/kxZnpoufqAFpz+fSO+W6nlzpfswGfPW3FZ2x7OWi/er23B0JNNXMnjIpwf2Y/ddU3cvv11uKvjLDaO69B59AJOVF1G6+xh+bkXuGO2HV2WuRb/3DCELSFTKFtiRJ997Zj/cxktGuxnnq5ttFwpnMAKbh1T07Wbna7jUc3BS+TMaBzLVt1Eo+Qj9At6iKEDT7Cn/nfoqP4YDz/uw9ziQ3gtkcXXUxvwyfoM7lW9275Kxt3xeZzuuhutL49Dkfowl/RwGremvo9myUEuI/9rzNmwGTV97dymdiltWi/AYlNmL++9/cT17/7CqyX/BwkfAAAK03pUWHRNT0wgcmRraXQgMjAyMy4wOS42AAB4nH1Zy44gtw2871f0D0xDJCWROviwD8MTJN4FkrX/Iff8P1JFzbQ4CBw/1rM0m+KjWKS0ny7+9c9vf//3f67nL/v26dN1tf/zz1rr+tNaa59+v/jD9eXX3/72/fr68/OXd8nXH398//mvS0Tx78RX+Puj9uefP35/l8j1/eq36loyrnabqPm67B6qoet8qdDT2y36cupN8enXvJtOa/PoGfTGLWtGGPQEen1efvfeQoteh96855AwpT2NTum9ho+BE66vr7/9Itefn//xS3//ZlxfoRG+3Mf10u4us60Ot4ZbixLmTMU2zJpScY4ZE7HefQ0bRdGhqHdo1xZUbGOMKVD0WL1ajFRco0tDUu5hbakgTa23tuabuy/bX3v/aOEj4Uc9UIvbwqMZPorRppXcSkvFGSLGIrjYRA6YRvgziqJAEXH7bEI3dInFgqKPPjWKoqZiiPWZgdmcjjLNe6GEUgITy8hwkvfMlTVUWFAxQVa6/UVo0vEZwtc21K4XuZvEVL8ClpzhngN21br0BUSwatPXEGh2nX3VNLBsDDrGWFk2wB1OrXsMm1rzwLpNWFrqu27dFC6tW9074fNoxtY0zYy/Q8bvFTFd/yo6Vg6YBdphQG7YXn0wJx7yoStaKs41YiUuXOacCAJYR5WKIisXN7pnGfuiuXtQUdECo1rUVBw6xTIJgKGuIJybmZTQ1PJsmWYrE2vaHS0GfNnqUkqgPZNg6LWWjbPQbdk4LmjzanMkGkAHrAw1Yw4IUQINJLtozg3uiLGEDbwErQAo32K6WimWemqCtCI6NQ1uIKvC0NqsNmODClWHdy/oa2mSeXdUIKrNlRHpCEVKoNnJGsYWddfaCbYrhMZellDxOQagAofQ6cT3oympqRKx2IVzdMvQAEd0cFHczQVmcyGpDBUBqcBxR9krF9pWhBAwUhBAG41MEI38WBRZoBfCJwz9rje4ciVekXNdJUW2CzRB08gBChho0kV0GFIgRXHXx5Hz0YmJOQJE3JHKjg4vipsDkZM+CEgB8/VORTeIimIkMsG9pBo0BcgsspW8eaUUWxvrRlPkFtDvzOaU9YGme/JeI4sIyEnYNYJOIYjaWLNAuEtyA45WRIGCw7TAdYwuoKgG3jU1wfROviHzaAgbCPwLwBRNSxBxnqAHlAWfjWgDpSP2otgTGQgcXiGXrk6GQkcFJkyBUB8ZeQMc2Py3c+i17GIDNurZc8OS+GXWWwwUmkylCLye7dkTEyBCeV44mDEdOguBrH8IZ3cPPAfHMfCltoKa4hiFpeQo/yvzIcvQ2yBvTBFyCBwGQKXUHLB+TWoBHp0KHc2JiAb63aTaRAleOZHAQbH7jATHOYLGG61ENDRPX05o79GBoi9WHatH1DYfBlWnUcxvsiZgj4SsG8vFrEgCxF8hdltO1NywaMjN4IJhFcToLiqCVWwlFUZgYrCoFiNmSSey/UqWWA25p5tIJnekhRwgcwWbmDevbEXtsgPq6FpESRgvgK+qRqoGKsKtAfMAjgxqKli8V82Vx3cFvzZqBhi4kYh7n72kc7ZUBBmFkdYUbCRJryF91jEwd4kwl3SvRdLYPMnEBGIJfmoaDWwNexSD3AF/7BVgZ/BMIblpGRIqKD2p2GN0dB7GJ6FSumOySIAyiylJ7xqcllhKUKU6tSerhCEqLZh7ZfR7YCkWgSiwmzMRgjUTKMrT0cwcEPeC0VV4YXpqAraL7UHQAxmDOcawXXV7ZJHmPQIJzMrrwKKRIwEAqCNjrtQ0jEilpaZYD7keTwzLVfLpLQ/PdbdntXTNkbOjtTpbfNcI2cjuuidgDFnGNVtJkevuTCygPZI7sb9hLWZ+UasCELc0CSi0xf/PVT/HEZYmWaWS3t8gJyEcQtjtsLzv0d0qkJzlYUeEbEsK14xOoDWwOhTNub0cI3c65SxUhoHsL/0Qj79pMvDFFGFOj8wllt1WSu6RMEL/GpYlJTU2bLGdO1Ofpdd9vTXbiJa7OkYGGtBu2UR/lvuWGVoAI5rVsNQYK4/ZwqWlHB2SGRoAhtAQfsBc5VzC0lu3uNDtowEumvTeMZC47qF36sIV9s4daCuOrT1RMK8xOLzkPHZxQOgIgr3RMWNarsvO8I4iixNc3CQpHXmSxikARpx1UMZMKsSSa3k9QZp741XNPG8NR9G3Rcy0SaZsyBLFSJiJ1KgjFbGvYlUlwEDHysGPlQwzviiuPBqLRvA2BtMB6iYNdhovZ69kN9CfYtdgwyzmlDw4u3Ur6F2yNTlNfG9bkqBM7mSujqpm27aOvW3vmR1rlLBOrXurx1vOKtxYjHMcmgN15naE+4TXJK2elNWw6SGxmL5oQ+CNJUf39OroJrcmuDUMamJjhSPsXNxzWvVzbnoBDwJu1JwOnufoxnZY+2z5jig67mPc5rDhI2KsE4s3lqK4uW2SIkgrKCHHimOfCrLyUVxZy+GOzZJrBdYNkAHWpcn2qZfXzW0glZH8M7BsoB+JZQdYq6a8AWnkXEKxkZm8uuAi+2HwYzJtwsRFmSME4MX1NckYN5ZpH1Rto26BpBPHg5e9XEEkyI1FtWf0ZA7h1oUBgDIQ8vC1IgTr5z6f1z22G5CEtYopBtnXBubiSpsKiDjv7JiTGNX4AYS9Pty0myeYBvcH5yqnm+X/t5Ew95I/8N+Gqcmtz4CmyXZdPKCqvrEcRjUK9QK8KfDPxQt8ph+sSksH2J4c1djBLHAF46jGXHKrqpKqdABRvXCfxjrPpgeR24f6i6YDXDHRg/mIIthm97V51hH86/dvH96e9mvUlx/fv53XqH6hy8+jExB3jfO2hOvGNc8TEmjt8vM6xN8+36LzL+T2vPTwXnGt86DT+Xs9bzWdzytSH2X6lW9o5fWl871Eno9Et+QIjMdQVh5KhG8g8rxiCR1lv5d3j85f5AkU8Oj8RZ5YJf2Nqz5RpECe+JXp4x79xAxKU4qfo9kfqfekQekw9epjgVJRTxUGdRDlEyYu/0rJMcz04ssnAt4Q09TjMFZjSVPPR1lvmjoFFkpwF3vitvQ46lfpMXaJI2GK0/NySxYas3odFhp7Ajd6DFv2BGV0mVz3uMyRRFvPR50ew5Q9Hnd6DFt2UEmPaet8RY9hy57CdHrceUC5OQqN6ZHQY24hTx06Xaaxo0OXeVl/XO50Gcbs+WpkVwFi5xmVPsOYPv7wskJj/TQTfYaxJ6xBl2nrGKbLtPVkkBeEtHW+osu0dQzTZRg78Bp0ma8ejx0+1dLYAQbfZGmsSOgyjB2oTLrMJfNI6PMkOsvFQmisSOgzjBUJfYaxk55Jn2GsSOgzjJ1yOX2GsSKhzzB2QOj0GcaKhD5zqXqywYWbxoqEPrsUpDp9hrEioc8wdrDr9Nkrmp0+c0I/QODdj8aKhD7DWJHQZxg7kA/6zD8jOBL6HK00QdDnkCqhz9yunzwHfYaxIqHPMFYk9BnGTuss+gxjRUKfYew006LPMFYk9BnGet3XhMZ63cuExnrdv4TGTgsu+gxjRUKfuRGeSUKfYexIuCIJrVURvYa507rceYT2qoh+w+BpXu4wkharbKSMnV/3EkmbVbanoH2QRcoQwRk2baVsVJnkeIHVKpOUeSGO/GOntFplORZhtceRZRycskWWcXDOFlnGAauPOa4RdWng79//AAw/f/ovZ40c+DS5grwAAAXRelRYdFNNSUxFUyByZGtpdCAyMDIzLjA5LjYAAHic1ZY9kpw5DkSvsmZPbHUF8UMCjApZ7Yw1F1DIkj8nmMPvS5RirTXkrlVfESQIJBIJfv/zx5d9+/r4/uePP7782/e//v3j25d9fbxX3sv/83M+fmvXT/v7+8/PHx8/P97/3lv+7y/w3zv7O5v+9c/HfVrufHyuZ3dV3kc8o3fXq5/bj40le/ntRz7XinAsnffGYz1XVfWxx3762XZe91lx62Ax89iNgfXYr3qeuy/ru+ycw3r7yX4d2VuuwjuDK+6ufdhvJ+Lq8vAsP4TVcdM4ERwwWS5Hbj78WeZ+sfhu98enP7P2Lgy3yuPFZbb6yNC9SdGfHocTusWWyVC3CIMTd9vlxD7ue05waRKXA8LmjpVmY9jZ6SSyVpa98nnqxly+d5IqJ05U4SpJqR+f+UzzdgzRJIRhHesrw0nlAYhrcQcGO/fIQPjmpH72DmVOHj7rp0z1iBtXaay971Sq2Em49oxN7bDcMLmyp0XkxpBUVVjtbhXXAO1cpwDn2s7giFu3SnX4e42PJqOpVWzfjo/lpvXTF1qotlND7ii/xe925dHksVezsXKRuIgU1d5Y7iI+LEnGlB2TX64doqwjPnjUOj6RO1SRt6tAjdsCPihZ6+wD9E31avzbWcOIXbEaC7yY5AI6lri93Tte1Hnn5LBjXTcMK5ex3p4+5AJSQhRiBQUa09qA/TZZ7xLKbgmYKdMkcfaBZ6aSQ3rqY8MJoQ97hcTagPIy3Z+tOjXcgNt7kSMOrMoGQbGMqMCrWA/6ZhwcuoodQ9hjL0WTK+fSrvQSWTKXzhw/EIytAXNlF4Bzh68r52ADu0gsaWqCpvDy7OSnKga1mRugb+ytje2gq8ZToV+fCt/GhZ81gOQK+kpnNqxRpxGFC2jCd8IadyYgAQWySivOWbYuCBNZqAmLztqSm6CRXwgScOqWc5xuTQlTHhXLtsCj5CaBgL/npXp1Kr0+xMDu8oSNRI3wtLrpoh/T0mR4ObDpBK3z25KmfWvc00mp/3Rpan9VZCqcoDPU+U3fqm2RPoXPRVeSl00W3N9KmFKzTRJJVKvkDqlh+QCdiA4aZ7/XgUNu0F9RFi6ceqA4CImoib2myWPlhQOoF2QVa9oiT4/pHNqGQxcvhslRzx7TL6ku0gvBkQXsohSdF608CvxCFQrqXDHtqfbZ0hlT5gtuDJ8XSMWjlcv0Bn22R5IBndJdStQ3yDLuOqM/cVELDDDEfoFF+p/qFS/oCUOIY+DqmRRoT86yGKVa12hRJq2ldahyh7Xv7Zzy1sXULF/D3pb2BJ2HCJlpitlLIlUjVu6iLutU/7zx8DUqhcjm+wSE50QD6giO0Vgw1uiEO8QEHhggzWeMbLBFr994kMxMiWLgNAYGpr/Um+c9iAh2isu4iSODuvXxKYFiAWRrMSmRkK2WksHpm8RwYYNcMdLW+0QgN0L2rpIrWGozPa4DuoqE6qQgJA8gh7DupaDoCHep5F1gqhPREyxyXqzritaBpP7yRE+PrjES75AqmXrbsCQPgJc4fAZ15io9JFdgqNGFxBxJxy+FLpX5lAyI+riSCvEUIREv9V+Ju0elXchCFobuOD239J6CHBVArCIuqEsC5RVvgjI2RBNkXxMyQ9N6FHuut6JpGCpFQYcNvWIgEe4aUPTKm4WUVrDznCHDK0ch4aZ5Ut0UwK3HEprS0jtEjUJeMWC9J526GKqgPfHf8jHNBLo0lQeTMiiHfS8CXtI40/vIJgFbelehQKMKOOWJpfXEo/bjVvFcqRJ+Tkqk5EfyIvG3mj6iQfR60qwcXabZFT8SWCoQYrak8ZtyiJhITGv/SSm7K24XQK3CASi4LYkugG8JDrAIz121JH78lNwA49ZUv/bW6s2otNHISr1P6DDr0VJGtq71fc9LfRE+yoNoumDmQOl1QNkkVXEoJ1HzJJmXH4+vcVPS/E0Zec788c9/AEyHk/Nx3J1FAAAAAElFTkSuQmCC\n", + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "for ref in refcell.refmoleclist:\n", + " if ref.iscomplex:\n", + " for lig in ref.ligands:\n", + " print(lig.natoms)\n", + " xyzblock = generate_xyzblock(lig)\n", + " raw_mol = Chem.MolFromXYZBlock(xyzblock)\n", + " print(raw_mol.GetNumAtoms(),raw_mol.GetNumBonds())\n", + " display(raw_mol)\n", + " conn_mol = Chem.Mol(raw_mol)\n", + " rdDetermineBonds.DetermineConnectivity(conn_mol)\n", + " display(conn_mol)\n", + " charge = -2\n", + " rdDetermineBonds.DetermineBonds(conn_mol,charge=charge)\n", + " display(conn_mol)" + ] + }, + { + "cell_type": "code", + "execution_count": 57, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 57, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "fconn_mol.GetAtoms()" + ] + }, + { + "cell_type": "code", + "execution_count": 64, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "N 1\n", + "C -1\n", + "C -1\n", + "C -1\n" + ] + } + ], + "source": [ + "for rdkit_atom in conn_mol.GetAtoms():\n", + " if rdkit_atom.GetFormalCharge() !=0 :\n", + " print(rdkit_atom.GetSymbol(), rdkit_atom.GetFormalCharge())\n", + " # if rdkit_atom.GetSymbol() == \"N\":\n", + " # print(rdkit_atom.GetSymbol(), rdkit_atom.GetFormalCharge())\n", + " # for b in rdkit_atom.GetBonds():\n", + " # bond_startatom = b.GetBeginAtomIdx()\n", + " # bond_endatom = b.GetEndAtomIdx()\n", + " # bond_order = b.GetBondTypeAsDouble()\n", + " # print(conn_mol.GetAtomWithIdx(bond_startatom).GetSymbol(), conn_mol.GetAtomWithIdx(bond_endatom).GetSymbol(), bond_order)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "conn_mol = Chem.Mol(raw_mol)\n", + "rdDetermineBonds.DetermineBonds(conn_mol,charge=0)\n", + "draw_with_spheres(conn_mol)" ] }, { @@ -2959,9 +3172,9 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "cell2mol", "language": "python", - "name": "python3" + "name": "cell2mol" }, "language_info": { "codemirror_mode": {