From 1907e1d3cba0c71e23082a8d17d9d4541f6e79ea Mon Sep 17 00:00:00 2001 From: xinzish Date: Thu, 30 Oct 2025 00:39:29 +0800 Subject: [PATCH] add api for data clean && sensor_placement --- api_ex/Fdataclean.py | 207 +++++++++++++++++++++++ api_ex/Pdataclean.py | 207 +++++++++++++++++++++++ api_ex/kmeans_sensor.cp312-win_amd64.pyd | Bin 0 -> 66560 bytes api_ex/remove_sb_columns.py | 32 ++++ influxdb_api.py | 143 ++++++++++++++++ online_Analysis.py | 82 +++++++++ 6 files changed, 671 insertions(+) create mode 100644 api_ex/Fdataclean.py create mode 100644 api_ex/Pdataclean.py create mode 100644 api_ex/kmeans_sensor.cp312-win_amd64.pyd create mode 100644 api_ex/remove_sb_columns.py diff --git a/api_ex/Fdataclean.py b/api_ex/Fdataclean.py new file mode 100644 index 0000000..deb2048 --- /dev/null +++ b/api_ex/Fdataclean.py @@ -0,0 +1,207 @@ +# ...existing code... +import pandas as pd +import numpy as np +import matplotlib.pyplot as plt +from pykalman import KalmanFilter +import os + + + +def clean_flow_data_kf(input_csv_path: str, show_plot: bool = False) -> str: + """ + 读取 input_csv_path 中的每列时间序列，使用一维 Kalman 滤波平滑并用预测值替换基于 3σ 检测出的异常点。 + 保存输出为：_cleaned.xlsx（与输入同目录），并返回输出文件的绝对路径。 + 仅保留输入文件路径作为参数（按要求）。 + """ + # 读取 CSV + data = pd.read_csv(input_csv_path, header=0, index_col=None, encoding="utf-8") + # 存储 Kalman 平滑结果 + data_kf = pd.DataFrame(index=data.index, columns=data.columns) + # 平滑每一列 + for col in data.columns: + observations = pd.Series(data[col].values).ffill().bfill() + if observations.isna().any(): + observations = observations.fillna(observations.mean()) + obs = observations.values.astype(float) + + kf = KalmanFilter( + transition_matrices=[1], + observation_matrices=[1], + initial_state_mean=float(obs[0]), + initial_state_covariance=1, + observation_covariance=1, + transition_covariance=0.01 + ) + # 跳过EM学习，使用固定参数以提高性能 + state_means, _ = kf.smooth(obs) + data_kf[col] = state_means.flatten() + + # 计算残差并用IQR检测异常（更稳健的方法） + residuals = data - data_kf + residual_thresholds = {} + for col in data.columns: + res_values = residuals[col].dropna().values # 移除NaN以计算IQR + q1 = np.percentile(res_values, 25) + q3 = np.percentile(res_values, 75) + iqr = q3 - q1 + lower_threshold = q1 - 1.5 * iqr + upper_threshold = q3 + 1.5 * iqr + residual_thresholds[col] = (lower_threshold, upper_threshold) + + cleaned_data = data.copy() + anomalies_info = {} + for col in data.columns: + lower, upper = residual_thresholds[col] + sensor_residuals = residuals[col] + anomaly_mask = (sensor_residuals < lower) | (sensor_residuals > upper) + anomaly_idx = data.index[anomaly_mask.fillna(False)] + anomalies_info[col] = pd.DataFrame({ + 'Observed': data.loc[anomaly_idx, col], + 'Kalman_Predicted': data_kf.loc[anomaly_idx, col], + 'Residual': sensor_residuals.loc[anomaly_idx] + }) + cleaned_data.loc[anomaly_idx, f'{col}_cleaned'] = data_kf.loc[anomaly_idx, col] + + # 构造输出文件名：在输入文件名基础上加后缀 _cleaned.xlsx + input_dir = os.path.dirname(os.path.abspath(input_csv_path)) + input_base = os.path.splitext(os.path.basename(input_csv_path))[0] + output_filename = f"{input_base}_cleaned.xlsx" + output_path = os.path.join(input_dir, output_filename) + + # 覆盖同名文件 + if os.path.exists(output_path): + os.remove(output_path) + cleaned_data.to_excel(output_path, index=False) + + # 可选可视化（第一个传感器） + plt.rcParams['font.sans-serif'] = ['SimHei'] + plt.rcParams['axes.unicode_minus'] = False + if show_plot and len(data.columns) > 0: + sensor_to_plot = data.columns[0] + plt.figure(figsize=(12, 6)) + plt.plot(data.index, data[sensor_to_plot], label="监测值", marker='o', markersize=3, alpha=0.7) + plt.plot(data.index, data_kf[sensor_to_plot], label="Kalman滤波预测值", linewidth=2) + anomaly_idx = anomalies_info[sensor_to_plot].index + if len(anomaly_idx) > 0: + plt.plot(anomaly_idx, data[sensor_to_plot].loc[anomaly_idx], 'ro', markersize=8, label="监测值异常点") + plt.plot(anomaly_idx, data_kf[sensor_to_plot].loc[anomaly_idx], 'go', markersize=8, label="Kalman修复值") + plt.xlabel("时间点（序号）") + plt.ylabel("监测值") + plt.title(f"{sensor_to_plot}：观测值与Kalman滤波预测值（异常点标记）") + plt.legend() + plt.show() + + # 返回输出文件的绝对路径 + return os.path.abspath(output_path) + +def clean_flow_data_dict(data_dict: dict, show_plot: bool = False) -> dict: + """ + 接收一个字典数据结构，其中键为列名，值为时间序列列表，使用一维 Kalman 滤波平滑并用预测值替换基于 3σ 检测出的异常点。 + 返回清洗后的字典数据结构。 + """ + # 将字典转换为 DataFrame + data = pd.DataFrame(data_dict) + # 存储 Kalman 平滑结果 + data_kf = pd.DataFrame(index=data.index, columns=data.columns) + # 平滑每一列 + for col in data.columns: + observations = pd.Series(data[col].values).ffill().bfill() + if observations.isna().any(): + observations = observations.fillna(observations.mean()) + obs = observations.values.astype(float) + + kf = KalmanFilter( + transition_matrices=[1], + observation_matrices=[1], + initial_state_mean=float(obs[0]), + initial_state_covariance=1, + observation_covariance=10, + transition_covariance=10 + ) + # 跳过EM学习，使用固定参数以提高性能 + state_means, _ = kf.smooth(obs) + data_kf[col] = state_means.flatten() + + # 计算残差并用IQR检测异常（更稳健的方法） + residuals = data - data_kf + residual_thresholds = {} + for col in data.columns: + res_values = residuals[col].dropna().values # 移除NaN以计算IQR + q1 = np.percentile(res_values, 25) + q3 = np.percentile(res_values, 75) + iqr = q3 - q1 + lower_threshold = q1 - 1.5 * iqr + upper_threshold = q3 + 1.5 * iqr + residual_thresholds[col] = (lower_threshold, upper_threshold) + + cleaned_data = data.copy() + anomalies_info = {} + for col in data.columns: + lower, upper = residual_thresholds[col] + sensor_residuals = residuals[col] + anomaly_mask = (sensor_residuals < lower) | (sensor_residuals > upper) + anomaly_idx = data.index[anomaly_mask.fillna(False)] + anomalies_info[col] = pd.DataFrame({ + 'Observed': data.loc[anomaly_idx, col], + 'Kalman_Predicted': data_kf.loc[anomaly_idx, col], + 'Residual': sensor_residuals.loc[anomaly_idx] + }) + cleaned_data.loc[anomaly_idx, f'{col}_cleaned'] = data_kf.loc[anomaly_idx, col] + + # 可选可视化（第一个传感器） + plt.rcParams['font.sans-serif'] = ['SimHei'] + plt.rcParams['axes.unicode_minus'] = False + if show_plot and len(data.columns) > 0: + sensor_to_plot = data.columns[0] + plt.figure(figsize=(12, 6)) + plt.plot(data.index, data[sensor_to_plot], label="监测值", marker='o', markersize=3, alpha=0.7) + plt.plot(data.index, data_kf[sensor_to_plot], label="Kalman滤波预测值", linewidth=2) + anomaly_idx = anomalies_info[sensor_to_plot].index + if len(anomaly_idx) > 0: + plt.plot(anomaly_idx, data[sensor_to_plot].loc[anomaly_idx], 'ro', markersize=8, label="监测值异常点") + plt.plot(anomaly_idx, data_kf[sensor_to_plot].loc[anomaly_idx], 'go', markersize=8, label="Kalman修复值") + plt.xlabel("时间点（序号）") + plt.ylabel("监测值") + plt.title(f"{sensor_to_plot}：观测值与Kalman滤波预测值（异常点标记）") + plt.legend() + plt.show() + + # 返回清洗后的字典 + return cleaned_data.to_dict(orient='list') + +# # 测试 +# if __name__ == "__main__": +# # 默认：脚本目录下同名 CSV 文件 +# script_dir = os.path.dirname(os.path.abspath(__file__)) +# default_csv = os.path.join(script_dir, "pipe_flow_data_to_clean2.0.csv") +# out = clean_flow_data_kf(default_csv) +# print("清洗后的数据已保存到:", out) + +# 测试 clean_flow_data_dict 函数 +if __name__ == "__main__": + import random + # 读取 szh_flow_scada.csv 文件 + script_dir = os.path.dirname(os.path.abspath(__file__)) + csv_path = os.path.join(script_dir, "szh_flow_scada.csv") + data = pd.read_csv(csv_path, header=0, index_col=None, encoding="utf-8") + + # 排除 Time 列，随机选择 5 列 + columns_to_exclude = ['Time'] + available_columns = [col for col in data.columns if col not in columns_to_exclude] + selected_columns = random.sample(available_columns, 1) + + # 将选中的列转换为字典 + data_dict = {col: data[col].tolist() for col in selected_columns} + + print("选中的列:", selected_columns) + print("原始数据长度:", len(data_dict[selected_columns[0]])) + + # 调用函数进行清洗 + cleaned_dict = clean_flow_data_dict(data_dict, show_plot=True) + # 将清洗后的字典写回 CSV + out_csv = os.path.join(script_dir, f"{selected_columns[0]}_clean.csv") + pd.DataFrame(cleaned_dict).to_csv(out_csv, index=False, encoding='utf-8-sig') + print("已保存清洗结果到:", out_csv) + print("清洗后的字典键:", list(cleaned_dict.keys())) + print("清洗后的数据长度:", len(cleaned_dict[selected_columns[0]])) + print("测试完成：函数运行正常") diff --git a/api_ex/Pdataclean.py b/api_ex/Pdataclean.py new file mode 100644 index 0000000..a5c0bd1 --- /dev/null +++ b/api_ex/Pdataclean.py @@ -0,0 +1,207 @@ +import pandas as pd +import numpy as np +import matplotlib.pyplot as plt +from sklearn.decomposition import PCA +from sklearn.cluster import KMeans +from sklearn.metrics import silhouette_score +import os + + + +def clean_pressure_data_km(input_csv_path: str, show_plot: bool = False) -> str: + """ + 读取输入 CSV，基于 KMeans 检测异常并用滚动平均修复。输出为 _cleaned.xlsx（同目录）。 + 原始数据在 sheet 'raw_pressure_data'，处理后数据在 sheet 'cleaned_pressusre_data'。 + 返回输出文件的绝对路径。 + """ + # 读取 CSV + input_csv_path = os.path.abspath(input_csv_path) + data = pd.read_csv(input_csv_path, header=0, index_col=None, encoding="utf-8") + # 标准化 + data_norm = (data - data.mean()) / data.std() + + # 聚类与异常检测 + k = 3 + kmeans = KMeans(n_clusters=k, init="k-means++", n_init=50, random_state=42) + clusters = kmeans.fit_predict(data_norm) + centers = kmeans.cluster_centers_ + + distances = np.linalg.norm(data_norm.values - centers[clusters], axis=1) + threshold = distances.mean() + 3 * distances.std() + + anomaly_pos = np.where(distances > threshold)[0] + anomaly_indices = data.index[anomaly_pos] + + anomaly_details = {} + for pos in anomaly_pos: + row_norm = data_norm.iloc[pos] + cluster_idx = clusters[pos] + center = centers[cluster_idx] + diff = abs(row_norm - center) + main_sensor = diff.idxmax() + anomaly_details[data.index[pos]] = main_sensor + + # 修复：滚动平均（窗口可调） + data_rolled = data.rolling(window=13, center=True, min_periods=1).mean() + data_repaired = data.copy() + for pos in anomaly_pos: + label = data.index[pos] + sensor = anomaly_details[label] + data_repaired.loc[label, sensor] = data_rolled.loc[label, sensor] + + # 可选可视化（使用位置作为 x 轴） + plt.rcParams['font.sans-serif'] = ['SimHei'] + plt.rcParams['axes.unicode_minus'] = False + + if show_plot and len(data.columns) > 0: + n = len(data) + time = np.arange(n) + plt.figure(figsize=(12, 8)) + for col in data.columns: + plt.plot(time, data[col].values, marker='o', markersize=3, label=col) + for pos in anomaly_pos: + sensor = anomaly_details[data.index[pos]] + plt.plot(pos, data.iloc[pos][sensor], 'ro', markersize=8) + plt.xlabel("时间点（序号）") + plt.ylabel("压力监测值") + plt.title("各传感器折线图（红色标记主要异常点）") + plt.legend() + plt.show() + + plt.figure(figsize=(12, 8)) + for col in data_repaired.columns: + plt.plot(time, data_repaired[col].values, marker='o', markersize=3, label=col) + for pos in anomaly_pos: + sensor = anomaly_details[data.index[pos]] + plt.plot(pos, data_repaired.iloc[pos][sensor], 'go', markersize=8) + plt.xlabel("时间点（序号）") + plt.ylabel("修复后压力监测值") + plt.title("修复后各传感器折线图（绿色标记修复值）") + plt.legend() + plt.show() + + # 保存到 Excel：两个 sheet + input_dir = os.path.dirname(os.path.abspath(input_csv_path)) + input_base = os.path.splitext(os.path.basename(input_csv_path))[0] + output_filename = f"{input_base}_cleaned.xlsx" + output_path = os.path.join(input_dir, output_filename) + + if os.path.exists(output_path): + os.remove(output_path) # 覆盖同名文件 + with pd.ExcelWriter(output_path, engine='openpyxl') as writer: + data.to_excel(writer, sheet_name='raw_pressure_data', index=False) + data_repaired.to_excel(writer, sheet_name='cleaned_pressusre_data', index=False) + + # 返回输出文件的绝对路径 + return os.path.abspath(output_path) + + +def clean_pressure_data_dict_km(data_dict: dict, show_plot: bool = False) -> dict: + """ + 接收一个字典数据结构，其中键为列名，值为时间序列列表，使用KMeans聚类检测异常并用滚动平均修复。 + 返回清洗后的字典数据结构。 + """ + # 将字典转换为 DataFrame + data = pd.DataFrame(data_dict) + # 填充NaN值 + data = data.ffill().bfill() + # 标准化 + data_norm = (data - data.mean()) / data.std() + + # 聚类与异常检测 + k = 3 + kmeans = KMeans(n_clusters=k, init="k-means++", n_init=50, random_state=42) + clusters = kmeans.fit_predict(data_norm) + centers = kmeans.cluster_centers_ + + distances = np.linalg.norm(data_norm.values - centers[clusters], axis=1) + threshold = distances.mean() + 3 * distances.std() + + anomaly_pos = np.where(distances > threshold)[0] + anomaly_indices = data.index[anomaly_pos] + + anomaly_details = {} + for pos in anomaly_pos: + row_norm = data_norm.iloc[pos] + cluster_idx = clusters[pos] + center = centers[cluster_idx] + diff = abs(row_norm - center) + main_sensor = diff.idxmax() + anomaly_details[data.index[pos]] = main_sensor + + # 修复：滚动平均（窗口可调） + data_rolled = data.rolling(window=13, center=True, min_periods=1).mean() + data_repaired = data.copy() + for pos in anomaly_pos: + label = data.index[pos] + sensor = anomaly_details[label] + data_repaired.loc[label, sensor] = data_rolled.loc[label, sensor] + + # 可选可视化（使用位置作为 x 轴） + plt.rcParams['font.sans-serif'] = ['SimHei'] + plt.rcParams['axes.unicode_minus'] = False + + if show_plot and len(data.columns) > 0: + n = len(data) + time = np.arange(n) + plt.figure(figsize=(12, 8)) + for col in data.columns: + plt.plot(time, data[col].values, marker='o', markersize=3, label=col) + for pos in anomaly_pos: + sensor = anomaly_details[data.index[pos]] + plt.plot(pos, data.iloc[pos][sensor], 'ro', markersize=8) + plt.xlabel("时间点（序号）") + plt.ylabel("压力监测值") + plt.title("各传感器折线图（红色标记主要异常点）") + plt.legend() + plt.show() + + plt.figure(figsize=(12, 8)) + for col in data_repaired.columns: + plt.plot(time, data_repaired[col].values, marker='o', markersize=3, label=col) + for pos in anomaly_pos: + sensor = anomaly_details[data.index[pos]] + plt.plot(pos, data_repaired.iloc[pos][sensor], 'go', markersize=8) + plt.xlabel("时间点（序号）") + plt.ylabel("修复后压力监测值") + plt.title("修复后各传感器折线图（绿色标记修复值）") + plt.legend() + plt.show() + + # 返回清洗后的字典 + return data_repaired.to_dict(orient='list') + + +# 测试 +# if __name__ == "__main__": +# # 默认使用脚本目录下的 pressure_raw_data.csv +# script_dir = os.path.dirname(os.path.abspath(__file__)) +# default_csv = os.path.join(script_dir, "pressure_raw_data.csv") +# out_path = clean_pressure_data_km(default_csv, show_plot=False) +# print("保存路径：", out_path) + +# 测试 clean_pressure_data_dict_km 函数 +if __name__ == "__main__": + import random + # 读取 szh_pressure_scada.csv 文件 + script_dir = os.path.dirname(os.path.abspath(__file__)) + csv_path = os.path.join(script_dir, "szh_pressure_scada.csv") + data = pd.read_csv(csv_path, header=0, index_col=None, encoding="utf-8") + + # 排除 Time 列，随机选择 5 列 + columns_to_exclude = ['Time'] + available_columns = [col for col in data.columns if col not in columns_to_exclude] + selected_columns = random.sample(available_columns, 5) + + # 将选中的列转换为字典 + data_dict = {col: data[col].tolist() for col in selected_columns} + + print("选中的列:", selected_columns) + print("原始数据长度:", len(data_dict[selected_columns[0]])) + + # 调用函数进行清洗 + cleaned_dict = clean_pressure_data_dict_km(data_dict, show_plot=True) + + print("清洗后的字典键:", list(cleaned_dict.keys())) + print("清洗后的数据长度:", len(cleaned_dict[selected_columns[0]])) + print("测试完成：函数运行正常") diff --git a/api_ex/kmeans_sensor.cp312-win_amd64.pyd b/api_ex/kmeans_sensor.cp312-win_amd64.pyd new file mode 100644 index 0000000000000000000000000000000000000000..8bbe8006cbedbfe5599fcb7575619961cef48dd5 GIT binary patch literal 66560 zcmd?Sd3aOR_BWnF0|g8TNVK9zgyJnCQmJ}dE091dCz?nV&^iEhK&zrwOtcC}p^Y{k zf(6BkUcFYmB93sKKon>Jp@>Swsft6@t2Z^Mb%Nq_-p^<4lc9t6_x*jJ=l$#D(e9kR zhQ0P$Ywx}GaQ5+EP;2RFu~^dazhPJ`>yh%GTmJs{KaMnurT@?k{VmV;dE=P%X}&j( znO1h?4X*OJ*IhC9;%i)&UVQDf*9Bac%y!KUUhBH@T36wesjh3TyKMFeHe25;ll2qN zZo25(x0~$o|NrykB>O%*Uvu7D{h2=cq)M+j@9q9yOa43kzmoLb{+}cL6?oqu{r2s* z?AIeTrZ@NhPSPeb9b;F|d*u0?D=#f$S$|fu;v$RXvfFxDZkyQS$~aD&Wtgj9kN(G5 zo=48PD!2AUqz?IA2{owSJuH@9l9fy?rGQeo_}4=o6p=}b+|JJ`r|l<;We%Vh{S0|X z>`t?^GtEh}RKyG6b8(vGBP4CD7Ryn5{*!-eqrxHyZGx>$xssX5Unc}+R|N2=J%>LO zk#+Bc@3L4*Pndhz#es`0mS;}^96T*~_&@GM_{o26{C9#vvy==#1|25@X&e4e`RBG+ z>Q5+Fh*C$?Rq&+ZNkrhePndhd+)I%qAn3*d-CHb3Q~tRHZ}yz)fJhxtSAy503*N*O zvj3-kYT=8rGPH0Y%WkpSgPgb^)xv87|DKAaW+2G(42^PZ@J)e?FA*A;Wwq z;R7v5NHHI1TRv!fps(=3%?CKMJmm4A#eB%wQJ9s@z6_}giN$cf*p zR9di)rtdPG*VThu-{FhgpH&AuuYRG!tKX36)xVgv_4SVbw9uZQO;4+;PYd)2r)he_ zR?Qf(4k|Ix+*7BYd!FZ9=&ll+G~N2-A4qbyQPc{r-g;K#q%B3cJH2|FPq*G=O|#g~ ztM`R{S;hQQ{-b>%YUGVn7DH~5l;U>@T$#ViypftF`R&Up(L&$YtAC(&^zB+`e_*Uu z6-^5alfQ>*Rqbi^(1(Ef^fMg7UehPStfITYw6Qy8)pR79zQAF!=r6=4k~lCu7WC<9 zUj02ye=Ntjt~}oUQBA(=3{?+pmkbCh|NJSrZ5O6iPrvmnJY-858qt6SLT- zHv%~@4ZO6FVPCqO>?76zTLsmJl5?~Oi-%ig`|t~{@dgEQMTcB#v4rYx?z`c=eilpX z@gf(LDOGDd1CNIDoCZh__I)?(CXTq5ZhL1lawAsfb|mnR9EgQ7ZNJ0I$r;61mu(@6 z^A~2(h3HC}-p6oefOf-MWcD($3};(C08+xM%v>dSwz;Z==S*_;T`b}6DJ6_iC9LdH z!m3RuVF#GSE;5|UO!l@78<}&4R`o~Ez{Q$wbEp!|7LI|e9B?0(!u>5)JmNgH3wMp& zzXUIzZcBfYNVXS%>sgO1OZ1kNXtnjRA>uLmV0w@^j!zfwc{@5%>gPja&V4<4{DYuh zw5<|(o|&^|&bY8sgQ2gv&eN*)xzHAFM@Osazj)zBI^3o~(`%*W>XS2kdM(Pe_(OXF zR*$hK*2{4AsR!%uoUAd4P zxD2?C)_rl4RdII3A&VSh;DsKjL-*4LH)7eP9m~FwWpk+AcmgT$#Hvn8En9lT?#tc< zBcrN^C=NN`5S!b%;G0me;q)YnJrq^`VHb+!sJ!u15}~bGccU9lsJBjZK#O%z?|G`s zafJ?;^|%j%4~2e^>o$c-37GO8Nb?`CK^ zJVMjy8|cI|t!kUg?5WN1Yfe_Cbi-BqioJUK8Ij|+pcl%JYJ890cTs)3BSoM5J$o2* zF1{?yAKGcJ{+lRz2fCs_rs_#`$?4(M&zGi95gY2)A=A5(p}!fu{R)T(CKc=jyqd8~ z3+=TpeF;?PMixUUz)f^x6c&umH&r)w?v-M!-W%?r$}iQ>)GA{04d*>H4S3o% zEF%C--$;|D?}X*heLVxEkR6-;#0LIWT~~jW8F5Zv7P~H(3^JhW@H9yZHY(!S<}C!co^I+%hEjh?}^s=hL3H!h!{l~zOdrWv=(Lf^hzdvR2&Vr zi4AuI>x-_W)Tg)jz%h>v^l(idUAWU?*>o7NJ)7S}wOSpBg(1ff$l=vrM|mqzcQ5n_ zf9=9`x44!x_I{qt@00#gD=LtGeAax2x9a#T1PI__IL}&7u9Gsf0_#+mHBgElpRdD@ z^AJCHckD3u&zUD)-4Yk=@y3K`pOW<$u-0@poVTV00&C$V5-e$G*q4?`d-O8arYGSi$$>Bl{zaI}k%vih%d4R1MDq3iGq)AEH%drv+h< z(@-+>{K9Y4?#CRD;~g*T*;UiCp*5f02G!WBkC!s6w+QK*pmU3fPVbIxH0U(Y>Fq-L z8TwjJ{?YU8-l`MCHg6ivXI~Ua<{|OxKYJia{^=RPe4l<}rXPV(jz{mM>9;yq<({5! zFTJ;_eDD1JR3+J<85+*3LI09!C#w2|5l@=5CoM89E1QGWZmprEr#7mIzUzv9jy4Om zv5Y}Q`6qOc>1ZVkyMgK+4>)+%90OSAx{QT+r|gFxMMxZvmT zxi$HGB%eP>J|BT+=(!aUD+cZP^m=rVE=_l?7l|7+ePHyByTnk{^d`eO27s=e*&++33DQi1E8h`02CcxHmAUvk|ug z&~1ZzS||)sbZlvN1K@-IMQnrNEQR`ZkN$r2ZxNdKd1-U4R&?FTIo^MZ9~uvk*Z+n1 zVe%}<^JnqHeYEudNBnTu#b*5QGBlhLKb%X?koe(60+R8=yVO_2dMz8N89(%3S8K)( z$Fj_b)8MezA%3udZo_u)Lj16^fk?K&!1b(Wpb(v~t-_5B+A?GE;n5oFoF5TUk*jGG-Mvz}yj4VZDnN*QV_ zYhf+J9{MXf1Jv&NpV>?}bIC&nYthXBnuX$gdKAVKk^S`$d8mkNXDK&3gw>4@&#!Ou z<`-qyL-zpe)n_0!z$g^PgWl7t2W24a5$>t?QsZE+qWmHYqn&FkGNQhLc?h_8)ib(SK?;diJ-%VCR$w!xA*Q~kN(Fvyz2C;+wxfe& zfC`mb^aB$t7A{x=HML12=+@_R;7}>e^n8|bj9B$zl+qeP5r(r$m^M?F4ikDQp+5tP zv2l^_cO|s=?n%ldveL0?FuOBFzKSgUm;!WltUW?xJriP4YkFH?rWqsF~rM0+pI-irch;&<0iSmhY&mH#*mQ5!k2r%q^t2Q&Btv3O^zp z(`RnKi2OO*N%-VEVqXmG&OYLh0Dzvk?pp{O z4&H>`nR|&*4dIQc4i>SDMGTX6yo^O`L=iCYRiMI>15Dg3J!iKO{}7QjF$|ubWOWx= z%>b*;^TXu`BEY_%ut(=#?EaIUzX7h6BIeiVDNHf%A@4OAQ$AcY7u-w=xfBpo>_wd}Uz`>Qy8qMs;Dvim#PjSzOud||g$8r$+r41Ivs;rU`X zwD`#R%uPtXDoIpFqN_nv#4#qt8Gpcfc~@8Z2)R(%bg<{FH#}n@RK_)QF2i&)Y^RIs zHKrKhIF>zJ-0T^aeXa;>oeIWYMs@FK`JRj*77%aFLzc5Op3>UqBwzbV7`YpE_z0**x9ps71`3FiszqIzkk*o28zNQq4uK zYf;O<6<9y|{h#~2vT_UbRB>i*4|{wrK< zv$+U!Ow#Or6po&4dn*$abI!ehzZLx=p|?IJV6M?{v*=Ge^#bTNz7%#G=sKJ9^G8VP zNALN*;oK+>s*q1at7nh?f6p6kCMV6y{t_pntEo0Qz^tZWWe8+U6w2NM! z1*CI(?*Sp;`<-QrJxuF;)OmP(2iedDXLDD;aF5Y6j{GB%mU&JPtj9w{)vt2<=TGA9 zPgy-vKfg;&tUvvSSsZn^P^{_QjyjK-YEdz}y?QBVJI(2i7cp#)3OnT>XP1znvR#oX z?KT0^K~_pL?Y1k+7b!8xQp=@3_ULNJO%IW^^>!rCS$uksk2ygm*G+QRq;87O^soTR z*2TroP7{QiKrjA)&g?O@-dWTe{aA0kfgQS9?L>nxqCZb#fe=U*WjcoSCn_n>GVRsZ zrF9*E2LP_=@4-6J{^%3q^~x7iOzl$U6<^f>4oi1Z?rvnm33hn(KM+@7XQq=6tP!|1 z`zvyfkFVU9`x8HSAI6G{@pmxi2)>n|CbF9F$aC0`~@w{@eoi+Y3hIP1c z7bml4nc2s9UO5{-Tl+zWz?eXczX+jUFoIa__$F|qr)q(v(=@gEK{5^^B!?wCZ?#=P zV0SBdX92)XF|_)%3*+RTAA3^;~Mr0CU( z)*xLIJ@gi_>6+ZVFd4EJEQDF@)dp(R+N!EmwkA=nJk~Y~D|M*t*~#ovm~BJ0=apRi z6dgc9ARP}R`;s62z@jefXb(mCi{`8P2S`4|6Oei*Ih_b{E>CPzkZ+y=w_`=4!vJZ0 zT4leN%sz(M_o(d6$?T!Rr>%9qf;^gp3?{*BmAy2XeYi<*rh)_j$*o6ALsxDJ;3VPk z3Ov0lcp$+x1wJJOT&TC13>?IW{mw0_qB>d_Bi3BJO1xtnf3}^cT zz#j<&sOaCA^u*;M#}T3v3pw5^9Ch)T`rr^CuLJVt;{@YIV4%Z5D{EaQMFlaqHq$-m z8#K2vyMIlGwd61Fr=uBkvYy&GPfT!t0yB;oUTY1p@2N zp)VynO#c^cTg{5h^33$i@LV)&=7lr1a+Nr_KFm`IatgQQzkJaLV#<)xI1= zV|guzF0k3Gf9=D5b2;{xWkdPiQWt;l4 zl;z*F`^d|Bi?C!YP=oIOPlCnKV5vk2IyQ}g!aSRgBHtT9(s^^_@q2~uFyYHpPIVBb z+f~kYMI~`jhms|CVa?$1F@d(56*-2`VSsAk5@>Uvh%cm$O*tac3B5vBbU$*qwJ?|B1i~ z{n`+bViz`>5VTf`KP=aT(aR_rjBF;OlfWqE2Sw~-5@MlnB6hr%yFp?&@9c)X1$k%!Cqrs1>v2cA9F@h~q!zV|zUL>@v{A0u$WM*EqHJV*Bz1+2+JXHiuJ^3C zlQ=={CdFe4c%XJokBRfymo(8_I7Vp?ZRan}(y|~|>v}wO_RfAHu=SmLDDcNzEW?Vn zJ@l#&O4w?zCY|0Bi`?Hh;E_MRhXrxgMcG0j`f@o21 z>nukVDB0JPBc2Ifj@tA1=_sKm>ey^N52gTVk~eAHiQBR{$>B|nZm++t&acL|;TG+y5&Kxbd>e%3 zDOT*nofSE4i=PLlF3T$S>QM}NPkl&HT5+0{vk{Dg#W=azLhpibJjL2h*W)1OkJPx3 z&r_)->b$FHw!e7dt2oyh@W?^b<5SPF?ia`9Nklo_c|nGzpGswG9G={8 z`i+fiw^zT>p`o&1n|i*rBv#?)`Bt>!H=~UnC7L@vZ7oM>s8jO?$T!}^`BqlMME<*p zj9!B6%U=C&Xb8wpCXiWE6Zu*b8T;pQ#Pxiz<2-384`<+P*i5uBpeD|Sxp{C(9dYFZ z4C-=4Pn{$kR^JVu!U0-O1fiyzM_eHTswk%t4d+#KYz*doH#8H`Jm@FAJKj4YgcO_cB^GvSE9v)!Re*u`Oi zB|t!)bA67_W(hZ`64sh1$r3i93*XQLW<2NmIGt0Xk?qVuH1d1Tz-1^QTa|F7aGWED z-KM8-KSFUorwjK|a{m^*5RDx52$)A~RyMEoc@QzW{2E3ht=LyZE#N)$DL9yHzuLxu z)LYOH_zOa|AiT@=8#17z*;EpBGYvAuZiM#sqF;;Usy@{A$1+9MA7p`@fWAtoQ`0>e zHE>N_coVNC)}!ofvRrjraX10{}z=UVR%~K+fMj zbpr%iw2jR)r71^U5IQ= zzk_a=a$J6~xMg?}7U?GQG~K&USg@ZulXYRZckDd0q?c+*x!YS%AME8Zic4ZgA#({b zlR0SM=VTks_XSgisKe2d=RJh!hZ5&RigSv@Zq5M^N1x)_-{1}o(;KpmIxRiedy@h2 zVux0KcGM=ra9Ao2W~}GiPWD2~SSF*FndgG+j!v@k*;`fbC}h)|l+9>yfebYX30|yA zNN^e@cph3U%rD^3m3Y~phX`$(_bzRDK%#XF+oiya4m~|ZucwM$y=0T7f=~@45);+t zNoU#Z>d~*zs0w!Fd}J0SQmoD-t6|+pRDLH(tgY`bL5y2m@h9tG4+`&g)AYH*7kaOM zkREUr*caxHw=cN>MW{))3P3UG?vD~$CrZ)dJH*02FO{M}NYnZCy%j{TmkVY`=*jvhC%! z$xG~)qv_t_=p}{B;Ftg&1!BDG$)iwgW!)QL6%DE4q>NyJrq9b{tT9&at?5(a`wn^! zwGfU85@O#E3-j}nE1f&RAz95HvOI^lr&OvrpU1PXPZGd0#*nh;BoD>Gik*Zr6FNjs z3&9CWNXQwYob3>~b6w9RBF=yBp|5_cCP0l+pKd2W3&oL46^s`Z+(p(Ns1a0kN^Y==Rgforsa(SA6QjBP&TU_jWoj58i| z{}JbMiP12`JOt9s@L0617nO|XQAWgm*%-+RK*2C42LT|L+BnpD^!K7KUvAcs>vhM) z4mD5s^^^mCg8qA#MZD{D(mbD4$r|H8Bzh>M&u3wVd!& z$L^L|1`2a^e%T-TH*rqOl2ys%`RZ@)pl4vc5br?ot&;cVfNx?2=Ax}#&gbEb^K>op zR8~%#VXS8W5xsOCNJh2TM?Rlbk7xVHSF+Y3EwqoU&Z<*CaNbt^xI3#7Kj4x5cOv0n z7%I0j!}jHFBoXTq=+Q2>UuT_(n&tU<$PzWl-O;+ zCa3}&{tXUVxTq2sh*68`ko7%03qMyQs~Bk`vc6JSO~@)iRvWUO2-U(Q$ek~NQb;TrSyai)aA>`P z`M?Cy7pin?y#+A6_VOJ?E0JAvUmAYZTheN# zZvnWuuL8FqGfME+Y2ZBFf=34*2a-p2mO{!wZXR<_B~lTJDk|mUB}A$?N|07mBDs7= zh)oZ$5FEsk1_LZyQHNgutVKcq<2TTS8j+%BLiJWOApu4+642PMv;B=~jw~WCJdSH* zy=h0*DcJ7{)~sOlijI>G$nr0U1}@YD z(_PClIm)VSDH}q$J$vR{6zd{yEGCWf{_5QqqcYg;8N*+(HacuBcV{ zb!Ps<%m?LqqzEP;B{I%53~w28-)*WeN9klUHU;jGsZ9u8Fi zA8Kw9W}%hxWW7VWhftk@i5G@L^$NCKnXXa6Xb5nLD%fkvbZrXOu3!$vU=izJ&e+1C zOnF*lfe!GLr(nw!%&lM<3dS%LbcZRPa%3*D$f00zG9+R>PRUuT=xBcMT&w8rQ9Rcx zm`lN$73^W# zg&1%s>sR!L$f6u55s#yHU>O0*%vQ<_ak`A&NuCzvDMc46*sBaK&^sxZTfr(6Y?(4% zrGjY+R;OS;E7)oUD^{>Z1#4A2n-r`>!P*pTx-xpZf|V**CZ~z$j}%WiV-Q(XreJOb zo2}}g$>tWqY7}g( zf~{7t8&!F06>PnNH7VG!O3r2lt5>jg1zW1@vcQnMM95(D8_Ej5Se1M0-1^J1mQazNoW z!#TIa&?!%gT9nFHD%jh~ZR!*(s$lgBHdUFiQNh|2EUI9W6|7Cc+7--!UK{28sLIQf zr$rVFASgXg!DOBr4!IR9L%~WFY!_yBz${ZRhl158*eu1fR>3kAY^{RjDwkTXU@iq~ zRqVN@!_ys_%30O`56uufMI~v$(t-@GuVr(D=t=y#Y zTg?2=nGec#q$VYYH1N(q!b1BZ70+_1{2ViXEb~FBAvKL$tT0MUj0=dtGRswdrJ28& z`Jh~h)S;xUQy8mFj3q?>g>g|Jly z6O5)xui5fsJ+xL}ISQ7iV8seHLD{iH!Q2W~p4)x zT)`?7tX9FE)J4vf3RVdiZ0wEo3fA_KY^i{7-tXz7J9T)<(HfJbD0my z8l-kbS*tKsni%&H1C$S}Rr&R1{ubtgvKc93A^32M!ibs}UlRkA5r$5X#esy<4n#`I zY^3bGK$)X3@=T0Vi2=%u#VWtl%)f;BpsYY@mRYGVYD|nfh(XFamA}@^e^%vRXEk-a z9?!7U`O1p*^0cT%X}m?j-c^>0Dp;+8WuWPR?q_8UhdeD>sbDz@7QBO13gszSor0Ap z*p~`cs$i=XtWv?Os=OKnTdQEJ6>OQJTdQE}6|70YzEQAd1*=!Eb_Kgv@noRQ0aL+T z@?;GuFU?l4CI!fsAy1216l|q}1r@AL!J-OQuVCNnLf5EZZ3-4uu#@i) zSet^iE0_ZsKzU~?o|*Erh;4?RsbCSslkY%7PENi6D^;+Aif5UEITWl$!Nw_`wF;K0 zU~3g@L5*-(uV5|(YgVum6xfGM`0_UaV+jDM zhLocgNQ)IliHUImG3XuTD!GH-@1ZD z9t*Lip-KhQ6s%6c`YG6I1uIstMg=?XcA;xhuo4ApQ?TI*)~;Zs3YLjs2jxwH$5Dk2 z0C8dgMY{xW5tltfF67b{!;n%hc}PvYxbYJhnu&1+F(i^GRe*96z?a-eT7%T=SZWo< zN)zKgVvuyjS_P;#0ek_3q|HbnG(Dh2VMI-guZaQiS8x!8G(6r(d3g|&l-Wql{vk(U z@boM(P9+8_?q)vvredby6{SdOR$PL=*e@OB(dHUkV<@Oss$hcgy-?}p}u5y0XL`xvUaQbe<4#{col`Vu|>HLP~49$4$aoCLlRRC#!!32`wJQmMdxHnG$;P@Uco){NIhCVM_tJ;c98@w|; z2d7wZyPO`81$%52JUNXUDEDV^6F^MOt{P5~$gD$KR;w5*as9c`aGtH8&$hHQh$?itz~k>ZD66g@@b#ic1J zxj+z?vKGIgE+yqvb45ztzXB;6<7BU)3Mo?NrAj%77K3DQDLtlylofb@LDoMPG*bzB z@C6aHEnXt^sOgJys7{W|P_2cf>WkH)Y9aUg7P}}rs*g?I#r0Aa_v7@%YniX=M6e6@ zwp}k3TaIE-F|jKB=6fKfRE%*_Lj6O8VQR%xi10tG*o0@q0xngt_ISk>s#>fC^SGxf z)qht`35T>jPYJ7<m6@5jkm(WPIX>hcJgSHFv0vk%G@Ul}dJPNNm>Qoe%CXjY|! zf0U2#mz#gELpu#V^1#y~;jJhUZHe8gw55;93T-Jx^_-1E^(t&ps#hJnQoT4j_)_)d zK!xgEv>&pd`omPk4p$XxQ56f^B9+1eYR5`Zw~(DnNr_73{Vyuj_!KK;d;B?;zbJ#1 zI?#M%{efnq(L9>76{kypqt5B1>}KLTdBk?*vmh`=48D$+NV<8G7JRGKaJHfnnyxu; zkBj`Wtz3(Ys)I&wtQNNFe88Eai5lJm+KAuiGsnvXY4+;As1YwovtCO;)wY}&TeX75 z0ABdDWuX)e=a{5N!6x$(JflBVI5<$}1~7v_b#MY_3+31kj=CGp=4+5|zS)J@BYvZY z{&fysnwfy#%^i*y97SAg8PCBMIr?>>=HMh|A-h8u!NC-0XLmmD6qOucgeW^g7@|lvfO= zjRHiRAEFP^@Is*TXFyaH11Bqio)ZuRx&zeR2$ZV?ve#S<3&)OziqDQ0Q!K^UtIwwK zoez_B#CdzXs5`n8RT?j9X!k`mUM+RE*WhMmRn%9aKh)?;fY`Zwh4(;ku3r_E7cXik zo_S6lCtvCV1vqx~Cg=iARly$dQXfaA9z_*9D10h8e{=KjRvoNSyBp#~&NYj?L5kdl zC%iV%IuB3D-TK#GPu%FTXkUdiQ!C#~mp%9DtWuTJxB>Ob{iY~)KaQ5-JQ}4O4c`v- z#`(O<4kcsTMJzVrjG#{y6TX8gNuz3=M1-)d3QW{Nvp&3+)pi4(d9M~*w}$J}q<;=S zpS50|U&&&pVK^gKkzXU4N^A}BGVl^o(7qWk%f`>ZxOMEnI9p!#8nIml0SxCRlMI65 z>;#YOl(A6jDq*5!0P##gT#-batq^Al;sS-}z{$3`{26hMlA76Gd5OJ~8|4r}}3jFnz=iA z)W;Zy?XjMLYxFH@ByET~>@V_s!m&4_Bk{3QPO-@nR434V#1u4?fweH`2j$ z&BZsM3-jw%gTeehFz~ZzXEBFXU_1B(8s8Y?G)c>|d1(+Gmq(x-KMh8FGq%4rvVa%U zq3_&-vUB(Hy;x+9H*c%kC<0kec^YMfOLFAk(qcP~!F=?QVL*x+Ub92J^FB1?;C%(R zrc4GyylLkF2JdpAKPZgoQ{+~|_z>!{SMyyq8Tk$8s|lp73aKKF^h4r^aDii%Z#)c9W$Ppm1GH4^ z3425?SN$rxHN*KSp0vo2%*TlY(`;7Ini+t~Yaz-T>Of5!oQ`bW`ZYYD=MH=Isdx%o zKf+V%7~$53+|+@=z8eR#^jB@Ev1atW&<1^viJFT8M_aW zjMtFCGH#2O{Xic^Uoo5_aM1dK7R1XcWmz|+mIbjb8~r4I5#*@(*@?1><7MrKTA?BPWWHW2NVS|H+8%f|D`qevCl0v}bCtBn^(J4%6qd?Asq>ssJnRDtaC zV5kff82t@hgJ~%HX0A$XoK51rAV$@0Oq661>7d3V<0T356-2gnEotK=LVpVA#bo6W zG@5fvN?E;F*7b2gQC8nXSwt4kSOZs-S_$$EM83CA7d4)(%6jcFQia#ER`(1`DGAb9 zHZF(@if-S*3vf9D3f13wjxTa`=H-`P9$!zv1r^{@p25?ZE_o{+!^iH6DM#?o_-#L% zba%0)S29NOBP!zH9{mMVBA>hsiNHt-S;RuP zZ-;K_dC+=6LGu+f2S%%y^4c=y;^N%tC{r{ngw-Vct0J7I2%{v7*xq;)35v_{;$meM zOP(Ty8xwiH2J_=7ERPW^#DJe#HZCQ%*4~ueT2xAjaY4ZA-hABk8yt+hBcxijC18U~ z%<<-8DTIj)!WxXh()8to7|yX!lZYrxYsvI(FqIdO`+G-DFt0p26D6R&a&Jb&wg4O; zVjHFJUhGmO_GzhNgSTbleJo+D5_{>au3}#z!pU8LO6%lscBPcJ$`*fDI|MW4U@|Cmd!L*S z1?Ntdds4Bs3rjjwj+j#x`3&cNIHZ_c2!0?z3rdAUWhDq?$`E&`S33~=NP;zr;Kv6^ z#}Zbm5|l&T-GSgI65Of?wtxW9k!TG!p8XvGVPbd`aegN9*+A~3tv(kpwCK>V_Mb0` zyACYw9T| zo=)(%NUCBxkC7qzsWR_+WbmJKPNjUWw=9~**k8ri6yw7yA%?K9lEo9i?$W6o3o$xX z@la$L&W}YtF}7`ECg?FWYkdza;5OI)LllbI=ZOVS`@=4xeeK}|z}Pq#=zjf6dA%{>C?{+*Ls1lHYEBc*^y2NxTU%Gj6ui*k-Eu!mNQz?hhg2S)3C z6wVHlU{;U!e;+_C!lgMI+p)#~bv;>zq-tF3F<=GIVO#G2f>xDh2^>*`H)Zf*CqI8N zeHM%2=P!EZlFarB@G7 zx-FQ=$gobdZ7nA+yv`%*;|H4gMi1@GH+r}lc4R$Y-C=8$m1^AnE$eaEQ#&soSc|q* z02F_zG`?t^xc|h#72-@T54-TWBYfY*r_Xh0k-3@t_>=wt>|tMegOs*iE6B80SL3(W zz)Ev{hQDC@{FtVX%R?1%zwzn4yamnn>PrQAQbyH!s@Yw zieAIK(9W9&dk$2_!Mv$yS9B(#ExAM8zTvaM_Nv!`=hsUdxjX&nrj|ZS1`&6r$OG5M zs(SYOm4Y zt=a-5?XoYq7OM{U7SCV-#xOrW2+FL?$_3p4=U^x~cE z*s%)m#kl!c5U<+i^32?-SbY*-qlM+NSz+@^MPGq+Ah-qxd1Uj4wyfHgjopd#b5SJr z2l7(y(Y_PDX!a0u#S;wYE0SgNOVY0C0}$Yp(xH8tjs=N=*8;deV`4_qa_NEI(Z_H@ zP|O;2oXq;*BUbv9Qf&u=2CxU|fDTs#ZJEYG;dZ>G2ldROgmL#&_w$(F8#<`zZ=rPz zKtpBWQ3teXikodXPo;=ZAEpYnPp0E{#Q7@Ey)9RRY>x;`=i0sq}T!!(KgK z9M(Dzh`D=tnPj29wNT&XfjPB&R1x=t7lj(_)t{gmKK*QmuizYqJv1IvKHZav&nsp6 z3hD#6c^g-b^~Hx~u)BQ>fQ7o(RfIQL^$T1ceXFM-+5`W2s#<#Ezrq6CRP!77dGm zv+V&DZNk?r1z&s;rHUu3thhI;`bkG^%a$W4upnD*Gx;A95Gmufu%)Fs3`^K&Bv9C+Buj?DLdh|5CZ$X1SR7svyzry$&ieP+%Jv*bP zPlPV^hw5=B{`S~KFh{|!s4|Q`jlbU(52z1LR%Vxrub(N7TOI=&ZrWRgfElB>+S^OK z86O1}RNZhkk%t@yjbmYBjI4mLciU%QvNX}Vm9Td!W$#u-7l;p2v$j@zim%;#V~KsJ z>FHZy@9>%t7morsI(%OMcb8Y~%fn@K9cWbSE~o>AK@UD8~hqR&*zVvR*%Oh zFAK&k2>MEn0gJxCpY-sIQosJD>d$$jLuOI_iSrNj@_Q|qMeQ`4L#4WYT>1re?7#|v zkocfW8mk!&a;9~;{|xBf4gR4Xt*mVRIkq{wTFU0@)jJRcVX)W&ayKwJ|WjIk?X@Lkc)Gheu&?c zT)}fhroDlIYWzG9tj8$uIEFgx@I;($j(OPE2cb_ds*qvf5In?=MDAVq0i0>TNdUh) ziv+Mo7n|D@MC%Iq(Qj^K#>V4zJk(N=J2~c|XHee&dJ{;5b|2-g1ynl4yhJmIUsspB zuTI5l=Ez9Pwrlzi2-GxvbtT;cZ^UEMv53{1ioo<8^yh~2TZoMjIUO+wq(dkC9^((} zax`N~}lxvpfSO?`4N|K3;>ic|@cOl}H0oVHM;A)VMGU!NB{J4kALjIUbGdhvHgots$;^)H<^JH)+!S_%hAJry97yNahslegwK^($xKCz58vQCova_{Aw zoaoV`M1MF0G0A~Cr3WtIYBg?lknk)ndFuhXI6mzboCH3XC_cZ6WG0{YjuSpvCZCJ| zf)>olx-eq&M;$nAKCwvPXoR&1!*wcFGg{0M;oO2E;VXl3_a>|zui(}M5Bbs+1^7v0 z`;O1uar-ZI8ga>)iI1_uhN_#lK7dTQ@B3UIKH^tID|irPZPDIp_Z2jI>yarwSm^b{XHE zh=N63wPXwv*wuFuv%~1{8J^*L52ib02vcOyoP@Q${u5n~SgN2-^E)~{CKA`fBR<4O zZ_X9b52}a;&!1?RD;%#uj}6fYRc&e59oM&TcO0#_C(ays9QLXm;$DRfnXTW*uWQRS z{l+q0$mQ0m7Gzig`QovU&_WU%n&0x;ZA{NVjDS-s9d7^FiX!o?UT)T*vik#9YgM;q z@-B_g-kWi$N-ruP*Jwh5MYf*J^2v<1jdL?&^Ps9*g&4ik3-|&3vtzr|aQ23}u>rzR z1&HU+%G+AO*^ar{BIqxo1BqgLMRDERNAiCL+NKUTkGvTHjD5pRhx)>A2w#UKed&^I_{! z>awA=Cmvoy(3=iu#!V-oC4L*k4tDh&)jHrkw-_z@tfLsXHhgik)~`8sL5D}}z}^_D zhZ6-BrJ>lV^3+#!E~}7zv>jWa9wTjC$Th~Z$mGPMLG;+w9OmQ5poA@TkNc= zmsX>)mIdfao$ag^wYt0!S<&xtmx%Jf-~_|@BTB=-j-xzIO}cd3Lnr~(fMwx8>3GQ9 z&n~fqpHGHWOOw{w3J4ZP@{-GHv{Med0aRpH^s$46A+Fsvf1}(CwN5-MT7<#(smp68 z++nfs|Iua|e?B_pIo>U}r=BazFK&F8%%i{M)f+wf2hfjb!&EvP#xrv!8a<{(r4+23 z3w#yahrqGS`B^!Oel?(gn=(8+yJCrD;u#s9jaZ|{UDTYs71>YUcT;*%WVTVLSGYX; zcW`of;C)Zke=zmfT?Y(r!*`f^{K2V5MIO#@B}Pg8J8#wY%y>NN$7&B&WElKb@TS1G z+qF^XD0bMF@YJ~n=VtT<%wsXv-H#d0fc*`aym%wUn7rUGCNEySIL}k?YH%=#5s|#z z!nN49TfC7&dPko-8PX4uPFU_?L_B~!4KMws=dVeL%ozpmm$NUC$W!0$V`#$f-Boa= ziRK4Y@k4l=`25KpnAE} z;ej`XWz4{)DGyIU0RZCLAaL-1inF(KsdI?92x2RY>#t5C3IZ|qC(?_g{A7_V4^v;x zch2oF`St6ij@Watvxw}ds<;63+tFgxzsPX;1z@qG4Cf6bV@tn ze4Vs8){}Ks4qVJW3WJYY0g zQ1i;mYrXnu7#gQ!nB$~<3D>R=OuSyrq}BlvO;?vAp{FXa{u`PJ9-de&Hk^B ztn*6o%JDF2qHetK?nJ~1ES5_W#cBM3-tKb}tGT1o3h4o?$pf%;Cji`%P4vwezvT+d zY~D2Y`S%#a8NEn+nuuP8qfkz(s@sI^b=0mtR(c! zt)uKi-FP`G@qOAeu_saZKKk_Zx03NBHWv!so&Py+4)o!=VFchD9VWXeeZ|f6snnk0 z&T%M+2Z^>}W#k+ClEXPoL<_k_s18C5iAFjG+82-r;qAHTmsU~Ih;xiY73h)9g_=4S zavbr&w98i|1^|E&vA;)M^NABnz2&aAi_js%g)!Euo8ObS5hqsjUw=Khe#ePb##A+` zv5PA4VQXtG)aKP&F)h#C>&3m%jI;92h#b;iYxt}u_7Nesxd-YroOL1y`u!5Dli*7K z83>6m4H?2yL!*8}OS2C}jnnNKPkLkKUG;Su_Gj2l+#{B1hu;KY4V-Hl*OirAG3mp{ zxU>C{llHM4;;ZZ!)%jH#I8M$+h>CYP`hgX@(x>sr9XZy(&+nS3YxaElcm7e|CeMrU zN<0@P0~y0tGjLU_biH|oT=L9i{f=?Bc}rll(^oo= zETB8=zZKd&L#m|Ai>C4H13emVrA}H0Bb2${qT_x%tW`ni2{O6V z<)nB3M%5{AcoNp!>Vv~MrgL&@{RX+JwH$?6ux|aBSOQG8Ass}Ll{z58v9hs5s^NRz z_Cn_>vgupZkcK1S*ca3=-zrwCdmU;t=2Ut(m}PHG$g7$1LZ@Db4C5zKli^7kCB=d~ zB~2AZ>naQcqdYzc6WS9Pk1eBQZ{om5K$Dk+dNpsj7@RwIfh(|+Yf1*&ehN9!jMqyL z3ftZZ#m845%q=vYE=utzPcgK)C8a(7Pvt3wJjGq*v38JWijt=|E{}odPVyx5={EHV zlh=X?6zTpGMKUl>u4g3l7d!@4cS=#V_sD^vRis0S9Te4OA;_a!+KaX%d`58M`uV_E z#>Q!`D|Osn!<02#wUWJuqiPFX8UE;TA+IhquY*U=_^Vl)gd)vKNxaDJ*Toef{UEm! zW@5KyzG&3uq?{4w=?I|AD_W`L-mGCy?l|DX!CH%5hr0?C=K$NBNx} zY!RfalZoQEq-DJqyc5cg``Dl9V~z-7qrXGG?T{};9|y=OMIW>S^dU(XeeA)c7EQQ= zKKg?}Cw;g?AHPF+DDT*~KAwgh|7ZG;!n^3>0j!&J(gzVmA7SwBs*gl_L%>_$oDH17 zX&nsbLPOMT)0CjjxS;7kjN2@Gt7INO zV0CB)l$g?8!1oF7B0cEib9El|Hj6tS~rB{|mL zA+plaFzK7V!E+2t5!QX!!Bqg&X?=uelX+6BCPy{2LrWawIp)R}KRLds^$K$l)P+u- zyX@X_nt6~GzPdbrWOC*^C1sMQ9;fE})dpXBdU zlvfl=T?z_+r=r}eP~we-T=+W`x!y0w{R^B-`XQ&oeow-m)TuN$l$XmW-f|zdJE6?}I7HJN4J}cg<|RU50FNq?Bx@$`%(&$^HWs1hzO&N;ca!vc+vuvYS-4I7~|RODbDj zB_;b|l`RgElFjudO;5N70#hY=Fo@hb0#VOo?kwcSEu75tswHVGVC3O~4KeR=x%>0Q zV7Iz7UCf$9cXTjl3k}B=ZLFDrJh8a}JK~60an}t?KXa$h5;1V+oO+5|IltOlVeQm zP>l}W8J}OB??~h%Md@N)%PUxlGz}dq(&RlidN=WDDcd8)iZz$aiK0KlV z+8I7|GxW^)!|@=A?JC0btv`G;EqH>a_w^d;kX|NN&9T~Mn@kNwoU0Kg!!&W$4<%Vc zlbBq3jZ3j8!TI3-P~P53-Y6Kwnfj7UyOno7WKi-B6nS$YFWe_1Av2XBLN5m=2z_Kn zq4$kS=mYy0ajBC%!v7EXb=w$Kgjr4ZVO;}m|#)E6gpl!y8q1h z@1!Wn@dhiSsT#$3lRS8Vasl{4mK9!bf_MRU0`7$fsQ42I?c-YAf1FwCoC~)m#TqPP zsksCIoi!;-Td=pZ!nZxfRLTpy4id`zqD%jX1+O+p^4II=k6wI{CrGtmiM<)!n!S2I znq_#s^k(6TY%N@d&XtqEigGRrVW`6OO35fKJjsO*-nj8Uw$Ifh4QHrujf)n(p%}SH z=V8d6RGK6K6hl}iE@k1DE7^gDU#^h^gG%^i#w|$dR3cN<@XKq}69*-H!h6wYclJP7 z+uAcwR@IsoxTvZ%J#bc4>n#EIqJtJd@UFaWn-2{`>3&QTq}{!P=T7VUk7pCMBQ{ox zj2nS@9KH$`8Mg$T0nX<~#%;w4I!|?t>nls?IC44E>0j0N(0IVT z7*X#hmjj*Oi@>cmJVwv@;5%qyjSo5MUWwg-hpr*dSbT3_I*L6G z#mY4TIH|jg9uDXsfO0Jx=RWl1YXzFGpg6LzHFj9%8No#;kAplff5?43 zTGR){;^Um2nX+t}bWzWA*)Ph&`N=|#*S1II5EI=WKI?=N@Y{GM4V!DITQN>av#MSl zJ9&sZvhn4q9Yv+02W;b@p0{SGCE^}fu+(Q}vu1r#=E+z`L!eN?X_m`yZjkERc>fAs z9ciLAp=ErAc48(@02@dY`Bu16~rU>^30KFjA6M zRPf>acL9tPr^Vs+vzoYTfzx#zIy~tV3(>wZ&HZb%JhxXr=+WEsuYLMwxEBQ~_vuks zB9nV#jnSLY#PM7ozybVv&jD^wV}E-H~_C4M(%#`!qGEw_jgSilf7Rn8u^;_voK+ zi*G2lATknts8^rW$%P@aw~z;Ut-qrZUVIrF-^s%)$8h?ySSj(z8RC;CNk$8s`eB~@ zUh^l*rlm(VE$fHFmrp6gdVD1$Ase2$c~%^*z564g5s<)7g5&9CT{efW0~GSBQ|x@F z)xKl_W;^`)vOUC?c(8}jIvXh~InbLwOF1cW4d)l2rIY%@3rg7ocJxoDl-o^_zYU{6 z-Ru_e-Ua9>{lLeu-soXjL%dD3VlH#dlg!x2oBR^r2Q{20Lru_;y_$IOFBpZbmjgXq zhIvJv-g}3&7*B4%cUULl7pfsQ_LSgrbJGg-lgA>jP=6KQPC~RA{SBjo@6vW0e8wrK zbWg!{dt?;|^5^8@Bv5*;U42PMC+iLirpH5UvwXDme_SCO#Mo8DkpN5!u%1Fe)xi(6 z@5XRBx`;XvzMIDT=`nwzob3R}sTB4b7^~ppXufXsh3Xkwb?cSaqcR0;!IKLk75Y@Oir7+(&(zyZF=K@rB`&#}OCxM}qsra*+TVXB)f)@862B?=@_prE5L$W(2k) z6Z?^Nho=o!HDZ6P%{~Crw1W2+?lYAYhK@W2&Tyyadi1z5P|jQBqHufl7Z{X1VZ5b` zG(LiGD#Dj(?uU}Sc;BS@F8sg_a|z4@!v>D`Vg85pFhBYn9E7Q7-_Ey)uHzB!s}e7i zZJGa{+ed> zPQImPujckO7hMyt%4Hr{yc(4F=MEf#Stwbw=D-isi*gyn=UE2M`fkB6%W#ODCwH=! zHiArxTmaQ?t@RaT*sI@ShI)4nudnf4zLV#`zc6g5;5X&{7)|e&BI$ziq?Ccnd{Fva zZ&HiLCf}pes6){g-laRM!@G2&y?R5^9~1UYzEAfiOdRWpfeCg#)QgQqyldbq_|(3n z9QBo>zo;m{98#qX<%mL+=&{dC?O>^*SB zdipit$VzWe-2vPVLmFH?McxA2y%*ugghOj+a9n<%!7zepIBmMq(Ry1?E{xFLn9<~O zRU?~!4&&4&#qQB-d!hHYMkgS#zo=$Bd@qm&wZkr$Q`%SHwMIN+t|?mMA|DQe7XI5`7X8 zVU%Aii4H&(Li5*LOl&VBiT1ZK77C7x`W3nh@cR^8Cp>xxd*h?Pw1)*QIo7yiSM=H5 z_*HELCq((Z0rr)ASr*gd#;74;1~9?j&M*C7ym?d8d*W)-!@$gAyk%HVz69F`6c8H* z^=33L)=j?nDAzqkp91e%SVZC+S?hP-aZtDZd@%)HXj8;)s1Iw@MkDTP0w8)OswtXr zRgJp}U3-sUZgduE7Uc(dq913Ve;aWQ@-3)F46^3#ao@2!YDH1@7n?<|2(ecw0Y*`P zTC9s<7E$8zVO|!-0#{hez;QiJGu)U==)b@`*jTgYJ1;`{3waOTmtY^6XkhkzlN&RU zYIqs@H!Gs6SO)I-!n5AYcFQ9z&j1<9QlB^GdI&&;0KwGkUlX8SkY9oL_Q5IOHx@7+ zWP3CEH8hQA0nQZNjSd@SI#ii~{&5_f&pjG?GXmBsqi6861lzx+l*o6StnUFWAuUd< zD9Lb+1o=wpSg+i@lwQrnDqsMdx|}srCJq>EQ{KPew&V9`6hUkB@31%qEUvM^A!G_i zqjPnv_3vtsNvtRHGPKD@(6O*)bP9eLPPaTZZs?P4VLNhNh}L@&^LWw(rpSBga45@g zW+M+nQE#~m?TwLe@*7KG_s0R`JQ;Ii>kODl!)Xr`j#sS0)^C7-dY{hhuyx%mh3Urc zmjT&fU4|SMyThvE7pj5073Qz?@F$CQbv;{Cgp@nd)9CfdeBCy^1o>fKx+|WEttne& zJQIX-&@V*|=%5xg%iD@`6z9_0;W}7ZADRhj{ke89%|(w`&uM6>tku~FL0$8;s>`zS zAOY^8)pS68oS*vmwaM+wEayWuMhn=?WjH^9qKq9?MvwXbqUC1cu(8?3GBGvXhk`M| zwbq>n<}Z(Z_|ws{gi^3I7|^O_J5Y)r)1lA7L*t4lzJq8uU&j+pK%U6o>(%dwKT|!| z@wo`E1o_+W5mc9!i{H6LZCY-F_7b{y`~#=qc*so_rP&J`p#34|u@Fprv;$3eJd2NL z>zJUT-=$MI@Xr2Uj)G!>H~{+HVP>xa9Gs7UTK$)Tz!)(e5}*%>*v8?xwLc7o%Mv~K zy)}`7UlfPc@Nz{)R`e_Nqld|du$R#fK+(s{#)+d5IkdR0WyGw!r*%+neVb+I(Qq?# z_A3;>IIlZ_SB*CzKr}GASUvQP-i#lyeJS)(M)Z;-(|8}+X>AHqHFjWVqpi_qbcN_r z4#Nwi!2!}|XgW^Y?=YM*2eZDoq#5fASXYsj+hOsgbq{Ezu6GF&y#5feWuTL1g<9W0 zUi3r7-5PyWP#F@bdTunITc4L77*BTwdvCWkpe?3F%X@R%sR&#lA z{?F`Isu~S^@ZRbH0(|cbBz$Cw$EirRTZjCg_Pz!>j^ax9j{biwS(c1#Y^M!YWOlKT zWrgr(kStk9Ae%v!f#oNy(M-!4Jb&r#u_Q0#iI@-@;jF`96DI*ye|GG^*?5zzaSTgj zn+<}6eKH}elkgG+SXPkiuEHb4Ax`|hTU9;NGuYT?^WL6)XK9_*S5ej7Wx2n4C z%;3LAFaXDk?_!F<9D{$N-!c6iKDhQ)eEj(G>0bQf zDBib#xMlE>X%ySC<1hJV*f1Bs`kksckDq|`s~oGhV@f;~L`lns!M%fjt4-KGVSQi( zv1Eq*E91gGHzDd`mEJCXMwA&5i46yLajjx*5I<=;Jq0~O8t(bkwMHbq{f>9n6dv2Y_Qc?Cw{JKZzpZ6!@!#ACX;sc#bE*Vi zner7xXn<(dc#encWve=KGRR)CD#oxNd;f1wVgaS~E^UEh2LGYGC8LKed_PY1>vv!& zxpl(Rv_1Pcg1^W`+W$lV2Vo6^zgjbRYz^LqZT&;;T}H!>hFs~;qTqTV(qT^;#anIw7k!aeF=GxLy{j{iu$oxziW<^#82E0F^&EAg?F9={Y9zn3=(7%*h=Y#0qv93ry`YeLm7j2US=a=BX3`_|OfiM(gCS0M;(rR~{Q zA|$?f;*PBo+Oye95ruN$!379?7gUpAn?D!abRl+_Q@;m;CvV_3K8LR?=L(KoafH2p zvfo7N*5Y5_e?91avVGD2LJ-vd2oSzSe~&=Pp*PR7idsqz&CJ82ZzH9lblHn3tm{8nk@DZCKc`t#iF+5Fu>E%!XP4AN@(&9RmXi~B{{#U=n9vbnHj>(9T} zvc0%V6-nzD3c&ycb^t!?FQvUofmKp3O^j?!Ev&1gMuX%Vfw@%HpWlprML^aRd_ST=l zu^pD!eYy6pTYp*FhKZvMV>)l`%Ad4h={T?LPbb%wtp0(~Hn^~|6kM*o@(0%3?uL6_ zBYcg**9gbU4d2h?)c8!>Ao4=AeAZQZ6zlC--1?8>|LN;HHRVos!byHkRhGg2iEqu< zS_a>~e=QX18SRd@0{DoX{4#x+tWCEgg%`sbx-mg7iI4uqA;vQG?G}9d1pS$NWj0L=2Iy`t7>5#(zm)BDMjt5*> zSrX?dtX0^ouuI{9!iN;@Rrrj;pDKJ;VSy^YUg1p&9btV|#rqU~RpB9pdtul4dQ_p4 z?=|&&Mq!3}zEI)S3SUvrZ&h(qp;OLNzn7us6~3nM9fetckk1z={1b)s3cU(HrO;G( zpTb=V|61Ww3ja~zw8E1LXJIeH*HVRb3fC#@RJcXq{R$sdctGJZ3SU+Dwn90l@NX?j zx&Dbl{64#=w?*OY3i}j3tnhIqmkAYrMIn9!UgUpY;arc5*D7pM*sicg;bw)~6b>o; zj=~9ruPQvDa54WJIj&lTO$s{{b}Q^txJ}_Mh2K#4w8B>uzOC@2!qQAxPmRJRg|{e- zD7;VMmlTdEJfiTY3LjPSyI*0S!XAYk3a?UFt#GcwJcaL_6wRmC6nfzo5`o z7*M!QVU5DY3JVm@{7#ntV}%n6#}q!GFsAT!h3gbvs_^3q^Ax_v+Z$Z3DSS!c(+Y-+p`QgIK?lPuGKypsnjpJeBUoCjI)gTw--P-A2RK?~F(M zb})+i{L!cx2u6H1{|=0}&Q#9X@>xd6@Y{Mg3a+ANFtSmPMS~IBDVEMu?%DEf4q91}+eScJuGu|7OAAK>@sOeUeGx4h z3H57Th8BwY0w91OqewG2kSHaFl454mZ&+3}k**6`lQlwymN+S{vK50_X*{Z^v?galKHBXey-Abe3s@5 znT9XG?J*)+FdU1TXmqvKZY*aN4x)E0QKAMJM54B4#bb&Fg?Gh+Av+jJ=#jtQ?ukZp zGaj*n+#|Z~P)MG2om;7Q`7FcY6s|%`VM=n${i#(|Zm4o<4H(_ndj_2L7=1p!9qK2i z;>zmUni^}_MVjr~WLR2d^=0_CC7{Kk7TOSv_(Gb`+!$w_s;pjNY0LQG#^7cn0`Cr7 zm7I4YS{II^Ju&q2W?v|dlDea&w#n#+Y6lV}2*xBQiFy?c)Ff`%tQ(fi!z{s{Q$DmO z9x?i2Qnk)lvew%Jy@0ZXLMTP^X85w!E7)B z+8YF!-iYjyPiQUoT{om%C3%KqS#7~jFcL8OmRlIaLY0k;7%B-OS)tPlP&dnO24l9b zDF* z4YZ?C3_H}w9g`4?6gKKz<4CsBa+gLJg5p{ib;dDkSRm1px*j&8I<&(KvSP|!ggJN9nm-!p4h8-86|fS9mQ0U!BV)iYp<@9}k~*s4)0Rb|k>$*^ z-@XU~s!KjcmXtD}eq2f@qy=qD6Jy35BmvO0o^v#2}9J%Z zk8`8tws)(!vC7f{5FU?g@D}#<3)6&p0+{GvzZ^NCT9_wj%W3lohPW%ZERjXmVQ739 zRS6|<^DgwjVcn@hpqw*UJvD#7QyFG=S82ULbR-W#e;i{V0_&yuy1T_-fVh>piCbb2 zCioQ=07gXDH}$&8ax^I+Rc?Tru`}_LQpz92k!-%M3&SF(oYQr88e=RPJrMOf!h?jO zeoO=+oz*}+aP%ihm(&wZljS(Q#`58clJlZ?zA+T-@`Z#6)^!z8B;vV06t&`}s8ZKe zB;|QDkmSdYNLu+$hkze59f90k(@h(E#RUOw;SZPQ)+*SbD&?Dt5ptUx{+hK#*$ps+ASw z-eg;V#f2_wlBMyy4z{f;Dg$ah2}EI3c@Tzuw!g=rB9gYv(0JA2%(a^h(?X$IWuUq; zz?N9805LbNz+$6HL*Z(nqovTI*wd|ZW4=hhrxpaLJ{F4FpWM>H^1%NWnw(ksJALol;bHx?n+GzS>Jwd+{1>ewlE6G!m|I{5f#(=Xa7J1LK z!?ZSqP-3Lo=>%QJo?Oj0VF)-`{)hEA(M`^N85U6Nm&NKGmc`kZW516vpyl0Yj zjX4r}$o(tq2R;%({2Ok3srM%-J3m-Ll1pwpxB4UgV6-|hM=FV=uL~SLrERiAJyFmm z$AdJ8@^uT37Z?Lv!hZT}tiZeCY%!_>UDYrb$$n%a*Tv^}?kme3K{z`t&gX&d3l146 z7hF+NuKbi3Q9)X~Ff9(pyp)R%o!r+fcLZS^SfyNi{NujL+!2K5q{ZQHlybpOBIR0; z5+k}WEnY8Q$EKzJUx=2<_j&0h9ei2NlM@QlpR3oCA4vXxl`!zO)1- zkNeuEh_3iM%nP{R{vOc@+i@DxO)#1J%Ho`M~EH-UfbBw z#xw)p!Z@{!K6`0c%z6WOekqoXZ&x;YxM!)#VrA8vhY$e-74z2#hun6={i%=BwY{Ew!xrP_ep+c-c|eBbN&5@~NOUC*~>3qq2@tlC()9w zBInO{Kj>q=;Td&SJ0 zso9eyK0~bo&_0&P)fmHDf%imgL2U(4?K=?+B3XM z7U$98H?nBHolUa?`IKFkM+10EaT0gSzf1Y2#vymm6X*%_ltw>)Ne-2GE}^`coT;qI z3?Y|J56aQKww>r3xL*}~ zLiR|fv~&-=we!vFsmxBakM;Ro)c-Q-eWrv8>{%FBMYMkqW#FFoE0z<+4QL271R4Sj zrO`;)KO6ND?pzX)4YyCnpcN#d-+V+6xt%>R|5IXg8U%Ic9;CV znNQz*qkzsIId7(7d`_&abGEm%20Z0uP~OYswDTd@mA?dAf=_%L!4((sL*IKE@~GjB zLb@tYK%a!ofi6o~L!Kp&C(By&Gh`oolxXn~(ZYA6ez@bF(j4d`Y!C93))q|(x!bvv z9mt_^#F$%(Wh`f#dUGtV~nva1$N%N^9{YbBlWP?x z?N0Vhnnb%@bk87d0BL=$w7YUVbP#E$yFT*Yl2ZA|e~UF3YYd^uSidX$$bU=9p6@KTLvBGcQyArcgaX03c|5NS4};q}5N zf@lRwX>?AnbWsR(aw+R|D==;&eBl=b{lU%%F{sh#uTg`U37=v@Yx?{K-s49jYx!^o?}pz{PjAMn+cnst)p$;9>urW_v(awg zY}W!S!tvrUl+!BT!b`xx3OtC~jLk+!%8%&p6^A$YkrPLLL_c$C7XowY|6EeW_#$l1(V(@w9^IQhrkTbLE;hbR5j%@Jn7UzZtB;?;@1X!*sLGO|a z1R6yjmY~kH;@nP%h4vS#-2DNqcRhz9GUKZlZ^hsX}s3AfM0<>3Nmt z=+-8DrY1@2G?munpcVQ$xLTp+I3L+x$E|1Y48h_zAWe5CIs84ANaqu?9MxJDpLqSN z+S{5a3l&PS(ThlP^q6&HO7#~DRx0z(6h5IEklfdu`qtrX6+Nl)H1lB@cz)b@DinV- znZG^g?`edS1&%J@;gjG-j*m-3+naHo?KFOjr6=_Z3`^z00u`yva51yhxT-|?@)?WB z(Yx84D*ILx|0Z$?ofk15*NvgAw_5AtT{aJLYDJE9@o<-6>YRrDwCM$xbgsj>1<}Py zwm0BFv^x~-m7@Ka;JwMf0XB~QTM{h;+dP(N_i0+f{fPZal6w%#w7Zd}?2i)Zaipv_ zJ{>o3BE*J5cHU=Pbf#h8lP}U%e85c5|433oTbpqD*_TXdiNv+{w?DTvaTZ)_Hu`;`#&_^6B1CV2Sp3EqvM zNyvyzh%hhmEc#})BTD$GhUjK^p~?RF>m)D9=7~r99L-0eqT^|glOo#MjGRixywP4B zNk4Y+CsdR6CM>bi9+xHl=ILsdYEG+?Ww)}!8`1C3&b6@Q&8!@Bscf~r-j2q26bqLY zdeM=WsQjps<+>>HMNP4drebztw9v(K^{qp|)* zIK{zCVZ714WZM$D*WuRXMmgH92c6PDuCgh(Ifx-~#3j#UpJ0!fELF_-21|!WrPP+w z=5l(nO0`021WZXYJIct`+~(#?nmsSCs$#;u4?EfF0317ioz-rK8htT1DH@GYNE36- z?>Vh48r>9+Irao{W=+k?#*>?ZrVW3)4H4{wvDSDqdqV`y&Vi5-NDTO9IP~(!8_R~E z&iJ@C18UtU9of9n%bA|Z?l4T8Dr2?gB+=f?R!c>qX0zdAUcES(hbsxivTlqUX1^CR z1<#DI-Hi|dSHALtrMCv?lde>`ltCgDWP%f^3LS9XX|T`$gF{&DHNkZm14R4}0y|1W;I_%=if@Dqti$E?YueYXXvwQE+? z)rbzoX_jz@fC$3$$z)hJaApt3j2?Ct$3vZae;&bsEYcm^2#LaB1537E9jE%6@##^S zQ#Q_dWJERx;n69a>TzNx>QNL}PW9nEi1O*Zfbf9_tYfj{PI)jF^)a0|yI&r*_+`p+ zziBU5E>O!uE0(WVE*o+sx6vPt;hnMwB|jG-cQYUT|6?WKd^^W~w44|3&|d0(S$^b& z%WyvY!n^JaQE0+LK}6+lTs#o}CVwW+v+-6`T$;lFC^+TgjAbCn1M%62h^Lo3rRv9_ zo_Ox$Pp3Om{NPVzzMWrE<$w9~^-P>2-K3(cSLN(HLwQr@NH?SCoc3v|ogZkAlEY5* zeDVx*rRPYeDLU_0PA?Cq{+e^7t5{@s1YhjJM&QL74G9+yp-i`WSpeHJLXCT=)s1Wk^Rq(Y=?V zPK0*?-@g*aTOmgd1H<19%0OhbO1N=hIuLBBf-8B zD@WjS1{%gD+q4wJ~%d(r17}yGiiFc+(yk?gai^g=c_UzmB}1-w7<- zC(9`XejI)gV@P-VJGlKL*l*%9kHK9DG>jv-Pa@3t3htu+E%A2Deeal?gV}Y_YlJD+c}AwWzP8Kcfd2!r-40x z3q+W4=K;v`U5PIu%wu5^nDvB|Gvhzu4uIx3u=OBx72!7E3EX`MGoJrlSyv0tufhS~ zBPu)ud|8DV-%{aU0?QA{yo?vA@WsFfa6b%w7@tz5Jzf`$>_l_0_$z{Q7AH)z<8Z69uy=i|WKBN%r` z*MJe+(+KYbK8Sk~;r&4NH93wjBm0g_LFSBasxTw_esDS?`*E=En?QuekoTkOzoG;@ zcpC&i7!0#i?%;oegqqW&!%Yr82huP0i>aAJ?5b<}!Xf8VZL6wQRM%8#I9-Y2+bA}! zs=_J9^2@65QC_?|^x@Ok#;U4*!>YRG>f(Z{@aZBW+=UOfqJW6CstTW=zQV#;W7uac z4+s5b)QWc7%kg=;D|}YCdh?1Zd=xJd?8ZAXbuc1pE)8TF3dG8!B08yy(kH99o9cXWJoVsvtJdUR&=ma*8_ zF7S9_jCL>EU9 None: delete_api: DeleteApi = client.delete_api() delete_api.delete(start=start_time, stop=stop_time, predicate=predicate, bucket=bucket) + + #2025/08/18 从文件导入scada数据，xkl +def import_data_from_file(file_path: str, bucket: str = "SCADA_data") -> None: + """ + 从指定的CSV文件导入数据到InfluxDB的指定bucket中。 + :param file_path: CSV文件的路径 + :param bucket: 数据存储的 bucket 名称，默认值为 "SCADA_data" + :return: + """ + client = get_new_client() + if not client.ping(): + print("{} -- Failed to connect to InfluxDB.".format(datetime.now().strftime('%Y-%m-%d %H:%M:%S'))) + + #清空指定bucket的数据 + # delete_api = DeleteApi(client) + # start = "1970-01-01T00:00:00Z" + # stop = "2100-01-01T00:00:00Z" + # delete_api.delete(start, stop, '', bucket="SCADA_data", org="TJWATERORG") + + df = pd.read_csv(file_path) + write_api = client.write_api(write_options=SYNCHRONOUS) + points_to_write = [] + for _, row in df.iterrows(): + scada_id = row['ScadaId'] + value = row['Value'] + time_str = row['Time'] + date_str = str(time_str)[:10] # 取前10位作为日期 + try: + raw_value = float(value) + except (ValueError, TypeError): + raw_value = 0.0 + point = Point("SCADA") \ + .tag("date", date_str) \ + .tag("description", None) \ + .tag("device_ID", scada_id) \ + .field("monitored_value", raw_value) \ + .field("datacleaning_value", 0.0) \ + .field("simulation_value", 0.0) \ + .time(time_str, write_precision='s') + points_to_write.append(point) + # 批量写入数据 + batch_size = 500 + for i in range(0, len(points_to_write), batch_size): + batch = points_to_write[i:i+batch_size] + write_api.write(bucket=bucket, record=batch) + print(f"Data imported from {file_path} to bucket {bucket} successfully.") + print(f"Total points written: {len(points_to_write)}") + write_api.close() + client.close() + +#2025/08/28 从多列格式文件导入SCADA数据，xkl +def import_multicolumn_data_from_file(file_path: str, raw: bool = True, bucket: str = "SCADA_data") -> None: + """ + 从指定的多列格式CSV文件导入数据到InfluxDB的指定bucket中。 + :param file_path: CSV文件的路径 + :param bucket: 数据存储的 bucket 名称，默认值为 "SCADA_data" + :return: + """ + client = get_new_client() + write_api = client.write_api(write_options=SYNCHRONOUS) + points_to_write = [] + if not client.ping(): + print("{} -- Failed to connect to InfluxDB.".format(datetime.now().strftime('%Y/%m/%d %H:%M'))) + def convert_to_iso(timestr): + # 假设原格式为 '2025/8/3 0:00' + dt = datetime.strptime(timestr, '%Y-%m-%d %H:%M:%S') + return dt.isoformat() + + with open(file_path, encoding='utf-8') as f: + reader = csv.reader(f) + header = next(reader) + device_ids = header[1:] # 第一列是time，后面是device_ID + if raw: + for row in reader: + time_str = row[0] + iso_time = convert_to_iso(time_str) + for idx, value in enumerate(row[1:]): + try: + raw_value = float(value) + except (ValueError, TypeError): + raw_value = 0.0 + scada_id = device_ids[idx] + # 如果是原始数据，直接使用Value列 + point = Point("SCADA") \ + .tag("date", iso_time.split('T')[0]) \ + .tag("description", None) \ + .tag("device_ID", scada_id) \ + .field("monitored_value", raw_value) \ + .field("datacleaning_value", 0.0) \ + .field("simulation_value", 0.0) \ + .time(iso_time, WritePrecision.S) + points_to_write.append(point) + else: + for row in reader: + time_str = row[0] + iso_time = convert_to_iso(time_str) + # 如果不是原始数据，直接使用datacleaning_value列 + for idx, value in enumerate(row[1:]): + scada_id = device_ids[idx] + try: + datacleaning_value = float(value) + except (ValueError, TypeError): + datacleaning_value = 0.0 + # 如果是清洗数据，直接使用datacleaning_value列 + point = Point("SCADA") \ + .tag("date", iso_time.split('T')[0]) \ + .tag("description", "None") \ + .tag("device_ID", scada_id) \ + .field("monitored_value", 0.0) \ + .field("datacleaning_value",datacleaning_value) \ + .field("simulation_value", 0.0) \ + .time(iso_time, WritePrecision.S) + points_to_write.append(point) + # 批量写入数据 + batch_size = 1000 + for i in range(0, len(points_to_write), batch_size): + batch = points_to_write[i:i+batch_size] + write_api.write(bucket=bucket, record=batch) + print(f"Data imported from {file_path} to bucket {bucket} successfully.") + print(f"Total points written: {len(points_to_write)}") + write_api.close() + client.close() + # 示例调用 @@ -3652,5 +3775,25 @@ if __name__ == "__main__": # end_time='2024-03-26T23:59:00+08:00') # print(result) + #示例：import_data_from_file + #import_data_from_file(file_path='data/Flow_Timedata.csv', bucket='SCADA_data') + # # 示例：query_all_records_by_type_date + #result = query_all__records_by_type__date(type='node', query_date='2025-08-04') + + #示例：query_all_records_by_date_hour + #result = query_all_records_by_date_hour(query_date='2025-08-04', query_hour=1) + + #示例：import_multicolumn_data_from_file + #import_multicolumn_data_from_file(file_path='data/selected_Flow_Timedata2025_new_format_cleaned.csv', raw=False, bucket='SCADA_data') + + # client = InfluxDBClient(url="http://localhost:8086", token=token, org=org_name) + # delete_api = client.delete_api() + + # start = "2025-08-02T00:00:00Z" # 要删除的起始时间 + # stop = "2025-08-11T00:00:00Z" # 结束时间（可设为未来） + # predicate = '_measurement="SCADA"' # 指定 measurement + + # delete_api.delete(start, stop, predicate, bucket="SCADA_data", org=org_name) + # client.close() diff --git a/online_Analysis.py b/online_Analysis.py index 52ea9b1..69d76f3 100644 --- a/online_Analysis.py +++ b/online_Analysis.py @@ -19,6 +19,7 @@ from sqlalchemy import create_engine import ast import sensitivity import project_info +import api_ex.kmeans_sensor ############################################################ # burst analysis 01 @@ -1064,6 +1065,85 @@ def pressure_sensor_placement_sensitivity(name: str, scheme_name: str, sensor_nu except Exception as e: print(f"存储方案信息时出错：{e}") +#2025/08/21 +# 基于kmeans聚类法进行压力监测点优化布置 +def pressure_sensor_placement_kmeans(name: str, scheme_name: str, sensor_number: int, + min_diameter: int, username: str) -> None: + """ + 基于聚类法进行压力监测点优化布置 + :param name: 数据库名称（注意，此处数据库名称也是inp文件名称，inp文件与pg库名要一样） + :param scheme_name: 监测优化布置方案名称 + :param sensor_number: 传感器数目 + :param min_diameter: 最小管径 + :param username: 用户名 + :return: + """ + sensor_location = api_ex.kmeans_sensor.kmeans_sensor_placement(name=name, sensor_num=sensor_number, min_diameter=min_diameter) + try: + conn_string = f"dbname={name} host=127.0.0.1" + with psycopg.connect(conn_string) as conn: + with conn.cursor() as cur: + sql = """ + INSERT INTO sensor_placement (scheme_name, sensor_number, min_diameter, username, sensor_location) + VALUES (%s, %s, %s, %s, %s) + """ + + cur.execute(sql, (scheme_name, sensor_number, min_diameter, username, sensor_location)) + conn.commit() + print("方案信息存储成功！") + except Exception as e: + print(f"存储方案信息时出错：{e}") + +############################################################ +# 流量监测数据清洗 ***卡尔曼滤波法*** +############################################################ +#2025/08/21 hxyan + +def flow_data_clean(input_csv_file: str) -> str: + """ + 读取 input_csv_path 中的每列时间序列，使用一维 Kalman 滤波平滑并用预测值替换基于 3σ 检测出的异常点。 + 保存输出为：_cleaned.xlsx（与输入同目录），并返回输出文件的绝对路径。如有同名文件存在，则覆盖。 + :param: input_csv_file: 输入的 CSV 文件明或路径 + :return: 输出文件的绝对路径 + """ + + + # 提供的 input_csv_path 绝对路径，以下为 默认脚本目录下同名 CSV 文件，构建绝对路径，可根据情况修改 + script_dir = os.path.dirname(os.path.abspath(__file__)) + input_csv_path= os.path.join(script_dir, input_csv_file) + + # 检查文件是否存在 + if not os.path.exists(input_csv_path): + raise FileNotFoundError(f"指定的文件不存在: {input_csv_path}") + # 调用 Fdataclean.clean_flow_data_kf 函数进行数据清洗 + out_xlsx_path = api_ex.Fdataclean.clean_flow_data_kf(input_csv_path) + print("清洗后的数据已保存到:", out_xlsx_path ) + + +############################################################ +# 压力监测数据清洗 ***kmean++法*** +############################################################ +#2025/08/21 hxyan + +def pressure_data_clean(input_csv_file: str) -> str: + """ + 读取 input_csv_path 中的每列时间序列，使用Kmean++清洗数据。 + 保存输出为：_cleaned.xlsx（与输入同目录），并返回输出文件的绝对路径。如有同名文件存在，则覆盖。 + 原始数据在 sheet 'raw_pressure_data'，处理后数据在 sheet 'cleaned_pressusre_data'。 + :param input_csv_path: 输入的 CSV 文件路径 + :return: 输出文件的绝对路径 + """ + + # 提供的 input_csv_path 绝对路径，以下为 默认脚本目录下同名 CSV 文件，构建绝对路径，可根据情况修改 + script_dir = os.path.dirname(os.path.abspath(__file__)) + input_csv_path= os.path.join(script_dir, input_csv_file) + + # 检查文件是否存在 + if not os.path.exists(input_csv_path): + raise FileNotFoundError(f"指定的文件不存在: {input_csv_path}") + # 调用 Fdataclean.clean_flow_data_kf 函数进行数据清洗 + out_xlsx_path = api_ex.Pdataclean.clean_pressure_data_km(input_csv_path) + print("清洗后的数据已保存到:", out_xlsx_path ) if __name__ == '__main__': # contaminant_simulation('bb_model','2024-06-24T00:00:00Z','ZBBDTZDP009034',30,1800) @@ -1120,3 +1200,5 @@ if __name__ == '__main__': # 示例：pressure_sensor_placement_sensitivity pressure_sensor_placement_sensitivity(name=project_info.name, scheme_name='20250517', sensor_number=10, min_diameter=300, username='admin') + # 示例：pressure_sensor_placement_kmeans + pressure_sensor_placement_kmeans(name=project_info.name, scheme_name='sensor_1027', sensor_number=35, min_diameter=300, username='admin') \ No newline at end of file