{"cells":[{"metadata":{"papermill":{"duration":0.026459,"end_time":"2020-09-16T09:09:10.277703","exception":false,"start_time":"2020-09-16T09:09:10.251244","status":"completed"},"tags":[]},"cell_type":"markdown","source":"## Background\nThis work is a fork from great notebook:\nhttps://www.kaggle.com/hfutybx/osic-feature-extract-from-ct\n\nI am just refactoring the code for me to make it easier to extract features from both train and test dicom. The output is just extracted features (volume, mean, skew and kurthosis) without merging tabular data.\nI also optionally remove / remark visualization with the purpose this notebook can be integrated later to model training and inference."},{"metadata":{"trusted":true},"cell_type":"code","source":"!dpkg -i ../input/python3gdcm/build_1-1_amd64.deb\n!apt-get install -f","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!cp /usr/local/lib/gdcm.py /opt/conda/lib/python3.7/site-packages/.\n!cp /usr/local/lib/gdcmswig.py /opt/conda/lib/python3.7/site-packages/.\n!cp /usr/local/lib/_gdcmswig.so /opt/conda/lib/python3.7/site-packages/.\n!cp /usr/local/lib/libgdcm* /opt/conda/lib/python3.7/site-packages/.\n!ldconfig","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","execution":{"iopub.execute_input":"2020-09-16T09:09:19.116200Z","iopub.status.busy":"2020-09-16T09:09:19.115191Z","iopub.status.idle":"2020-09-16T09:09:21.470015Z","shell.execute_reply":"2020-09-16T09:09:21.468835Z"},"papermill":{"duration":2.398479,"end_time":"2020-09-16T09:09:21.470137","exception":false,"start_time":"2020-09-16T09:09:19.071658","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"import os\nimport cv2\nimport sys\nimport random\nimport pickle\nimport glob\nimport gc\nfrom tqdm.notebook import tqdm\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom tqdm.auto import tqdm\nimport gdcm\nimport pydicom\n\nimport scipy.ndimage as ndimage\nfrom scipy.ndimage import zoom\nfrom scipy.stats import kurtosis\nfrom scipy.stats import skew\n\nimport torch\nimport torch.nn as nn\nfrom torch.utils.data import DataLoader\nfrom torch.utils.data.dataset import Dataset\n\nimport warnings\nwarnings.filterwarnings(\"ignore\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sys.path.append('../input/efficientnet-pytorch/EfficientNet-PyTorch-master/')\nsys.path.append('../input/pretrainedmodels/pretrainedmodels-0.7.4/')\nsys.path.append('../input/segmentation-models-pytorch/')","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:21.543435Z","iopub.status.busy":"2020-09-16T09:09:21.542717Z","iopub.status.idle":"2020-09-16T09:09:23.221752Z","shell.execute_reply":"2020-09-16T09:09:23.221125Z"},"papermill":{"duration":1.718,"end_time":"2020-09-16T09:09:23.221875","exception":false,"start_time":"2020-09-16T09:09:21.503875","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"import segmentation_models_pytorch as smp","execution_count":null,"outputs":[]},{"metadata":{"papermill":{"duration":0.037682,"end_time":"2020-09-16T09:09:23.297837","exception":false,"start_time":"2020-09-16T09:09:23.260155","status":"completed"},"tags":[]},"cell_type":"markdown","source":"## CT Image Preprocessing\n\nThe size of the images are not uniform, need to unify them to 512*512"},{"metadata":{"trusted":true},"cell_type":"code","source":"#INPUT = './input' #local\nINPUT = '/kaggle/input/osic-pulmonary-fibrosis-progression' #kaggle\n#WORKING = ./output #local\nWORKING = '/kaggle/working' #kaggle","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_dicom_df(path: str)->pd.DataFrame:\n    Patient_ids = os.listdir(path)\n    n_dicom_dict = {\"Patient\":[],\"n_dicom\":[],\"list_dicom\":[], \"height\": [], \"width\": []}\n\n    for Patient_id in Patient_ids:\n        dicom_id_path = glob.glob(path + Patient_id + \"/*\")\n        n_dicom_dict[\"n_dicom\"].append(len(dicom_id_path))\n        n_dicom_dict[\"Patient\"].append(Patient_id)\n        list_dicom_id = sorted([int(i.split(\"/\")[-1][:-4]) for i in dicom_id_path])\n        n_dicom_dict[\"list_dicom\"].append(list_dicom_id)\n        \n        rows = []\n        cols = []\n        for patient_dicom_id_path in dicom_id_path:\n            dicom = pydicom.dcmread(patient_dicom_id_path)\n            #print(dicom.Rows, dicom.Columns)\n            rows = dicom.Rows\n            cols = dicom.Columns\n            break\n        n_dicom_dict[\"height\"].append(rows)\n        n_dicom_dict[\"width\"].append(cols)\n\n    dicom_df = pd.DataFrame(n_dicom_dict)\n    return dicom_df","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:23.382140Z","iopub.status.busy":"2020-09-16T09:09:23.381373Z","iopub.status.idle":"2020-09-16T09:09:24.235921Z","shell.execute_reply":"2020-09-16T09:09:24.235392Z"},"papermill":{"duration":0.901675,"end_time":"2020-09-16T09:09:24.236032","exception":false,"start_time":"2020-09-16T09:09:23.334357","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"train_dicom_path = f'{INPUT}/train/'\ntest_dicom_path = f'{INPUT}/test/'","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_dicom_df =  get_dicom_df(train_dicom_path)\n#train_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_dicom_df =  get_dicom_df(test_dicom_path)\n#test_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:24.312029Z","iopub.status.busy":"2020-09-16T09:09:24.311355Z","iopub.status.idle":"2020-09-16T09:09:24.528806Z","shell.execute_reply":"2020-09-16T09:09:24.529391Z"},"papermill":{"duration":0.258194,"end_time":"2020-09-16T09:09:24.529534","exception":false,"start_time":"2020-09-16T09:09:24.271340","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"if 0:    \n    print(f\"num dicom: {len(train_dicom_df)}\\n\\\n    min dicom number is {min(train_dicom_df['n_dicom'])}\\n\\\n    max dicom number is {max(train_dicom_df['n_dicom'])}\")\n\n    plt.hist(train_dicom_df['n_dicom'], bins=20)\n    plt.title('Number of dicom per patient');","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"if 0:\n    print(f\"num dicom: {len(test_dicom_df)}\\n\\\n    min dicom number is {min(test_dicom_df['n_dicom'])}\\n\\\n    max dicom number is {max(test_dicom_df['n_dicom'])}\")\n\n    plt.hist(test_dicom_df['n_dicom'], bins=20)\n    plt.title('Number of dicom per patient');","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def mark_reshape(df: pd.DataFrame, size: int)->pd.DataFrame:\n    reshape_df = df.loc[(df.height!=size) | (df.width!=size),:]\n    reshape_df = reshape_df.reset_index(drop=True)\n    \n    crop_id = list(df[df.height!=df.width][\"Patient\"])\n    reshape_df['resize_type'] = 'resize'\n    reshape_df.loc[reshape_df.Patient.isin(crop_id),'resize_type'] = 'crop'\n    \n    return reshape_df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reshape_dicom_pd = mark_reshape(train_dicom_df, 512)\n#reshape_dicom_pd.head(10)","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:28.593120Z","iopub.status.busy":"2020-09-16T09:09:28.592487Z","iopub.status.idle":"2020-09-16T09:09:28.598679Z","shell.execute_reply":"2020-09-16T09:09:28.599122Z"},"papermill":{"duration":0.098745,"end_time":"2020-09-16T09:09:28.599256","exception":false,"start_time":"2020-09-16T09:09:28.500511","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"train_dicom_df['resize_type'] = 'no'\nfor idx,i in enumerate(reshape_dicom_pd['Patient']):\n    train_dicom_df.loc[train_dicom_df.Patient==i,'resize_type'] = reshape_dicom_pd.loc[idx,'resize_type']\n#train_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reshape_dicom_pd = mark_reshape(test_dicom_df, 512)\n#reshape_dicom_pd.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_dicom_df['resize_type'] = 'no'\nfor idx,i in enumerate(reshape_dicom_pd['Patient']):\n    test_dicom_df.loc[test_dicom_df.Patient==i,'resize_type'] = reshape_dicom_pd.loc[idx,'resize_type']\n#test_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"del reshape_dicom_pd","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:28.915633Z","iopub.status.busy":"2020-09-16T09:09:28.914725Z","iopub.status.idle":"2020-09-16T09:09:30.258773Z","shell.execute_reply":"2020-09-16T09:09:30.258246Z"},"papermill":{"duration":1.42923,"end_time":"2020-09-16T09:09:30.258880","exception":false,"start_time":"2020-09-16T09:09:28.829650","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"#merge/not merge with tabular data\nif 0:\n    train_df = pd.read_csv(f'{INPUT}/train.csv')\n    temp_df = pd.DataFrame(columns=train_df.columns)\n    for i in range(len(train_dicom_df)):\n        patient_df = train_df[train_df.Patient==train_dicom_df.iloc[i].Patient]\n        zeroweek = patient_df['Weeks'].min()\n        #if sum(patient_pd.Weeks==zeroweek)>1:\n        #    print(pd.unique(patient_pd.Patient))\n        temp_df = temp_df.append(patient_df[patient_df.Weeks==zeroweek].iloc[0])\n    train_dicom_df = pd.merge(train_dicom_df, temp_df, on=['Patient'])\n    train_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:30.486922Z","iopub.status.busy":"2020-09-16T09:09:30.486247Z","iopub.status.idle":"2020-09-16T09:09:30.489813Z","shell.execute_reply":"2020-09-16T09:09:30.489308Z"},"papermill":{"duration":0.059238,"end_time":"2020-09-16T09:09:30.489930","exception":false,"start_time":"2020-09-16T09:09:30.430692","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"def load_scan(path,resize_type='no'):\n    \"\"\"\n    Loads scans from a folder and into a list.\n    \n    Parameters: path (Folder path)\n    \n    Returns: slices (List of slices)\n    \"\"\"\n    slices = [pydicom.read_file(path + '/' + s) for s in os.listdir(path)]\n    slices.sort(key = lambda x: int(x.InstanceNumber))\n    \n    try:\n        slice_thickness = abs(slices[-1].ImagePositionPatient[2] - slices[0].ImagePositionPatient[2])/(len(slices))\n    except:\n        try:\n            slice_thickness = abs(slices[-1].SliceLocation - slices[0].SliceLocation)/(len(slices))\n        except:\n            slice_thickness = slices[0].SliceThickness\n        \n    for s in slices:\n        s.SliceThickness = slice_thickness\n        if resize_type == 'resize':\n            s.PixelSpacing = s.PixelSpacing*(s.Rows/512)  \n    return slices","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:30.593251Z","iopub.status.busy":"2020-09-16T09:09:30.592032Z","iopub.status.idle":"2020-09-16T09:09:30.594739Z","shell.execute_reply":"2020-09-16T09:09:30.595240Z"},"papermill":{"duration":0.059417,"end_time":"2020-09-16T09:09:30.595352","exception":false,"start_time":"2020-09-16T09:09:30.535935","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"def transform_to_hu(slices):\n    \"\"\"\n    transform dicom.pixel_array to Hounsfield.\n    Parameters: list dicoms\n    Returns:numpy Hounsfield\n    \"\"\"\n    \n    images = np.stack([file.pixel_array for file in slices])\n    images = images.astype(np.int16)\n\n    # convert ouside pixel-values to air:\n    # I'm using <= -1000 to be sure that other defaults are captured as well\n    #images[images <= -1000] = 0\n    \n    # convert to HU\n    for n in range(len(slices)):\n        \n        intercept = slices[n].RescaleIntercept\n        slope = slices[n].RescaleSlope\n        \n        if slope != 1:\n            images[n] = slope * images[n].astype(np.float64)\n            images[n] = images[n].astype(np.int16)\n            \n        images[n] += np.int16(intercept)\n    \n    return np.array(images, dtype=np.int16)","execution_count":null,"outputs":[]},{"metadata":{"papermill":{"duration":0.048416,"end_time":"2020-09-16T09:09:30.689473","exception":false,"start_time":"2020-09-16T09:09:30.641057","status":"completed"},"tags":[]},"cell_type":"markdown","source":"![image.png](attachment:image.png)","attachments":{"image.png":{"image/png":"iVBORw0KGgoAAAANSUhEUgAABFAAAAH6CAYAAAAgOLgTAAAgAElEQVR4Aeydi5IEKY5s+/9/etbYaZ/yUgkQEJHPg9lcCXA9OBmdROX07P3nPwwIQAACEIAABCAAAQhAAAIQgAAEIACBIYF/hrtsQgACEIAABCAAAQhAAAIQgAAEIAABCPyHH1B4CCAAAQhAAAIQgAAEIAABCEAAAhCAwIQAP6BMALENAQhAAAIQgAAEIAABCEAAAhCAAAT4AYVnAAIQgAAEIAABCEAAAhCAAAQgAAEITAjwA8oEENsQgAAEIAABCEAAAhCAAAQgAAEIQIAfUHgGIAABCEAAAhCAAAQgAAEIQAACEIDAhAA/oEwAsQ0BCEAAAhCAAAQgAAEIQAACEIAABPgBhWcAAhCAAAQgAAEIQAACEIAABCAAAQhMCPADygQQ2xCAAAQgAAEIQAACEIAABCAAAQhAYOkHlH/++ec/+k9DJ1824uytR53PRzFtrze0p3jsz2cFC1jwDPAM8AzwDPAM8AzwDNz/DLT31BHnuK/32lGM9qSV1Xq02s+sa7XfW2v72pNWtreufbdNq+Fx7ms/Wmlk477mbX91xJyVHDFmtWZFX+mjkgcNBCBwD4Hyt42+MLD3X74whjHPAM8AzwDPAM8AzwDPAM8AzwDPAM8AzwDPwH3PwM5PLPyAYv9WDQ/nfQ8nbGHLM8AzwDPAM8AzwDPAM8AzwDPAM8AzwDPwSs/A6o8oWz+geBEdPlvT3o5t+RTnflvT0L6vac+tdL6GDwEIQAACEIDA+xPQHT97F3j/k3KCGYHsGdCanpOWQ2uzfL6veFnfcz/mHunjnsdme77fakqTWe+pp+2tZ/lU2/dUw9c+yY984vyTzspZ+EHjm58BfZdV7c+vEZMIhzqRsg0BCEAAAhCAAAQeQoD3k4dgpggEIAABCEDgYwicvDvwA8rHPAYcBAIQgAAEIPB9BE5egr6PFieGAAQgAAEIQODk3YEfUHh+IAABCEAAAhB4WwInL0Fve2gahwAEIAABCEBgm8DJuwM/oGxjJxACEIAABCAAgWcTOHkJenbv1IcABCAAAQhA4PEETt4d+AHl8Z8XFSEAAQhAAAIQuIjAyUvQRS2QBgIQgAAEIACBNyJw8u7ADyhv9EHTKgQgAAEIQAACvwmcvAT9zsQMAhCAAAQgAIFvIHDy7sAPKN/whHBGCEAAAhCAwIcSOHkJ+lAkHAsCEIAABCAAgQGBk3cHfkAZgGULAhCAAAQgAIHXJnDyEvTaJ6M7CEAAAhCAAATuIHDy7sAPKHd8IuSEAAQgAAEIQOAhBE5egh7SIEUgAAEIQAACEHgpAifvDvyA8lIfJc1AAAIQgAAEILBC4OQlaKUOWghAAAIQgAAEPoPAybsDP6BsPgOCvhlO2AGBV2LfemFAAAIQgMDzCOhOuPP7+CR3FputPY8glSEAAQhAAALfRaDdw/rP6snLf/2pwF2Xvudf8VcP3NOv1HRtLx/r1xNw7u5fX6mW0XsY+bVsqCAAAQhAYIeAf//uxM9iPP+JrzonOTxW+bDvR8A/x8yPJ6poYsyV86x+tnZSM8s3WjuppdhR/rbXxkwjnXLKVuJco7hnW++pd7Zn90h9CFxBwJ/11Xxv/wPKFf9wO8BVfxU4+j0Cs89lL+tZ1Kwn3z+rRDQEIAABCPQI3P1d6/lPfPV/ksNjlQ/7PgT885v5fqqe1jV3+r36cf2kh5hrNj+ppdhZjeq+8rmtxkrnsc/y1Uu0z+qHuhC4k4A/56t1XuYHFG/cD1TxPXbHr9SQZif/o2LeocddFjpbz+7mvSKu15OvX1GHHBCAAAQg8JfAo79rvd7I/9tpvjLK4Xt5NKvvQsA/y+a34Tbuz/a0/6jzZ/3Ftat6iXmz+VW1lCerka01fVt3qxzRZvEeq/0Y94y5eunZZ/RETQjcScCf9dU6L/kDig7hBxv50u/aUe64t1vjEXHq9RG1Hl1DZ+vZR/cT6/X60nrUM4cABCAAgWsI6Hu22UcNr5n5q31kOeLaak70r0Mgfpaae4daW7Ee/wh/1tvVPYzqXV1L+UY1tSftilWs7ErsI7Tqq2dXe1Ce1Tj0EHgUAT2jza6OcsRJkdWmpPeaI1/6HTvKG/d28j8qxnt9VM1H1vHzZf4je8lqZT1pLdOzBgEIQAAC1xDQd22zjxpeM/NX+8hy+NpqPvSvRcA/S/fVpa+t+Ip/pB31d3Ufj6yl3kc1297J8Nwnee6I9d4yf7WmcqzGoYfAowjoGW12dZQjToqsNiW915z5ilm1s7y+v5r7kfp36fOEiZ/R/ZOcV8V6P9G/qgZ5IAABCEDgLwH/zv27e8+K18z81apZDl9bzYf+dQj45xh97zLutblGtuf70j3C9nq5o59H1hK7Uc3TMyq3ar2SVW89u9JrzLESixYCjyLgz+lqzZ9v50nkSZFJ6u6216z43USdDeVs2/JHtpPm6cux56c3dGMDr3rW2JfPb8RBaghAAAJfT+AZ37deM/NXP5Qsh6+t5kP/OgT8c4y+d7my17TPGrFPn1/dk+eO/tW1lC/WiXPpdqxy7cQ+Ikb9Rbta2+NXY9FD4FEETp7T8jfwSZFdEKqpeM17Vrqq9TyK8bXoS/NKNvbY5ozHE8g+B609vhsqQgACEPgeAvqufeT95zUzf5V+lsPXVvOhfx0C/jlG/3W6rHcSz+Dzepaa0nNHv5ZhXRXrxPl6xp8I5fpZeU3vpE/Fyr7mCekKAr//5YlVHuW/tvUPQrOPGqqpepr3rHRVqzyu11rPuvbZ/jv0+GxGj6rf+yzaOgMCEIAABO4j4N+/91X5ndlrZv5v9XyW5fC1eQYUr0rAP8fov2rPo77iGXw+itvZ89zR38lXiYl14rySo6fxXD3NO6/7+eS/83no/bMJ6BltdnWUI06KrDYlvWrGudYzK+3MeqxrfT3zXftMP+tNa8/s61tri31mv5UJ54YABCDwCAL+vfuIeq2G18z81T6yHL62mg/96xDwzzH6r9NlvZN4Bp/Xs9SUnjv6tQzrqlgnztcz/kR4rp/Vz/D8bO5/xuk4xScSOHlO3+oHlPbh+WEzv/oBe6zH+Hrmu3bVvypflsfXVvuKes/V/OpY0cacqhnX75qrXrS79WIen+/mJA4CEIAABOYEnvF96zUzf971b0WWw9d+q/dmV+fb6+KzohrT2XDu0Z/FvuJ+PIPPr+7Xc0f/6lrKF+vEuXSr1vOsxr663s8W/Vfvnf6+l4A/q6sU5t/8/2Y8KbLalPSqqXmzWhtZ12e+x/q+r/d818/8Xo5sfZZL+1lsZU3xPXuaoxLfNHH04qTr7cd16Wc2xo3ms1y+f1Uez4kPAQhAAAJzAv79O1dfo/Camb9aJcvha6v5pPccPV/ake3FxnXPEfd6c4+R39PGdembjXs+d13F99joK17rmrtte21IU7XKsapX3MjOco5ifW+Ux3Uzv5fH43qatn7XGNU8qet5e727Rr5rszXta8+t9kZ2pB/tKadrMl+6ZrXva/gQeBYBPY/Nro5yxEmR1aakV03Nm9XayLo+8z3W9309810787P4ylovbyV2pslyz2LifpajrUVdb+7xPU1b1xhpfE/6nnXtit/LF9dHOaOWOQQgAAEIXEfAv3+vyzrO5DUzfxz9dzfL4Wt/I8YrHlv1Rxl3cuzEqIfV2IpeuWe2kss1MZ/vrfrKVY2Tfmar+ZpuNka5ZrFtfxSvPeXRPLPSXG2zWr52Uk95ejm0P7Mev6L1uObPYn3/JNbzyI/5mEPg0QT0LDa7OsoRJ0VWm5JeNTWX1XrPSpdZj4n7vtfzY0w2j7EVjces6j2256/kdG2Wz/ejn+l9TXpf6/nSNtvTaN210Zem2d5wTfR7Mb4eY3zuOnwIQAACELiWwDO+b71m5q+eMMvhayv5PE5+jNd6tFEX51Ef51Hf5lET51mM1qL2dK68PRvzR13cb/M4Mk11bTVX1GfzrLZ0oz1pos1itBa1cS5dtNLF9dFcMVfbUc22tzs87yiH60Z+yzHa972r6imP5971lQsLgWcR8Gd3tYfyN8FJkdWmpFdNzWW1PrLSRusxoz3XuR9j4ty18qOmzbXXs1mMr/Xi2npl9OKz2Eyb6dpapvW1ikZ6r6G1nnWt+5ne993PtFpzXeZLl9lMzxoEIAABCFxDwL93r8k4z+I1M3+e4bciy+Frv9X9mce4n0X4vvuZ1tdcG33XuR91Pndd9F13hR/z+zzm9z33qzrFRL3m2h9ZaXt2JzbGZLmjxueZXmuui7400UZdm0dNNs/irljLavnaSQ3lqeSQ9io7qzmrE+Nn+tF+zMUcAs8g4M/oav3aX9vhy2y1yK5eB8vitdezWUxbc33U+F7PjzE+78W09Wys6j3HSWzL04v3GvIzrfYym+l31jz3LN618nsx2o+2p2/rs3ESO8vNPgQgAAEI9An4929fde2O18z8WC1qZvszfYzXPMZprn232ovWNZkf9T7P9G3NNdHvxWg96rO5tLNaLbY3Znk9zrW+nvmujX6mj2sxxudR63PXue+a5vue+1GnuWuiL020UefzqNXcNZkv3dU2q+Vru/VWc7h+5Hs/VZ3HyB/Ftr3RGMWO4tiDwDMJ+HO72sf4nwjLdlLE0iy5qpkFaW9kY5xr416b+37mZzG+lsVozXXytZdZaXo2i9FaL0br0mVWGreZrq31Rk+v9RindbcVzUjf4n0/+jH/jt5zxPw+dx0+BCAAAQhcS+AZ37de8xF+hVivj1HslTEtV2/06oxilGsU24sfxShvtFlM1GjuWq31rGuj34vx9Rjjc9e57xr3XeO+a+T7vvvaz6zr3M+0bW00ejFaH8We7Cl/z+7m9nyVHK7v+TFPT9fWZ2MUO4sfxc7qsg+BZxHw53a1h/k/Uf9mPCmy2pT0s5q+n/nKI+sarbn1/Z7v+uj3YrS+oo/aOFfOzEZtnGcxba03rtIrT1ZHe70+fD/zZzljzBV6zxHz+9x1+BCAAAQgcC2BZ3zfes1H+BVivT5GsVfGtFy90aszilGundjVmJG+16Ni1GfPSpfZXoyvZ3Fac5372o/WNe5HXZv3RqbVWhajvcxmeq1lel+T7mrrNTJ/t57nquRwfeZnOTKdr2Uxvuba6Lsu+lHr86hlDoFXIXDynPa/IcPpToqEVOXprKbvZ34s5Jq41+a+3/OzOK31Ytp6Nlb1nuOOWM/vfq+Wa9zv6bXu2qqv2J7N8vS0bT0bq3rPcRLrefAhAAEIQGCNgH//rkXuq71m5q9mznL42iyfa6M/io1an/fiXBP9K2OUK9bwuTTRuib6UdvmUZPNY5w0cT3Opcts1GbzLE5rV+hbDuWL9or8MafPs/y+5trou+5KP9aJ85NaylXJIW3PZjl6Wq1nMVqTpmely2wvpq0zIPCqBPy5Xe2x/GSfFFltSnrV1Dyz0vSsYnxfa9G6JvOjPpuvxGVarWW5fU26zLou+pn+ZC3mb/NZvixmtrabsxeX1etp2/psnMTOcrMPAQhAAAJ9Av7921ddu+M1M3+1WpbD12b5XBv9UWzU+rwX55roXxmjXLGGz6WJ1jXRj1rNo643l75qe3naemWsxr+73pmsnsVjd/1RzbZ3MpS7kkPans1y9LRaz2K0Jk3PSpfZXkxbZ0DgVQn4c7vaY/nJPimy2pT0qql5ZqXpWcX4vtaidU3mR/3uPMsd12a5o97no1jXXeFntUZ5M31lbZSz7Z2MWe5K/lGOk96IhQAEIACBMQH//h0rr9v1mpm/WinL4WuzfK6N/ig2an3ei3NN9K+MUa5Yw+fSROua6Eet5lE3miumYk/zrMaP9Kt72flGOa7Qe47VWh67649qtr2TodyVHNL2bJajp9V6FqM1aXpWusz2Yto6AwKvSsCf29Uey0/2SZHVpqRXTc0zK03VZjm0Nssh3a6d5ff9WQ3XRn8UG7Vxrti4ns2ljTbTai1qq3PF92w1j+t6ubJ1j8v8LEZrmZ41CEAAAhC4hoC+a5t91PCamb/aR5bD12b5XBv9UWzU+rwX55roXxmjXLGGz6WJ1jXRj1qfR+1o7nEj/zTHavxIv7qXnWuU4wq951it5bG7/qhm2zsZyl3JIW3PZjl6Wq1nMVqTpmely2wvpq0zIPCqBPy5Xe2x/GSfFFltSnrV1LxnpavYXo62PosfxY72Yl5p47rPpcms6zI/i9FapteaNKdW+TK7mzvL5WsreT2u+Rpx3efS9Kxro9+LYR0CEIAABM4J+HfuebZaBq+Z+bUsP6osh6/9KP96rsv8vxE/K5leaz+q3572M/tb+TPLtFr7UeWedJnNI8bvc70YrWd1RmuK69mT2JZzNX5V3+u7t76af1XvdU9iPc+KP6rZ9k6GcldySJvZXnym9bVeXFt3Xebvxo7i2IPAMwn4c77aR/mb4KTIalPSq6bmPSvdzPbitX4arzyyWT7tNZvta8110ZemZ6Pe572Ytn7VuKPGKGe19yyHnznb15rrMl+6zGZ61iAAAQhA4BoC/r17TcZ5Fq+Z+fMMvxVZDl/7rf47c230/6p/VqLW5z+q355rov9b+TOLOp//qHLPtdHPI/bfr1q+WKMy7/UxyzeK096ovjRuV/UeW/FX8o+0bW82RvGz2N39Uc1Kz726nren8XXXZ75r5Wc6X5Mus67L/CxGa5lea9JgIfBqBPSMNrs6yhEnRVabkl41Ne9Z6Wa2F6/1WXzbr44sV4zNNFqL2jiXLrNR6/NM72uu3fU9X/TvyNlqzEbsI4vJNFrbyV+NneVmHwIQgAAE+gT0Xdvso4bXzPzVPrIcvjbL59roj2Kj1ue9ONdE/8oY5Yo1fC5NtK6JftT25jFuNO/laOu7ccq5Gr+qV52qXck/0ra92RjFz2J390c1Kz336nrensbXXZ/5rpWf6XxNusy6LvOzGK1leq1Jg4XAqxHQM9rs6ihHnBRZbUp61dR8ZKXt2VGs9nqxvi7tyLre/Rjje9GP2jiPep9Hrc9dl/mu3fWzvFq7I2fL3Ruqm9kYk2m0FrVxLl1mo5Y5BCAAAQhcR8C/d6/LOs7kNTN/HP13N8vha38jfq+4Nvq/lb9nUevz38qfmWui/6P67UWdz38r/85cG/2/6v+uRJ3PezHZusfN/Cy+rY3iejG+vho/0re90zHKn+Ve1XuOk1jPs+KPara93eF5Kzlcn/lZjkzna1mM1lyX+dJlNtNrLdOzBoFXIKBntNnVUY44KbLalPSqqfnIStuzo1jt9WK1Lt3IShttFhM1Ps/0vuba6Lsu+lGbzWNMnHtM3Gtz349+pq+sxTxx3ssRdZpneu1lNtP7WhajNdfhQwACEIDAtQT0Xdvso4bXzPzVPrIcvjbL59roj2Kj1ue9ONdE/8oY5Yo1fC5NtK6JftRq3nS9EXNk853YXoyvZ7W05jr52utZ6UZWsZlGe5m9Qu85shpac92VvvL37G4tz1fJ4fro9+KjLs57cW09auN8N3YUxx4EnknAn/HVPvq3Rch0UiSkKk9Xa7re/WpBj8n8WZ4sRmtZrPYym+m1lul9TbqedW3m9+K07jFac+v70Xfdih/zxHmWK2p8foXec3ju6LsOHwIQgAAEriXg37nXZu5n85qZ34/Md7IcvpZH/ay6Nvo/qr9e1Gr+V/mzIk1mf1S/vUyrtd/KvzPpMvtX/d+VTKu1LEZ7zY6G66Lfi4s6n/difN310Xed/KjJ5tL2rGKyfe1l9gq958hqaM11V/rK37O7tTxfJYfro9+Lj7o478W19aiN893YURx7EHgmAX/GV/sY3xSW7aSIpVlyVbMaJH20p/HKN8sjXWaz2EyntUyvNWl6VrpmmyaOXpyvxxjNVzWub/7uiHniPMsbNT6/Qu85PHf0XYcPAQhAAALXEvDv3Gsz97N5zczvR+Y7WQ5fy6N+r7re/d+qn5lrov+j+utFrc//qs/+MGv5PH/0s3o7MZ63l1PrrnVf+9G6JvpRm81jjM8zfVtzTeb34mJspsvyaW1V3+JGQ3kzO4o72ctqxbXV/DvxMcbnvfquiX4vRutRH+fSZTZqfZ7pWYPAKxA4eU7H31x2upMilmbJVc2VIMW4rcZ7TObP8mQxvubxvp75rs38LEZr0mverA9fH/ke0/yojfuaR53PpVm1niPzs3yZztc8xtcz37WZn8VoLdOzBgEIQAAC1xDQd22zjxpeM/NX+8hy+Foln+vd78W6xv2eXuuujb40bqMmzl2b+VHv81V9i83GLKfHuFa+70dfmsxGbTbP4rSW6dua9kc2i3V9tj/L/aiY1uedwzlk/mptz1GN9Zjo93JEXZz34tp61Mb5bqzHKaev4UPgWQT0PDa7OsoRJ0VWm5JeNTWvWMXIVmKkUUzPStezvbjT9azeSs7T+KxWllNrmV5r0qxaxfdslq+nPV3fqZXFsAYBCEAAAucE/Dv9PFstg9fM/FqWH1WWw9d+lGPPY9zPonzf/Uzra67NfGmzvWxN+p7NYrSWxWivZysxmcbXYm7fi37U+jxqs7nro5/p21rU7cx3cl8VU+23V++K9VkP1RqepxrTdB4X/V6eqIvzXtysXsszGrGOzxWXrWkPC4FnEDh5Jsf/RNhpTopYmrLr9ZpfHVfFxTyaj/qQZmaVY6bTvvRutVexHud+JTbTeI7oZ/q4FmNm8xifzbMcma631uJ7e3E91or72TzGMIcABCAAgWsI+HfuNRnHWbzeyB9n+dkd5fC9n4ix5zHyswjtuc10cc31V/gxv89n+V0r/6oY5Ys25o/7Po/aOHdt5kd9nGcxWovalblyRDvLEfU+n8Xu7nuNK/1KP6N6WfxIH/eyeF+L+jb3/czPYrSW6X1Nusy6buZn8axB4BkE/FldrV/+ZeKkSLUprzHyZ/k8dqR13a7v+Ss5TvQe2/zVejG+msPrZDl28rScleG1K77nfAV97MH7w4cABCAAgXMC/j17nu1vBs9/4ivzSQ6PVb7Muq7qZ3l6a9WcKzqv9Wpx6i3rS3uymaaypnjZSow0iolW+ys25tD80TmuqKfeT+xKHyPtSg+jPL7nOX195HuM/JHe96SP1jUzP8Yyh8CzCPizutpD7a/Y8Mf6apGK3g9R8Uc5FV/RSLtrY41Rnqht81V9zHEaP+vB88faPndd1ff4zK/mcV3M43vRj9o2jxqfR73vVf2YgzkEIAABCJwR8O/fs0x5tOc/8ZX9JIfHKl/Punbm93KM1ldyzrTaVz3N77Sq1exuHc8hfzdXi9PYyaHYaFdyxVjNV3JIq9hotT+yLWa073sx/x1zr7fqr/SzmntXr55O45Wn2Uou1+ND4NkE/Jld7eXn23oSeVJkkvr/tz1/xR/lVHxF07Q7YxRXqR9r7sTEHJ84H3G+87zPqnvnmcgNAQhA4BMJ6P6863vb86/U8Djn3lt3Tc/fjfU4+b0aq+vKJ7sb73GVXNI0q6E1zaPVvsc0jdZHemlibCUmarym8rpGa6NarhnpenkV7/s9X1rZqNO6bNzP5tK6jTrtxfVnztXTzO72OMsb91uduDabe28zbdz32OhHreZRxxwCr0BAz2ezq6MccVJktSn0EIAABCAAAQhAoEKA95MKJTQQgAAEIAABCIjAybsDP6CIIhYCEIAABCAAgbcjcPIS9HaHpWEIQAACEIAABI4JnLw78APKMX4SQAACEIAABCDwLAInL0HP6pm6EIAABCAAAQg8j8DJuwM/oDzvc6MyBCAAAQhAAAKHBE5egg5LEw4BCEAAAhCAwBsSOHl34AeUN/zAaRkCEIAABCAAgf8SOHkJgiEEIAABCEAAAt9H4OTdgR9Qvu954cQQgAAEIACBjyFw8hL0MRA4CAQgAAEIQAACZQIn7w78gFLGjBACEIAABCAAgVcjcPIS9GpnoR8IQAACEIAABO4ncPLuwA8o938+VIAABCAAAQhA4CYCJy9BN7VEWghAAAIQgAAEXpjAybsDP6C88AdLaxCAAAQgAAEIjAnoJWisYhcCEIAABCAAAQj8l4DeHZpdHeWIkyKrTaGHAAQgAAEIQAACFQK8n1QooYEABCAAAQhAQARO3h22fkDRLzVeGP+f/8AABjwDPAM8AzwDPAM8AzwDPAM8A9/zDLQ/yPzzHs2zP948tvk+4p7mrun5URvnvbgr1rNabW00ZvujWPYgsEoge0arOcZPsmVREez3XAh81nzWPAM8AzwDPAM8AzwDPAM8AzwDPAM8AzwDn/oM2E8eJZcfUP7hH4ZP/YeBc/Fs8wzwDPAM8AzwDPAM8AzwDPAM8AzwDPAM9J6B0q8mJtr6AcXicSEAAQhAAAIQgMDTCOiF6GkNUBgCEPj//wnLDEP7Z3U0sv1szXP0/vnvrbfYWU7PX/GzWlpz67nUg6zvydeecmi92WxN+9pz6zG+/mh/1OOje6EeP6joGdBzWbXjbzLLogLNMiAAAQhAAAIQgMArEND7ySv0Qg8QgAAEIAABCLw+Ab077Py2Uf415KTI6yOkQwhAAAIQgAAE3pGA3k/esXd6hgAEIAABCEDg8QT07sAPKI9nT0UIQAACEIAABJ5I4OQl6IltUxoCEIAABCAAgScROHl34N9AedKHRlkIQAACEIAABM4JnLwEnVcnAwQgAAEIQAAC70bg5N2BH1De7dOmXwhAAAIQgAAE/kdAL0H/W8CBAAQgAAEIQAACAwJ6d2h2dZQjToqsNoUeAhCAAAQgAAEIVAjwflKhhAYCEIAABCAAARE4eXfgBxRRxEIAAhCAAAQg8HYETl6C3u6wNAwBCEAAAhCAwDGBk3cHfkA5xk8CCEAAAhCAAASeReDkJehZPVMXAhCAAAQgAIHnETh5d+AHlIXP7QT0QhmkEIAABCAAAQgUCXA3F0EhgwAEIAABCEDg/wno3WEHx0f9gNJA3DUE2e1dtcgLAQhAAAIQgECNwCfeyzpTjUBdpbxu69E1peeWX4tcUym37DWt1DgAACAASURBVFo0aghAAAIQ+GYCJ3dH+RcHFWn2Fced/Xnu6L8iC3qCAAQgAAEIfAsBv5ff/cx+luZfNWLebH5aK8sZ196hxmmPxEMAAhCAwOsT8Ptptdvy7XxSZLWpHf2d/Xnu6O/0SgwEIAABCNxDoH1HM76LgN/L73pyP4P7V5zH88383XqzvL7/yjV2eyMOAhCAAATei8DJvVR+0zwp8gicd/bnuaP/iLNRAwIQgAAEagTadzTjuwj4vfxOJ/e+e/7pebK8MWfUxP3KfJZjtj+rEePbPI6oifvMIQABCEAAAiLgd4bWqvbvDdSJPCnSSXnZsvcm/7Lk/yZSXrdX1yAfBCAAAQicE2jf04zvIfCO97L3PPJPP0Xlbnma3xvSyfZ02bpiZDNNW9O+bE+XrStGNtOc1ujlZB0CEIAABD6PgO6TZldHOeKkyGpTq3rvTf5qDvQQgAAEIPAZBHYuw884+XeeQvf+O33uWc++Jv/kE1UO2VEuadyO9NpzffNnY1Xf8q3ErGhnvbIPAQhAAAKfS8Dvi9VTzm+7fzOeFFltakXvfUV/JQ9aCEAAAhCAAATej4Df/e/X/U/Hfg75P7trnuLdzjK4tvmzEfU7Ma9QY9YD+xCAAAQg8HkE/A5bPd38hvw340mR1aZW9N5X9FfyoIUABCAAAQhA4P0I+N3/ft3/dOznkP+zu+YpXrYSLa3bUZzrml8ZqzGr+tbDTkyldzQQgAAEIPA5BPyuWD1V7cYLF9Jqkbv0fvCef1dt8kIAAhCAAAQg8HwCfv8/v5v9Dvwc8neyKdZtNY/HNH80VrSepxoXdbN+VGM3TvFYCEAAAhD4fAJ+V6yednw7WraTIpbmUjf25HP5lxYkGQQgAAEIQAACL0VA932z7zz8HPJ3zqNYt9U8HtP83oi6kTbmiLFxX/Ooq9bYjVNdLAQgAAEIfD4B3RU7J+3fjiGbilQvsBB++TTrx9fcv7w4CSEAAQhAAAIQeAkCn3Lf+znk7wBWrNtqHo+Rn8Vqz22my9Y8Rv6J7urYLB9rEIAABCDwWQRG98/spB/1A0o7rGC4nUFY2d/JO4oZ7a30hRYCEIAABCDwjQQ+5R71c8jf+TwVK7uSQzFus3jfb/7KiLG9+Kouq30Sm+VjDQIQgAAEPoeA3xE7pyrfeqeFdpobxfT68XX5ozzVPeVyO4t1rXzFaB6t9rEQgAAEIAABCMwJ+D06V7+uws8hf7VbxbldyeFx8rN47clmmt6aYtxmWt9v/sqIsavxK7XQQgACEIDA+xHQPbHTeflGUpFXuIRGvfie+1tw/vkn/TdalDfLqb2eVUxvv60zIAABCEAAAhCoEfD7tBbxmio/h/zVThXndiWHx8mP8Vp3GzWjucfJj3qtu42a0dzj5I/07EEAAhCAwPcQ0L3Q7M4oR50W2mmuFzPrxffl93L11hU3slnsSN/22qhostysQQACEIAABCDwm4Dfqb933mvm55C/egLFuV3J4XHyY7zW3UbNaO5x8qNe626jZjT3OPkjPXsQgAAEIPA9BHQvNLszylGnhXaay2IqfbjG/SzfaK3Fange+dpzO4tR7Mh6PnwIQAACEIAABPoE/D7tq15/x88hf7VrxbldyeFx8mO81t1GzWjucfKjXutuo2Y09zj5Iz17EIAABCDwPQR0LzS7M8pRp4V2mstiKn24xv0sX3XN88ivxErr1uN8vfkMCEAAAhCAAARqBD7pDo1n2XknOM2Rxcc+Mk3t0/qvqhJf0YxqnsaPcrMHAQhAAALvTcDviJ2TlP9iPy2001yMWenBtfJjvpW5critxLtefhbX9hgQgAAEIAABCKwR0N367veon0P+Gon8fyK8kkN1o/Ucca/NV0YlvqIZ1TyNH+VmDwIQgAAE3p+A7omdk5RvPRVp9lljpQfXur/bu+eQX8klrWwlBg0EIAABCEAAAjUCn3K/6hxuawR+VB4r/2d37inGbYzyPflRM5orxm3U+578qBnNFeN2pGcPAhCAAAS+i4Duh51Tl38NUZFmnzXUQ7W+9G6rsVHnOeRHTZxL5zZqmEMAAhCAAAQgsE/gU+5YP4f8VSqKc7uSw+Pkx3itu42a2dxjmx9H3M80Mcbnp/GeCx8CEIAABD6PgN8Tq6f7e2t1MpwU6aRcWt6p7zHuLxX+V+zx8md5pHM7i2EfAhCAAAQgAIE6gUfdsa3OncPPIX+1nuLcruTwOPkxXutuo2Y299jmxxH3M02M8flpvOfChwAEIACBzyOge2LnZH9vrU4WFVm9xDrplpd363uc/OXinf+vh2d5VM/tLIZ9CEAAAhCAAATqBB5xx3qNql8/wX+VWd5H56j0UNGM+q7EVzSnNUbx7EEAAhCAwOcS8Dtm55Rv8QOKH/IqfxVWVneWYydmlpN9CEAAAhCAAAR+COiu/Vm53lONVbvSSZZ7Jb5pT3NU4iuaUd+V+IrmtMYonj0IQAACEPhcAvGOWT0pP6AUiUXQbT4bOzGznOxDAAIQgAAEIPBDwO/an9VrPa+x4q90keVdiW/a0xyV+Ipm1HclvqI5rTGKZw8CEIAABD6bgN8zqyed/wrwb8aTIqtNuf6Kup5Dvteo+IpzO4tzrfxZDPsQgAAEIAABCNQJ6H5t9q7hNar+ai9Z3kfnqPYQdSt9xtg2z0bUZZre2klsLyfrEIAABCDwOQR0T+ycKL+1kkwq0uwjxxV1PYf7K+fwOPmzeOnczmLYhwAEIAABCECgTuARd6zXkB871Lps3J/NFed2FpPte3zzV0Y1NupW6sTYXn9VXRZ/EpvlYw0CEIAABD6LgO6JnVOVb1YVafZR48qankv+yjkU43YW71r5sxj2IQABCEAAAhCoE9D92uw7Dz+H/J3zKFZ2JYdiZHux2nfb08Z1j2l+b0TdSBtzxNi4zxwCEIAABL6bgO6JHQr9mytkU5GVCyykWJ5eWdNzuV9tymPkV2Klla3EoIEABCAAAQhAoEZA92uz7zz8HPJ3zqNYt5U8rpffi9O+257W110v3/fd175b3+/5rpff07IOAQhAAALfSeDkfii/bahIs48aV9f0fPKrZ5He7SzWtfJnMexDAAIQgAAEIFAnoPu12Xcefg75u+dRvGwlj7SysxjpZGf6ti+t7CxGOtmZfqdGJScaCEAAAhD4LAK6V5pdHeWIkyKrTTX9HfU8p/xqb9K7ncW6Vv4shn0IQAACEIAABGoEdLfK1qJeU6UzuN3t1HM0vzJWY1b1rYfVmFX9To0KGzQQgAAEIPA5BPxu2TlV7VYNl95OodWY04Nl9Tyn+5k2rrleftTEuXRuo4Y5BCAAAQhAAAL7BD7ljvVzyN+nsvZjherJVutKLzuKk0Z2pPU96WV9L/rSyMZ95hCAAAQgAAHdEc3ujHLUaaGV5rzW7sGyejGvzzO9r7lWvu9nvnRuMx1rEIAABCAAAQjsEfiUO9bPIX+PyH+jlEN2lEsa2ZE27imm2dFw3UzreVbiVrReAx8CEIAABL6LgO6LnVOPbzvLqCLN3j281pX1Yt44H50rajW/OmaUjz0IQAACEIAABH4T0H3c7DsPP4f80/Moj2zMp3W3UTObe2zzsyFN2+tpsjitKV5W6259b6eG58KHAAQgAIHPJaD7YveuyG+6hNdpoSTlryXPP/J/BU0mozyjPaUdaeLeTkzLwYAABCAAAQhAYJ+A38f7WR4f6X1X/Z0uq7mbbnd8So3d8xMHAQhAAALvQ8DvrJ2uy7flaaFZc55/5M/y+P4oz2hPOUaauLcT03IwIAABCEAAAhDYJ+D38X6Wx0d631V/t8tK/t3civuUGjoPFgIQgAAEPpOA31c7Jyz/BX9aaNac529+G3FN67Ncvq8cvrbiV2pWNLHmTkzMwRwCEIAABCDw7QR0z7/jvareR59hRTOK9z3lcuv7V/ieW/4VeT2H8rr1fXwIQAACEIBAj8Dp3fEyP6D0Dsg6BCAAAQhAAAIQ6BE4fRHq5WUdAhCAAAQgAIHPI3D63sAPKJ/3THAiCEAAAhCAwNcQOH0R+hpQHBQCEIAABC4lwP1zKc6HJtNnt1OUH1B2qBEDAQhAAAIQgMBLEDh5CXqJA9AEBCAAAQi8JYF2/7Qh+5aH+NKmT94d+AHlSx8ajg0BCEAAAhB4dwJ6AeLl9d0/SfqHAAQg8H4EdAe9X+d0fPLZ8QMKzw8EIAABCEAAAm9JQC9A/IDylh8fTUMAAhCAAAQeTuD03YEfUB7+kVEQAhCAAAQgAIGrCOhF6Kp85IEABCAAAQhA4HMJ6L1h97984QeUz302OBkEIAABCEDg4wnoRejjD8oBIQABCEAAAhC4hMDJuwM/oFzyEZAEAhCAAAQgAIFnEDh5CXpGv9SEAAQgAAEIQOC5BE7eHfgB5bmfHdUhAAEIQAACEDggcPISdFCWUAhAAAIQgAAE3pTAybsDP6C86YdO2xCAAAQgAAEI/Pf/+8j2IsSAAAQgAAEIQAACFQL8gFKhhAYCENgmwB8n2+gIhAAEbiZw8hJ0c2ukhwAEIAABCEDgBQmcvDuU/ysbFcH+8x8YwIBngGeAZ4BngGeAZ4BngGeAZ4Bn4NWeAf2tmvXV9rQuXbTal/UYX4vryiNNs1Hja76XxXoexWUxiu1Zj+1plLe3X83Ri2f99Qjo+drpjB9Q/uGLXw8QlmeBZ4BngGeAZ4BngGeAZ4BngGeAZ4BngGfgG54BfkDhx5D//ar8LQ/8M8/Z/oF7Zn1q/1xsfBY/LHguYPEtzwD/3POsf8uzzjl51nkGeAZ4Bu57BlZ/RNn6N1BWi6CHAAQgAAEIQAACdxDQS+UduckJAQhAYEagfQe1ITvT37Gf1fbvxtG+dLLen8dpf8Uql2LiXOvNtuHW9+THeMVof8UqVy+H9ldyor3vR4672OpzXrH8gLJCCy0EIAABCEAAAi9FQC9VL9UUzUAAAhCAAAQg8LIETt4d+AHlZT9WGoMABCAAAQhAYEbg5CVolpt9CEAAAhCAAAQ+j8DJuwM/oHze88CJIAABCEAAAl9BQC9AzTIgAAEIQAACEIBAhYDeHyraqCm/cagILykRIXMIQAACEIAABJ5BgHeTZ1CnJgQgAAEIQOB9CZy+O/ADyvt+9nQOAQhAAAIQ+HoCehH6ehAAgAAEIAABCEBgSkDvDc3ujHLUaaGd5oiBAAQgAAEIQAACPQK8m/TIsA4BCEAAAhCAQEbg9N2BH1AyqqxBAAIQgAAEIPDyBE5fgl7+gDQIAQhAAAIQgMClBE7fHfgB5dKPg2QQgAAEIAABCDyKwOlL0KP6pA4EIAABCEAAAq9B4PTdgR9QXuNzpAsIQAACEIAABBYJnL4ELZZDDgEIQAACEIDAmxM4fXd4mR9Q/CDN742oG2l7OViHAAQgAAEIQOAzCOi94N1PU3mf0Vkr2oyHx8vPdCdryhvtSU5iryXgn821mX9n8zrNv2t4nbtqkBcCEPgcAqffGeVvs9NCM+Sef9Wf5WYfAhCAAAQgAIHPI+DvC+98Oj9H1V85byXnSr5M+4gaWV3W5gQqn40082xzhXL17DzDXNHLrfV5BhQQgMC3EtD3RLM7oxx1WmjWnOcf+S1P3J/lZh8CEIAABCAAgc8j4O8D73o6P0PVXzlrNWfT7Y5H1Njt7dvjVj4baU+YKcfMvnqNk/6IhQAEXpuAfz/tdFq+LU8LrTbn9dxfzYMeAhCAAAQgAIHPJaB3hHc9ofpfsdWzxpxZXEWTxWmtEl/RKB/2OgKR+8p8p4uYP+aY7Ud9Np/liPttzoAABCDgBPx7wterfvlb5bRQtSHXeU35vo8PAQhAAAIQgMB3E3jn9wP1vmKrn3bMOYpb0XqelbgVrdfY9VVvN/7d43T+zLazZetxbYVBNbaqy2pXY6u6rAZrEIDA5xPw74id0/IDyg41YiAAAQhAAAIQeDqB05egZx/grv49r/zRWaWRHWl9T3pZ34u+NLJx/+r5o+pc3fdV+XR+2VFeaTI7itNejNN6Zle0Hr8at6r3WvgQgMBnE/Dvh52T8gPKDjViIAABCEAAAhB4OoHTl6BnHuDO3j138ytjNWZV33rYian0HjVeJ+59w1znb2dtfmUoJtqd2FnMq9aY9c0+BCDwGQT8O2jnRLVv1XDp7RTaifHDyd/JQwwEIAABCEAAAp9HQO8Gzb7bUO939K3cspUa0rodxbmu+ZURY6pxldyu8Tq+/i3+7vk9Tv6MmXRuV2NW9a1WZXhP1ZhKXjQQgMB7E/Dvhp2T1L6B+AFlhy0xEIAABCAAAQjcSOD0JejG1oap7+zbc8sfNvPvprRue3Gukd/T+rq0bn3/Ct9zN//bhp9/9ewe6/4oj+uaXx0rcStar78b5znwIQCBzyPg3w07pyt/050W2mrun38e9q977vRHDAQgAAEIQAACzyPwjHeTK057Z9+eu/kroxobdSt1YuxKfxWt56/oP01zen6Pl99jpH23PW1c95jm90bUjbQxx0lszMUcAhD4HAL+3bBzqv43Vsh2WiikK029pvxSICIIQAACEIAABL6CwDu+H6hnt1d9WJ6z+SujGlvVZbVPYrN8vhZzt/k3jpNzrzBc0cbPIcbGfc2jrs1XRoxfiUULAQh8JgH/Xtg5Yflb6LTQVnMP+DdQ/FzRH/WcabU2itOetM3G4Xvyo2Z1rjyZXc2FHgIQgAAEIPAKBPxOe4V+Kj14zz2/kifTZPkyXW+tGh91vXzZeoxt89OR5dSacse51kdWMW5H+k/Y87PK751L+2572rjuMfKjps215zbT9dY8rvkMCEAAAv69sEOj/E1yWmiruZt/QPEz9fys755W66sxrleOnnVtxe/lWV2v1EIDAQhAAAIQeCQBv8seWXe3lvdb8VfrZDlXclTiK5pRzdN4z53lqqx5juhX4qWJsZ8w19lkR2eSRnakjXuKcRs1be778jNdb00xbnta1iEAge8gcPp98NU/oLRHxAFmfvYYZTpfizG+l/nSZ3vZmvQz67GZ1vdnfhbPGgQgAAEIQOCZBPzuemYf1dre74p/kr8a23RZTzG+ookxPj+Nn+XK8sc1z+F+1FXnnuPd/Xjm3nmirs1XRiW+opnVvCLHrAb7EIDAexHw74WdzsvfdqeFtpq7+d9AiT35GeVHTZxL57aiiXqfz/yYP5vHHJlGa1Hb5m34urRYCEAAAhCAwKsQeLd7yvtd9SvMs5yVOGkq8RWN8mX2ND7LqbVRbu1JG632ZeO+5tp3q71PsH6u5vdG1I20WY5KfEWT5fa1K3J4PnwIQOD9Cfj3ws5p+t+MIdtpoZCuNPWa8kuBmyLVcFtJ5frmj0bUxnkWGzWzGi1HjMnyai1qK/kVi4UABCAAAQg8i4DfX8/q4aq6fpaeP6uVxc1ifL8SX9F4zuifxsd8Po+5fW/kr8ZFfZt/wojnGp0palcZVOIrmlGPbe+KHLMa7EMAAu9FwL8Xdjovf+OfFtpq7g3+DZR2LmfT/NGIWs1XY1b0I6321Ies1rEQgAAEIACBVyWgO6vZTxl+pujPzhj1q1wq8RXNqM/T+JXcI632Yj9aH9kYo/ko5h32dA7ZUc/SuB3p457Hyd/RxJg4V263UcMcAhD4LgKn3wflN47TQjsfi9eUv5OnGqMabiuxrm/+aETtTN9yrcZE/agf7cWYSl+KxUIAAhCAAASeQcDvrmfUv7Omn839UU3XyR/p455i3O5oYozPPbd83z/xlc/tLJ9rm18dMW4ltlrjkbp4nlntqF89fyW+orm7z1l+9iEAgfcj4N8tO92Xb4rTQlvNfcm/gVJh4/ybPxorWuWJMbMaisNCAAIQgAAEnkXA765n9XBnXT+f/FE9adyO9HHP4+TvaGKMz5XXre+f+J5T/iifNG5Het/zGPdd807+6hlcL3/lvIpxG+N9T37UzOaKczuLYR8CEPhsAqffB+O/xI3daSFLVXa9pvxy8IZQNdxW0ri++aMRtTO9csU4rWd2Rav4GFPtS/FYCEAAAhCAwKMJ+N316NqPqOfnc79X2zXye9psXTFuo8735EfNaK4YtyP9yp7nlD+Kl8btSB/3PE5+1LzDXL3LVnqW1m0lThqPk689Wa271V7Veqz8aiw6CEDgMwnou6DZnVGOOi201Rz/Bsr/sDn/2Ye9olWBGDOroTgsBCAAAQhA4FkE/O66swfVqdS4+v5Ubbe9Plwjv6fN1hXjNup8T37UjOaKcTvSr+x5TvmjeGncjvRxz+Pcj7pXn+/07jHyV86pGLcx3vfkR81srji3sxj2IQCBzyZw+n3ADyj2fDhM+bbddaWV7Qo3/u+ZKJdyy2o9s9LIZpq4Jq3bqGEOAQhAAAIQeCUCj7qzvE7Vv4pTVq+Xe0Wb5ajEVzRZbq2dxitPZldyZ9q2tjKuyBHr9XLO1mOe6tzzVmOazuPkXx2vvG5XalzR52o99BCAwOsTOPlOaacr3xSnhXZQek35O3mqMarhthLr+uaPRtTO9C3XasyqPqsxOgN7EIAABCAAgVcg4Pfdnf14nap/ZT+xZi931LX5yqjEVzSjmqfxK7lXtKusWu7sLDt5vM9eztm656j6nrMaI53HytdexSrGbYzzPflRM5srzu0shn0IQOCzCZx+H5Rv1tNCOx+D15S/k6caoxpuK7Gub/5oRO1Mr1wxTuuZjdpKjRiT5WUNAhCAAAQg8EoE/O66qy+vseJf2U+s28sddW2+MirxFc2o5mn8Su4V7Sor5b76PFm+ypr6qVrlbPrmrw7Fu13J4XHyY7zW3UbNbO6x8mcx7EMAAp9NQN8Fze6MctRpoa3m+L+B8j9szr/yYUf9KCZq/1cUBwIQgAAEIPDCBPz+urNNr1P1r+wn1uzljro2XxkxPouNmjtqZHUra7G3UUzUaj6KyfYU5zbTVdc8z4pfzd90Me9KrGtP8sTYNs9G1GWa3lqM7dXoxbMOAQh8JgH/btg5Yf5tlWQ6LZSknC55TfnToIlglEd7bifp/n/b9c0fjaid6ZUrxml9ZGNMrDXbH+VmDwIQgAAEIPBsAn6P3dmL14l3qepWNNKu2mruqOv12qsf4091MT7mX+0v5ovzmD/u+zxqNXdNxVec20rcSOO5qv4on+/FfL636sdcbV4d1diqLqt7EpvlYw0CEPgMAv7dsHOi8jfdaaGt5m74N1B0jqwf7bnNdHHN9c0fjaid6ZUrxml9ZmPcaD7LxT4EIAABCEDglQj4nfZKfV3dy8o5Xdv86ohxo9gVrdffjfMcIz/mH2nbXtS3+eq4IsdqzV29em3xO2eNdZXPbdT05h4z6iXqRtpY6yQ25mIOAQh8DgH/btg5VfmmOC201dzFP6DMzuD78it9Sys7ipHG7UivPdc3vzoU1/Tyo63mQgcBCEAAAhB4JQJ+n71SX1f3snJO1za/OmLcKHZF6/V34zzHyI/5R9q2F/VtvjJO41dqnWrVa8uzes5ebeV029P6uuvl+7772nfr+yPfY5rPgAAEINAI+HfDDpHyt8lpoa3mwg8oOzk8ZnYG35fv8T1fWtmerq1L43ak157rmz8brp9p2YcABCAAAQi8I4FvuOv8jJX7v32Oj4h5RI3VZ3K1p6jXvFpXerfV2EfrTnpUbNaz9mQzTVyTVjbux7l0snG/N5detqdjHQIQ+C4C+k5odmeUo04LbTUXfkDZPWSrXenfNfJnfUvndhTjOvkjvfakldV6ZqWRzTSsQQACEIAABN6dgO65Zj917JzRYypson4nZsZ/p8YsZ9yPNeJ+No8xlbMrT4zV+qtZ9bnTl2J7XHxf/qyOdLKr+hY3G8rtdhbDPgQg8B0ETr8X5t9A/3I8LbTzcXhN+Tt5Wozim+0N18jvaWNO6WV7cdp329P6uuubPxor2lEe9iAAAQhAAAKvTMDvu1fus/Xmvc7ucZ3FY7RWsR5XqbWqVw8rcSta5V+1sUabxxHXspioiTk0j7FafyXrPa72VY11XfNHY0XreVbjVvVeCx8CEPhsAv79sHPS8becZTwtZKlKrtdzvxQcRB7f/NGIWs0Vo3nVKk42i9PeyMa4FW2MrcxH+dmDAAQgAAEIvAIBv89eoZ9RD96r/J5e+2572t66xza/N6q6LL4aG3WjfrI61bVRHd+L+XxPftTEuXSycf/Zc/V1lR2dJ9a4Sut5Tmq0WAYEIAABEfDvE62t2PI3ymmhlaaa1utFf5Yr6n1+Eut5mq8R1+O8opMmszGf5pm2rWn/CturwToEIAABCEDg2QT8nnt2L7P63uuKP8s72o91ojbut/nqiDlifNzfqRFz9uZZrWwtxmeaUZ9RH/M9ex77u2I+O1OskekrmixOa5X4qGlzBgQgAAEn4N8Tvl71y98qp4VmDXn+O/1ZH22/Ut/zVPRVjfJW9dIpTlbrV1jlxEIAAhCAAAReiYDfca/UV9aL91r1szyra9VaTbc7HlGj2lull16uSmzU9HI9az32d8W8epaVWtWcUbdSo2kZEIAABCIB/x6Je5V5+ZvltNCsGc9/pz/rQ/ujHqSRda3W3Pr+zFfcTBf3FSfr+23N56u+cmIhAAEIQAACr0TA77NX6qvXi/c783s5dtYfUesRNapn7/VSie/FZuuVfI/WZH2erq2coVJrJV+mrdRoGgYEIACBjIB/h2T7s7Xyt8tpoVkjbV81KtpX07zaF/UJS8W6fTXe9AMBCEAAAhB453vqGb0/oqbXaP6zhvex04PHu7+T65ExV/WqPLu9K97tbq5enOd2v6dnHQIQgEAjcPp9Ub7ZTgvxcT2OwFWf1VV5HndyKkEAAhCAwDcR4J76pk+bs0IAAhCAAATOCZy+O/ADyvln8FIZTh8IP8yVuTwvPgQgAAEIQOAKAtxTV1AkBwQgAAEIQOB7CJy+O/ADyoc9K6cPhONQLl/DhwAEIAABCLwKAd1TzTIgAAEIQAACEIDAjMDpsiMrPwAAIABJREFUu0P5jeO00Owg7J8T8M/o9GXSc513RgYIQAACEIDA9QS4q65nSkYIQAACEIDAJxM4fXfgB5QPejr8YZC/czzFyu7kIAYCEIAABCBwNwHdU80yIAABCEAAAhCAwIzA6btD+Y3jtNDsIOyfE/DPKPNHFTJ9W2NAAAIQgAAEXpWA312v2iN9QQACEIAABCDwOgRO3x3KfyGfFnodZJ/diX9Op/5nk+J0EIAABCDw7gT8nnv3s9A/BCAAAQhAAAL3Ezh9d+AHlPs/o6dU8Adj1X9KwxSFAAQgAAEILBLw+20xFDkEIAABCEAAAl9I4PTdgR9QPvyh8Qdk5H84Bo4HAQhAAAIfSMDvtQ88HkeCAAQgAAEIQOBiAqfvDss/oFzcP+keTMAfmOYzIAABCEAAAu9MQPfaO5+B3iEAAQhAAAIQeAwBvTfs/i1c/gv6tNBjcKxVqUKr6taqP1btn5/7u13czcTzy/e+s7XsLIrJ9lbWsjyxh5hPMbJtXzHyfa5412stWsXJxn3NZ7lm8crz6lbnlG39un9F/6esKvFX93zFuckBgVcnwD83r/4J0R8EIAABCEDgdQjovaHybp51vfUDihfF/+d/f6jBAhY8AzwDPAM8AzwDPAM8AzwDPAM8AzwD689A+2P11bn5H9SzXv08ilOM5pltGg3pfU17K/Y0Pta6Ol/Mf+f8lOnPpzPp0gvhr38hwAxmPAM8AzwDPAM8AzwDPAM8AzwDPAM8AzwDPAOv8QxMfgJJt7d+QFEmPvi1D75xeyYzfW7eh/vqTTrN3WZ631/xd3Opv934lR53terR49Vvtue6K3zVirl661FXnV+dr1r33XX+DLw6w1ft766+7sr77s/sq/XfPqc2Yl//Lv9Zjzrma+8v8IIXzwDPAM8Az8AnPgN6b1ixRz+grBRCCwEIQAACEIAABK4k4C9zV+YlFwQgAAEIvA8B3QWx47aukWl8XzrZ0Z40bl0ffc2brf5Huav6qq7lrWq/QSfOK/bnqZpEOcCJlG0IQAACEIAABCBwOwHeTW5HTAEIQAACEIDARxE4fXfgB5SPehw4DAQgAAEIQOB7CJy+BH0PKU4KAQhAAAIQgEAjcPruwA8oPEcQgAAEIAABCLwlgdOXoLc8NE1DAAIQgAAEILBN4PTdgR9QttETCAEIQAACEIDAMwmcvgQ9s3dqQwACEIAABCDweAKn7w78gPL4z4yKEIAABCAAAQhcQOD0JeiCFkgBAQhAAAIQgMAbETh9d+AHlDf6sGkVAhCAAAQgAIEfAqcvQT+Z8CAAAQhAAAIQ+AYCp+8O/IDyDU8JZ4QABCAAAQh8IIHTl6APRMKRIAABCEAAAhAYEDh9d+AHlAFctiAAAQhAAAIQeF0Cpy9Br3syOoMABCAAAQhA4A4Cp+8O/IByx6dCTghAAAIQgAAEbidw+hJ0e4MUgAAEIAABCEDgpQicvju81A8ofhj5oq35rlUeLAQgAAEIQAACn0HA3wne9UR+BvevPo/ndv/KOp7X/StrkOucQPtsZkOf30zX21e8bE/HOgQgAIFHE9D3UrM7oxx1WmjWnOe/y5/1wD4EIAABCEAAAu9DwN8X3qfr//zH+575V5xrVqPtn45H1Djtkfi1Z0+f6Qo3xYzsSj60EIAABK4m4N9PO7nLN+ZpoVlznv8uf9YD+xCAAAQgAAEIvA8Bf194l6695xV/93yfUmP3/MT9JrDyPEj7O0N/Jn3F9rOwAwEIQOBeAv4dtVPpZX5Aic37weRHzWyuONmZ/q79Vp8BAQhAAAIQgMC1BHS/v8s96/3u+Kv0Yo0svqLJ4rRWia9olA97H4H4OVTm1W6yXDE2auI+cwhAAAKPIODfRTv1yn/ZnxZabc7ryV/N0fSKbfYZQ3XVxzN6oCYEIAABCEDgEwnobtVd+8pn9F7lZ/1qL7OZvrcW43u6tr6i9TwrcStar7Hrq95u/CfGicmKrXKIOXtxVV0vnnUIQAACpwT8e2gnV/lXhdNCW83988+vS/00x078FTHPYHdF3+SAAAQgAAEIvDKBd7lfvc/mV0aM0Xwndhaj3LIzfduX1u0oznXNv3N4rTvrvFPuHhP/LFzj67NzrsStaGd12YcABCCwQ8C/h7biq0Gnhap1XOc1m787lGc3/jRO9U/OcNoD8RCAAAQgAIFPI/Au9+tunx4nv/IZSit7R4xyy95Ro5Iz06inZhn/JSAmV/NQXtlKfmllKzFoIAABCFxFQN89ze6MctRpoa3mwr+Bsn3If/Ps9HBFzDPYXdE3OSAAAQhAAAKvTOBd7tfdPj3O/dFn4jr5I732pJXVemalcZvp4prrm3/H8Bp35H/HnHcy8dzVz3Qn5h250zMEIPCaBPw7aKfD8u11WmiruQ/4AcW5NZ8BAQhAAAIQgMA1BPyOvSbj9VlOe/R4+aMupZEdaX1Pere+775rml8dMW4ldqdGNebTdeJ+9TmV122lhuvlV+LQQAACELiCgL53mt0ZS1GnxbYaDD+i7ORQzC4kxe/aZ3Db7ZU4CEAAAhCAwDsReIc79rRHj5c/+oykkR1pfU96t77vvmuaXx0xbiV2p0Y15tN14n71OZXXbaWG6+VX4tBAAAIQuIKAvneaXR1LESeFVhuT3mvuHFB5nmXfvf9ncaMuBCAAAQhAoELA79mK/h01fkb5vXNo321Pm617XPOzETU93dWxWb64FnuL+984v5NJzH3yLHzjZ8OZIQCB5xDw767VDvKbsZPlpFAn5XTZa1a+lF0/Td4ReI7M74T9bzmLGa39LxAHAhCAAAQgAIElAn6/LgW+kdjPKL/Xvvbd9rTZusfJjzqtu42a0dzjmn86Yr44V36taz6z0kc7i3u1/dh/nJ/0e5LrJPakZ2IhAAEI+PfPKo2lW+uk0Gpj0nvN5leGx1T00nhcxVec20pcpvEc+BCAAAQgAAEI1Aj4nVqLeD+Vn7H5oxG1M33MVYmvaGJen5/Gz3Jl+eOa54h+1I7mMfbV5qPe495q7zG+zVfGafxKLbQQgAAEnIB///h6xV/6pjspVGkm03jN5leGx1T0TeMx8j1Wa9G6ppcnxmTzmIc5BCAAAQhAAAJzAn6nztXvqfAzNn80onamj7kq8RVNzOvz0/hZrix/XPMc7kddZe7xr+RXes801TOcxLYap/HVPtFBAAIQiAT8+yfuzebjWzhEnxQKqcpTr9n8yvCYVb1ie3Had9vTtnXXNZ8BAQhAAAIQgMB1BPyevS7ra2XyM87eJaJ2po8nrcRXNDGvz0/jPVf0Y+64P5pXY6OuzV91ZL1W1irnyfJU4qQ5jVceLAQgAIFVAv79sxy7EnBSaKWOa71m82djVd/yrcRE7aynqJ/1zz4EIAABCEAAAnUCfs/Wo95H6eebvXO0U0V9JcZpVOKjxuMrfoxf7XFUI+YeaX1PcW2t0o/0bj3fq/ved8+fnSGLm8X4/mm858KHAAQgsELAv39W4pp2/ouEZTwpZGmWXK+548+KZTlHMaf6UW72IAABCEAAAhBYI+D38lrke6j9fM0fjaid6bNclRwVTZZba6fxypPZmDvTxLUrYpQj5n71ufrO7Kz3nRjPeRrvufAhAAEIrBDw75+VuKYd38Qh20mhkKo89Zo7/qxQlnM1ZqSP+Uda9iAAAQhAAAIQWCPg9+xa5Our/WzNn42or8TEnJUcFU3M6/PTeM8V/Zg77mfznZiWJ8a1+buO7Cyz82QxK+c/jV+phRYCEICAE/DvH1+v+Evf9CeFKs1kGq8pX7o296F9t77f810fc2YxK/qoreTParIGAQhAAAIQgMBfAn7P/t1975XVs7le/ioBxbmNOXyv+asjxu/k6NWMuXs6rUf9Si9Z7Eq8engVm51n1NuqPuY6jY/5mEMAAhCoEvDvn2qMdEu33kkhFVy1XrP5leExFX3TKGaml87tKMZ18kd69iAAAQhAAAIQqBPQ3drsJw0/V/VsMaYa59wqOaLG4yt+jN/ps1dnNXfU9/L21mP8lWfp1bxzfeU8K9qs59P4LCdrEIAABCoE/PunonfN0tvGSSEvuuJ7zeZXhsdU9CON5+r5q/EjPXsQgAAEIAABCNQJ+N1cj1pXqs4ssqqr5mn5qkO13VZjpfNY+dqT1bqs1qtWcW6rsTOd55Tfi9G+2562t+6x8nvad1jXGdz2+naN/J42W1eM20zHGgQgAIGrCZx879RvZfu3NFrBRw0/XLWux+z26Tlm/qhGjB1p2YMABCAAAQhAYI2A37NrkXW111jx6xV+K73G7535zGObvzpifJYjau6osZpT+pXeojY7q/L2bJZjJ4/n7+WcrXuOEz/W6eWKujZfGafxK7XQQgACEHAC/v3j6xV/6ZvupFClmUzjNZtfGR5T0bvGY90faXwv+p6j+QwIQAACEIAABK4j4PfsdVl/Z/IaK/7vLLWZ569F/FZ5fPNXR4zPclQ0o7qn8Su5V7TZWUfxbS87y04er9PLOVv3HCd+rNPLFXVtvjJO41dqoYUABCDgBPz7x9cr/tI33UmhSjOZxms2/64R62ie1dOebKbRmjRutYeFAAQgAAEIQOCMwCPuV6+x4q+eLOZejW/6mKPNV0aMz2Kj5o4aWd3KWuxtFBO1q+doubMcO3m8z17O2brnOPGzOlm+qi6LbWun8b28rEMAAhCYEfDvn5k27i/dqieFYuHq3Gs2/66xUudEe+cZ7mJDXghAAAIQgMCrEvA7+a4evcaKv9JPzLsS69qYp81XRozPYqPmjhpZ3cpa7G0UE7Waj2LinmKijbqVecxVna/UmGm9Zk/rGvk9bbauGNlMwxoEIACBOwjoe6fZ1bEUcVJotTHpvebOAZVnZFdrrOhXtKMe2YMABCAAAQhA4C8Bv2f/7l63ojqzjFWd51GMrO+t+srhtprDY+Rnsdpzm+myNY+Rn+l215RTdpRHmmhHMXEvxrb5FSPLO1q7oqbn8Fq+Hn3XNX9lnMSu1EELAQhAIBLw75+4N5svfdOdFJo10tv3ms2/euzkX4lZ0V59NvJBAAIQgAAEPp2A37PveFbvv/mnI+ZbybkSG7XVvmPcSn+VGjH/LCbqV/s5jZ/194z9eKZRDytazxPj2pwBAQhA4FEE/DtotebSt9VJodXGpPead3y57uRfiVnR6sxYCEAAAhCAAARqBPyerUW8juqu3j1v86tjJW5F6/V34zzHyF/NH/VtvjJO41dqPUobzzSqG7VVfrtxo17YgwAEIFAl4N9B1Rjplm6Jk0IquGq9ZvOvHDF3Jf9qTNRf2T+5IAABCEAAAt9OwO/Zd2Jx2rfiszNrz22mi2uub/5srOpbvp2YWR++v5o/6jX3nCNfetmR9l32dJZmK8P1d8ZUekEDAQhAoELAv7cqetfUvhn/jTgp5EVXfK9Z/VKu5o+5K/lXY1b11d7RQQACEIAABCDw+w/yd+Hh7wY7PVfiXdP82VjVt3yrMav6Wc/Z/mqNqNc8y52tSS+bad5pTeeQrfQurewsRjq3sxj2IQABCFxJ4OT7Z36jWqcnhSxN2fV68svBBaFyRtsLjTrNXd/WfEjjtrfv6/gQgAAEIAABCMwJ9O7XeeRzFN5v81fGSuyKtvWwqlffK3ErWuVftbMa2ve8WnPr+z3f9c1/teH9VXqTvmlXz6NY2VE9aWRHWvYgAAEI3EFA3z/Nro6liJNCq401vddzfydXL8bzXuV7rSyn9uOe1rEQgAAEIAABCNQI+F1ai3ieynu9wp+dJNbo6au6LL4aW9VlNVbWYp021/A9rcn6nnzt9ax0sj3dM9bVU7S9XqRr+81fHYqX7cVr321PyzoEIACBuwicfActfUOeFNo5vNdzfydXL8bzzvyWY6Zp+3HsxMQczCEAAQhAAAIQ+EvA79i/u6+z4n1e5VdOF2vFmLjf5qsj5ojxcX+nRszZm2e1srUYn2m0FrVtrj3ZTPPMNfW1Y3f7jrVinrjf5gwIQAACzyDg30er9Ze+uU4KVRrz/Ct+JfdIU6nl8SO96+SP9G2PAQEIQAACEIDAHgG/Y/cy3B/lPV7lV7teqVfNGXUrNZr27jHrZ1R/Fpvtj/I9cy/rdbZ22u8sv++f1iIeAhCAwC6Bk++ipVvspFDlcJ5/xa/knmlG9bLYqM80vhb1mrsGHwIQgAAEIACBNQK6T5t91eE9XuWvnnVWdzVfpn9EjaxuXBv1EbXZfBTve1nsK615rxX/qt4fWeuqnskDAQh8FwH/nlo9+dLbxkmhamOtRnWsaKs50UEAAhCAAAQg8D4EHvFu8j405p0+ipfqzDu6T3Hag+Ize1/X12fO+m9rd4+s7t01yQ8BCECgQsC/nyp61yx9e54U8qL4EIAABCAAAQhA4AoCvJtcQZEcEIAABCAAge8hcPLuwA8o3/OccFIIQAACEIDAxxE4eQn6OBgcCAIQgAAEIACBKYGTdwd+QJniRQABCEAAAhCAwKsSOHkJetUz0RcEIAABCEAAAvcROHl34AeU+z4XMkMAAhCAAAQgcDOBk5egm1sjPQQgAAEIQAACL0jg5N2BH1Be8AOlJQhAAAIQgAAEagROXoJqFVBBAAIQgAAEIPBJBE7eHfgB5ZOeBM4CAQhAAAIQ+DICJy9BX4aK40IAAhCAAAQg8J///Ofk3YEfUHiEIAABCEAAAhB4WwInL0Fve2gahwAEIAABCEBgm8DJuwM/oGxjJxACEIAABCAAgWcTOHkJenbv1IcABCAAAQhA4PEETt4dtn5AefwRqQgBCEAAAhCAAAT+Ejh5CfqbjRUIQAACEIAABD6dwMm7Az+gfPrTwfkgAAEIQAACH0zg5CXog7FwNAhAAAIQgAAEOgRO3h22fkBpBT9pxPMI6OiMMWakvWIv66m31uqN+hvtKTbmVozWR1bnzTSj/FEvrWzcj/NR3ao21orzmKfNfYz2sz2tVXNIpzi3cS+bj/S9vZhHut56249DMdFG3WjueWOe2bzldU0211rsQXU9XmtR25tnsXEtm8d8sa7mslHf5p7X90frrrvCH/WnHmd1ZjkUH3WVeUWjPl3rvupnut6ax8hvOXt5penZk9hezndY17l3ub3DGekRAhCAAAQgAIHrCJy8O/z9S2fQlxfC/+9LLhzgwDPAM8AzwDPAM8AzwDPAM8AzwDPAM8AzMHoG2p/Zo/3VPf+zvcWO8se9OFf8/ycJ/89oz/OEsD/TWZ4YsKqP8aN5y63/jHTZHj+gGDxBxP48ULCABc8AzwDPAM8AzwDPAM8AzwDPAM8AzwDPwCc+A9mPJKO17R9QlPQTIX7Lmdpn+GlnfYXnUj24nXG+87M4yX0SOzvzlfti3cv5Lufo9c86Lyw8A7+fAf0zH//Z1jq8fvOCBzx4BngGeAZ4BngG8mdA7w5Vu/UDSjU5OghAAAIQgAAEIHAnAX8hvLMOuSEAAQh8CoH2vfmo8chaOlOsGefSXW1bnVEt7fc0cV362KfWo226uJbNq7os9hPXIt/ZfOmfHgGbJWUfAhCAAAQgAAEIPIKA3k2aZUAAAhCAAAQgAIEZgZN3h6W3jZNCs0OwDwEIQAACEIAABFYJ8G6ySgw9BCAAAQhA4LsJnLw78APKdz87nB4CEIAABCDw1gROXoLe+uA0DwEIQAACEIDAFoGTdwd+QNlCThAEIAABCEAAAq9A4OQl6BX6pwcIQAACEIAABB5L4OTdgR9QHvtZUQ0CEIAABCAAgQsJnLwEXdgGqSAAAQhAAAIQeBMCJ+8O/IDyJh8ybUIAAhCAAAQg8JfAyUvQ32ysQAACEIAABCDw6QRO3h34AeXTnw7OBwEIQAACEPhgAicvQR+MhaNBAAIQgAAEINAhcPLuwA8oHagsQwACEIAABCDw+gROXoJe/3R0CAEIQAACEIDA1QRO3h34AeXqT4N8EIAABCAAAQg8jMDJS9DDmqQQBCAAAQhAAAIvQ+Dk3eHlfkDxwzS/jbjmc30Kvlb1Fduz1TxNx4AABCAAAQhA4PEE/K5+fPX/VvQe3H9WP9SFAAQgAAEIQKBP4OSuXvrL/6RQv/2fHc9f9RVd1Ued4jMbtaN5Fs8aBCAAAQhAAAL3EvC7+d5KeXav3/PzSFYhAAEIQAACEHgGAb+vV+u/1A8orXk/zMyPh53pfT/GzuYe2/w2tDaLZR8CEIAABCAAgXsI6C7W3XxPlTyr15bflPLd5hlYhQAEIAABCEDg0QRO7ueX+wHF4fnB3HdN5rs2+pm+uua5qjHoIAABCEAAAhC4j8Cz7mav2/xsVDRZHGsQgAAEIAABCNxHwO/n1Sr5jd/JclKok3K47PXcHwb9u+l69yuxPc1VeXr5q+uv0ke1X3QQgAAEIACBuwg84070ms0fjRXtKA97EIAABCAAAQhcQ8Dv5tWM41s/ZDspFFKVp15TfjVYerfV2KjzHM1/9lA/z+6D+hCAAAQgAIFnEtB9+Ki72etVaq7qn8mS2hCAAAQgAIFvIOB38+p5l34JOCm02pj0XlO+9mZW+mhncdm+58j2H7n2Sr088tzUggAEIAABCEQCj74TvV7zK2MnppIXDQQgAAEIQAAC6wT8Xl6Nrt38/2Y9KbTamPRe033tj6zr3R/F9PZO43t5d9ZfqZed/omBAAQgAAEIXEXgkXei15JfOYe0spUYNBCAAAQgAAEI3ENA93Gzq2Mp4qTQamOu97ryfX/kS+92pO/tncb38u6sv1IvO/0TAwEIQAACELiKwCPvRK8lv3IOad1W4tBAAAIQgAAEIHA9gZP7+Ct/QGnAVsYJ4JU6Fa33snqOSn40EIAABCAAgXci4Pfi3X17rZU7OMatxN59JvJDAAIQgAAEvo2A38urZ1/6JeGk0Gpjrve67rtm5HuM/JE+7imm2WcO70P+M/uhNgQgAAEIQODZBHQf3n1Hex35K2dXjOxKLFoIQAACEIAABK4joLu42dWxFHFSaLWxqPfa8qOmN5febU8b13diPIfH93zXR78XM1qPOeJ8FNv2Vofn81hf38nrufAhAAEIQAACGQG/a7L9q9a8jvyV3IqRXYlFCwEIQAACEIDAdQR0Fze7OpYiTgqtNhb1Xlt+1PTm0kfb0/u6x/j6zPe4qh9zVuOiLubRPOpmc8VFO4qTtqfRPhYCEIAABCBwBQG/b67I18vhdeT3tNm6YtxmOtYgAAEIQAACELiXwMld/PE/oDic6Fc+Fo+p6JvGY+THWK1H67q415vHmp5DfhY72pNeGlmt92xFJw0WAhCAAAQgcErA76PTXKN4ryN/pI97inEbNcwhAAEIQAACELifwMld/DY/oDSMflD5M7zSZXYldqb1/VjL99yPujYfjVW9csU4rUcbdZq7Tmsn1vPhQwACEIAABE4I+H10kmcW63Wavzpi/E6O1ZroIQABCEAAAhD4S8Dv5L+745WlN4CTQuM2arteX/4oUpqRrcaPdL6X1fL96Ed93I/zVX2L95iYL85dKz9qNNd+Zkca7WEhAAEIQAACpwT8DjrN1Yv3GvJ72t664tz2tKxDAAIQgAAEIHAfgZO7+O1/QGmH740Ixufye7FtXZpmq8Nj5I9ipZFd0baY2VDeila5PEa+9txqL1rXNN/34x5zCEAAAhCAwAmBR9wxXkP+as+Kc7uaAz0EIAABCEAAAucEdBfvZJr/BW5ZVajZZw3vQX6vl7ivudtZbNOuDM8tfxQvjeyKtsXMhvJWtMrlMfK151Z70boGHwIQgAAEIHAnAb+D7qrjNeSv1lKc29Uc6CEAAQhAAAIQOCdwchfP/wK3/k4KWZoj13uQ30uY7WvNbRY/289i2prHNX82VvRRO8vf02frs7XsHFlMpmMNAhCAAAQgcBcBv4seUUP1Vmspzu1qDvQQgAAEIAABCJwTOLmL53/hW38nhSzNset9yI9Jtd6sD1+X7/vN13qMjbreXPG9fa1L51Z7mXWd/EynNWmussorm+XVHhYCEIAABCDwCAJ+F91Vz2vIX62lOLerOdBDAAIQgAAEIHBO4OQu/v3rwqSXk0KT1Evb3of8mEDrzfrwdfm+33ytx9ioW5173p4/ypnFrOqzHJW1rE4Wl+lYgwAEIAABCNxFwO+iu2q0vF6n+asjxu/kWK2JHgIQgAAEIACBvwT8Tv67O15ZegNQoXHK+3fVh9tYtbrXdHGMYqO2Mvd8M3+UL4td0Y+0O3ur/ezUIAYCEIAABCAwIuB30Uh3uud1mr86YvxOjtWa6CEAAQhAAAIQ+EvA7+S/u+OVpTcAFRqnfMyuenGrytma9pr1ffnZvq/t+Modreca7bmu+VHb5qMR9SPtzl7MP+tnpwYxEIAABCAAgREBv4tGutM9ryN/Jadi3K7Eo4UABCAAAQhA4BoCJ3fx+C/w0N9JoZDqeOq9yFdSzZvNhu/Ll07zXqx0M+t55Gcx2pPNNFqTxq32Muu65l89Yv47alzdM/kgAAEIQOCzCPhddOfJvI78aj3p3VZj0UEAAhCAAAQgcC2Bk/t46a/qk0LXHnn8b2PM+vR9+epP82Z3h+eY5TnRXpl756yx91k/OzWIgQAEIAABCIwI+F000p3ueR351ZzSu63GooMABCAAAQhA4FoCJ/fx0q8EJ4WuPfJ/s3k/zdfQuuaZlaZns5jKWsw3i1nRR22bj8aqvperV+eq/L26rEMAAhCAAARmBPwummlP9r2O/Go+6d1WY9FBAAIQgAAEIHAtAd3HO1nHf4GHjCrU7CsM7yfzRz1mel8bxfb2PF5+T6t16WS1nllp3GY6rblOvvaqVnHNxuF78qOGOQQgAAEIQOBOArp/snvq6rpea6VejFuJvfoM5IMABCAAAQh8OwHdyzsc/v5VPMiiQq9y8Xs/mT84yv9vZTEnZ4v5ZvXb/kpM1M56zfSzmNizcsT1rPfV3FlO1iAAAQhAAAIrBHRPPeIO8lryK71KK1sD92FGAAAgAElEQVSJQQMBCEAAAhCAwD0ETu7jt/4BpeHU4aOtoI4xmldiM43i3WY6X3Nt80cjamf6lmsnRj14rNbc+r5838eHAAQgAAEI3E1A90+zjxher1pzJ+YRZ6EGBCAAAQhA4BsJ+L28ev6lt42TQquNVfXek/uVeNfLr8T1NMrhtqfVumubPxpRO9O3XFmM1lZqZVrlcZvpWIMABCAAAQjcReDRd5DXa/5srOpn+diHAAQgAAEIQOCMgO7mnSzzm9+yqlDlhcHCbnW9J/crRV0vvxLX0yiH2xVti/ORzT133G+x2o95tJ7Zmdb33Z/lci0+BCAAAQhA4A4CfhfdkT/m9HrNn41V/Swf+xCAAAQgAAEInBHQ3byTZX7zW1YVqrwwWNjtrve12ttJbDxYzKW5dJqvWMXKxthsXWuyMWZ1rjzRZnmihjkEIAABCEDgTgJ+F91Zx3N7zeb3hnRtf6TrxbMOAQhAAAIQgMD1BPx+Xs3ev/WTTCr0ai8B3tdqbyexCaL//RsgMW82b/HZuq/FGr7X82NMpc4dubI+WIMABCAAAQhcScDvryvzznJ53ebHMduPeuYQgAAEIAABCDyGgN/RqxX/3viDDCeFBmmPt7yv5q8Mj12J62k9X8/32J6mrWdjpO/FKM8s1vcVE61rKn6MZw4BCEAAAhC4koDfRVfmneXyujN/lot9CEAAAhCAAAQeQyDe2atV87/SO1m8WEfylGXvq/mrQ/GrcT298mU2i4m6TONrUa+5a0a+9D17EhtzjnKxBwEIQAACEDgl4PfOaa6deK+f+Ts5iYEABCAAAQhA4B4C8a5erbL0a4MXWy10t1697dRpsQwIQAACEIAABN6PgO7/V7jL1cv7UaRjCEAAAhCAwPcQ0H298+6w9MvBSaHv+Tg4KQQgAAEIQAACjyLAu8mjSFMHAhCAAAQg8BkETt4d+AHlM54BTgEBCEAAAhD4SgInL0FfCYxDQwACEIAABL6cgN4ddjDwA8oONWIgAAEIQAACEHgJAnoJapYBAQhAAAIQgAAEZgT07jDTZftLbxsqxEtKhpI1CEAAAhCAAAQeTYB3k0cTpx4EIAABCEDgvQno3WHnFFs/oOwUIgYCEIAABCAAAQhcTUAvQfyXO1eTJR8EIAABCEDgMwmcvDvwA8pnPhOcCgIQgAAEIPAVBE5egr4CEIeEAAQgAAEIQOAXgZN3B35A+YWSCQQgAAEIQAAC70Tg5CXonc5JrxCAAAQgAAEIXENA7w472fgBZYcaMRCAAAQgAAEIvAQBvQQ1y4AABCAAAQhAAAIzAnp3mOmy/aW3jZNCWXHWIAABCEAAAhCAwAkBvZvwA8oJRWIhAAEIQAAC30NA7w47J+YHlB1qxEAAAhCAAAQg8BIETl6CXuIANAEBCEAAAhCAwEMJnLw78APKQz+q1yy2+t/a6YGLcXFePW0vzteb7/OWO5u7bra/25/XUI5sTT1qb2SVx23T+4hz7fXWfV+1tRZtL4evK4db5XGd1qJ1jfuuU25fm/mK8Zxaa7Hyfb+XU9refrZeiXGNfFnl1FxW6822teqoaEea0V61h3fTjc482nu3c35qv+0z0n8+9YycCwIQgAAEIACB6wicvDfU38rDHyIqiv15cYMFLHgGeAZ4BngGeAZ4BngGeAZ4BngGeAZ4BirPQPtJwHVx7nvuj3Rxr82zoXxtr/nZkKa3n8W8w9rJuXJSnVN7IXy+FHgGeAZ4BngGeAZ4BngGeAZ4BngGeAZ4BngGeAbe8Rno/OwxXN7+AaVlfUdI9HztP9x6ula4PuLZeUSN6pnFqNdTb72af1W3U08xsqs179Dv9LITc0fvnvORPfVq9da9z3fzP/FM7/YZ3N1v9hm3tTburk3+a98l4AlPngGeAZ4BnoFnPAP/vjYsma0fUJYqIIYABCAAAQhAAAI3EdAL103pSQsBCEAAAm9AoN0Fzx5ZD35HZb7WqradsapFN/9RaueZWXrS9CHsFCIGAhCAAAQgAAEIXE2Ad5OriZIPAhCAAAQg8NkETt4d+AHls58NTgcBCEAAAhD4aAJ6CWqWAQEIQAACEIAABGYE9O4w02X7S28bJ4Wy4qxBAAIQgAAEIACBEwJ6N+EHlBOKxEIAAhCAAAS+h4DeHXZOzA8oO9SIgQAEIAABCEDgJQjoJYgfUF7i46AJCEAAAhCAwMsT0LvDTqP8gLJDjRgIQAACEIAABF6CgF6C+AHlJT4OmoAABCAAAQi8PAG9O+w0yg8oO9SIgQAEIAABCEDgJQjoJYgfUF7i46AJCEAAAhCAwMsT0LvDTqNbP6DwkrKDmhgIQAACEIAABK4moJcg3k2uJks+CEAAAhCAwGcS0LvDzun4AWWHGjElArzMljAhggAEIACBAwJ6CeLOOYBIKAQgAAEIQOCLCOjdYefI/IBSoCbAvJwVYP3nP/9xXjCrMUMFAQhAAAJ7BPzO2ctAFAQgAAEIQAAC30RA7w47Z36pH1B0ELd+KF/v+a7f9Xu5tb6b91PjxKVnP/XcnAsCEIAABJ5PwO+eR3fjtaP/6F5aPe/hGfWpCQEIQAACEHgHArovd3rd+gFlp1AlRge52lZqSzOrLR32vwRmvNo+AwIQgAAEIHAXAb+H7qoR83rNmR9jd+azGr39nVrEQAACEIAABD6ZgN+ZO+dc+utWxXYKVWKU/y57RQ+VHN+kqXxW38SDs0IAAhCAwGMJ+D30iMpeb8U/6W2ljrQn9YiFAAQgAAEIfCoB3ZOyq+d8qR9QYvM6lNuoiXPXZn7Uj+an8aPcn7oXmX3qOTkXBCAAAQi8BgG/d+7uyGut+ie9rdZqegYEIAABCEAAAn8J+J36d3e+snTDnhabt/Nb4fXk/1b0Z9Jnth/1dyfG/1Ww4gQirzZnQAACEIAABO4i4PfOXTVaXq/jfqzpe9GP2so85qjMK3nRQAACEIAABL6VgN+lqwzKf916keY/YsSaq3WzeK1V+5dethr3rTpxkv1WDpwbAhCAAAQeQ0D3TbN3jtU6rnd/pceVOGlX8qOFAAQgAAEIfBsB3ZfN7oylKBXbKbQTo3puV/N4rPvVPB7TfMaYQOQFszEvdiEAAQhA4IyA3ztnmfrRuzU8Tn6/yt8dxTTLgAAEIAABCEDgnMDp3Vq+kU8L7RzVa8pfzaO4aKt5duOq+T9NF3m1OQMCEIAABCBwFwG/d16thvfmfrVPxVT16CAAAQhAAAIQGBPQ3drszihHeaHdYqsNxpo7dbMcWqv0I61sJeabNeIk+80sODsEIAABCNxPQPdNs3cN1djJr1i3lTyr+kpONBCAAAQgAIFvJ3B6v5bfNrxQ8x8xYs3dulmeaq4Y+4hzv3MNeL3zp0fvEIAABN6PgN87d3R/mt/j5Vf6lLZZBgQgAAEIQAAC1xA4vV/Lt7IXetRlHmvu1s3yVHPF2Gs+ts/NEnlVOX8uEU4GAQhAAAJ3EvB75646J3eZ9ye/0qe0spUYNBCAAAQgAAEIjAnoXm12ZyxFqdhOoZ0Y1XO7msdjo1/JtRMT88YccR71O/OYM853csaYmNPn0vqafO1hIQABCEAAAlcT0F3T7CsO70/+rE/pRnaWg30IQAACEIAABP4S8Lv17+58ZeltQ8Xmaa9RqJ7b1cwe6341j8c0f3XE+NF8NXfTj/Jle4+o4XV36hEDAQhAAAIQqBJ49TvH+2v+bET9bD7Lxz4EIAABCEAAAj8E/F79Wa1785v831xeqPmPGLHmTt0sx0qeGF89d4yrzqv5m66aM+pOayg+5s3m0mIhAAEIQAACdxDwu+eO/Kc5vb/mz0bUV+ezvOxDAAIQgAAEIPD7b+gdHvOb/N+s8QLfKbYaE2u2+crI4k9zVOvH2r24qKv2F+N6+dt61O7WGMVlNUb6Ub/sQQACEIAABKoE/P6pxjxS5/1V7sWoX5k/8lzUggAEIAABCLwjAb9Xd/ov/yLhhZr/iBFrrtTNYlfidb6YR+sjq5imqdSU3m0lf1WvPlzf/Nk41VdqzHpgHwIQgAAEIDAi4HfVSPeMPe/tijsx5svmzzgnNSEAAQhAAALvQsDvzp2e539F/5v1tNBWc//8k/7bE95L1d+p32Ji/lmeVX1Wo+UYjatqjOpcVWN0DvYgAAEIQAACpwT8vjrNdXW89za6c1frxrw+X82FHgIQgAAEIPBNBE7vzPFf6kbytJClKrte89QvFw3CWDds/5mu6pUgxrV5Nqq609hYJ8uXre3GZblYgwAEIAABCMwI+L0z0z5y3/tq/tUj5vf51bXIBwEIQAACEPgUAqf3ZflGPy20A9xryvc8viZ/Zj2+4sd8o5iobfPqqMZGXTV/08XYXn9Rd1JjJRYtBCAAAQhAYJWA31mrsXfqH9WX15F/57nIDQEIQAACEHhnArorm90Z5ajTQlvNJf8Tnmoe7zfzd/OM4mKdkTbuxVjNZ7q4P5orZ7QxZrYf9T4/ifU8+BCAAAQgAIEKAb93KvpHaLyn5t85Yi3N76xJbghAAAIQgMC7EtA9uXs/l2/100I7gL2m/NU8istsJVeMG8WsaLM8Mb7N44iauD+bx/ira8T8s37YhwAEIAABCJwQ8HvnJM9Vsd5PdsdeVcfzxJqPqus94EMAAhCAAATegYDfmTv9/v0LvZPltFAn7XDZa8ofBnQ2FZvZTsj/lmPM/zYSZ0WbhE//JzYxf5uvjlmO2f6sXoyf6dmHAAQgAAEInBDwe+ckz1Wxz+jHa8q/6jzkgQAEIAABCHwSAd2Tze6MctRpoa3mDv4nPLGe9+9+1MW5a5vfG1E30u7meIcascfeWVmHAAQgAAEIXEHA750r8p3k8F6a/8jxzNqPPCe1IAABCEAAAicE/L7cyVO+3U8LbTX3gB9QZi84fu6RNupG2h6LWY7Zfi+vr89yzPY9V+bH+EzDGgQgAAEIQOAqAn7vXJVzJ4/30fxHj2fXf/R5qQcBCEAAAhDYIeD35VZ8Nei0ULWO67ymfN9f9ZUj2lGeqjbq2nx1zHJk+6t1shze52zftZkf4zMNaxCAwP+1dyZKjuPIEuz//+ldw3RFV1QWjsRBiZIcZm8zAUQecLFJNnd2HgQgAIFTBPy5cyrnbB71UOKK/4yhHmSf0QM1IQABCEAAAncnoOfk6vM6/ZT3QqvFZmHGmrt1a/lGOWNM6wxRp3lLX1tXjFvX+br7rhn5HiffY7Tm1vdHvscVnwEBCEAAAhC4koA/d66s08rt9Z/53LtLHy1OrEMAAhCAAATuQMCflyv9TP0Nd7fYbINeT/5sDtcrR7Suif6OtsTOjFirFp/R9GqO4kf7vdxlL8aP9OxDAAIQgAAEdgj4c2cnz0qsapfY4j9zqJdn9/FMBtSGAAQgAAEIjAjsPi+nnva7xUaHifteT37UZOeKr9lejqif0ZbYmZGpFTV3qxH7mzk/WghAAAIQgMAsAX/uzMbu6FW35Jh9Fu/UbcWqnzv00uqRdQhAAAIQgMCzCew+L9N/w/dCj3o4x5o7dWu5MvliXO8Hj1rNezG+J72s78nXnlvtZazHFb82oqaly8TWNKxBAAIQgAAEThHwZ9apnKM8qjnS1fYVO/NsreXxNc95Mq/XwIcABCAAAQi8AwF/Zq6cp/436Eqm3UKVlMMlryl/GNQQKD7ahvzf8q6+xGeH1+rFuE5+T689ad1qz63vu++alu/64jMgAAEIQAACVxF4xjNHNVfPpPiTz8grcq6ejzgIQAACEIDAnQnsPjPTf8PdLbQC0WvKP5Gn5Cj5MkN1ZUcx0kU7G9fTx9xlnhkxrhcTtZqfjunlYw8CEIAABCAwIqDnU7FXD6+1Us/ja736fja/x9RysgYBCEAAAhCAwDeB3edm+m1jt9B3y3nPa8rPR//+F5oqR7HZ4THZuBgziov6TG8xJluj5B5ppanV0Jp61DxrFYeFAAQgAAEInCDgz58T+Vo5vM4Jv1anlremK2sz2lYO1iEAAQhAAAKfRsCfnytnT39J8ELFf8SINTVv1db+yLbi43orT9TV5rXYHd1ObOyllqu2FuNOzGt1WIMABCAAAQisEvBn02qOUZzXOOXXau7kruVjDQIQgAAEIACBnwT8WftzJzdLfwnxQsW/YsQaV8xHfc/W7OWbzVX0s2O2xlX5S95RL7O10UMAAhCAAARGBPzZM9Ku7Hv+U36rj5X8rVysQwACEIAABCDwm4A/a3/vjlem/sauYuO0awrlP2lnO5mtPcqfzTfK09t/dg3vLfbie/gQgAAEIACB0wT8uXM6d8nn+U/5vT5navTysAcBCEAAAhCAwG8C/pz9vTteue0HlHHr1ykK1CuG/1iPqnFVnSv4kBMCEIAABCAwS8CfrbOxWf2znqV+NvnZntFBAAIQgAAEIPCbgJ6nq8/2qS8FKva7DVYgAAEIQAACEIDA4wno3WT1RejxHVMRAhCAAAQgAIFnEdh9b+ADyrN+OepCAAIQgAAEILBNYPdFaLsBEkAAAhCAAAQg8DIEdt8b0h9Qdgu9DFEahQAEIAABCEDgZQjwfvIyPxWNQgACEIAABJ5OYPe9gQ8oT/8JaQACEIAABCAAgVUCuy9Cq3WJgwAEIAABCEDg9QjsvjfwAeX1fnM6hgAEIAABCEDgi8DuixAgIQABCEAAAhD4HAK77w18QPmca4WTQgACEIAABN6OwO6L0NsB4UAQgAAEIAABCDQJ7L438AGliZYNCEAAAhCAAATuTmD3Reju56M/CEAAAhCAAATOEdh9b+ADyrnfgkwQgAAEIAABCDyYgF6EHlyWchCAAAQgAAEIvCABvTcUuzKmolRspRAxEIAABCAAAQhA4DQB3k1OEyUfBCAAAQhA4H0J6L2BDyjv+xtzMghAAAIQgAAEGgT0ItTYZhkCEIAABCAAAQj8I6D3hss/oOwW+tcxDgQgAAEIQAACEDhEQO8nh9KRBgIQgAAEIACBNyag94aHfkBRMS+O/+d/MIAB1wDXANcA1wDXwLXXgL/TRdbai+sz85JjRo/22t/7Sr6931rXkttaL9r3vdpa3C/zMuK6x7ovvdZkFe9z97Xfipf20VZ9Pbou9SAAAQjo/rN6X0z/O1C8EP7rvizw2/HbcQ1wDXANcA1wDXANcA1wDXANcA1wDXANfNI1oE9Hfmatzdj0B5SS1Ivh8weOa4BrgGuAa4BrgGuAa4BrgGuAa4BrgGuAa4Br4BWvgZkPJ9IufUBRcMbWQJY4X9fc85X91ujttWKuXL+inytyXsmA3BCAAAQg8BoERs+Xsu8a92dOuBIXY+K8Vr9o9H+1/ZW1mC/Tx0odxfTy9/YU71a9y2pP82jLvq9JL6s9zaPVvmzZl19snGvN1/8Tff2HYkf70tXy+Z77yimrWNeMfI9txUfNKCf733+GYQELrgGugUdeA+V+PTvaXylCJj9I2GIKAQhAAAIQgAAEnkJA7ydPKU5RCEAAAhCAAASOESjP9KuH3htWa6U73C10NQjyQwACEIAABCDweQT0fvJ5J+fEEIAABCAAAQjMEtB7Ax9QZsmhhwAEIAABCEDg5QnoRejlD8IBIAABCEAAAhC4nIDeG/iAcjlqCkAAAhCAAAQgcDcCuy9CdzsP/UAAAhCAAAQgcB2B3fcG/ic81/02ZIYABCAAAQhA4GICuy9CF7dHeghAAAIQgAAEbkRg972BDyg3+jFpBQIQgAAEIACBOQK7L0Jz1VBDAAIQgAAEIPDKBHbfG/iA8sq/Pr1DAAIQgAAEPpyAXoQ+HAPHhwAEIAABCEAgQUDvDcWujHTUbqGV5oiBAAQgAAEIQAACPQJ6P+lp2IMABCAAAQhAAAKFgN4b+IDC9XBLArpAb9kcTUEAAhCAwEsT0DNm9SXopQ9P8xCAAAQgAAEITBPYfXfgn0CZRM5LWg6YX5jyc5GoIAABCEAAAnkCesbwfM4zQ5knwPWVZ4USAhCAwCsQ2L2v3+oDih+m5V/9o7TqxvWr+3jF/JFRnL/imegZAhCAAATuTcCfNY/s1Ou2/JP9tGq01k/W/uRcke8ui5hP8928xEMAAhCAQI6A7rvFrox0lBdaLdZrMObvzXt5dvd6dX1vt867xTubmv9u5+U8EIAABCBwDwL+zHlER14v6+/2la3jut2axP8l4EyLvzpintZ8NT9xEIAABCCQI+D331zET1X6SbBb6GfZ9szrtPx29N5Oq15ZZ/QJ9NjBr8+OXQhAAAIQWCfgz5/1LLlIrzXr5yrUVbO1ip6xT6DGfSWr5ynxmruvtWIZEIAABCBwHYHd+236Lr1baAWB13R/JVcmxmu4n4lF85OA8ys+AwIQgAAEIHAFAX/eXJFfOb3Oqq9cM3al1kx+tG0CNfZtdX0nmyPq6tlYhQAEIACBXQJ+v13Jlf6b7W6hpeb+/Pn3ld7rF//0iPl9frrWJ+Rzflf8Xp/AkDNCAAIQgMCYgD9vxup1hdfpPdeizucr1T0+66/UycZ4D9mYV9T5Od2fOYvHyW/Fa99tS8s6BCAAAQisE9i9z6a/ROwWWjmi14z+Sr5eTMzv814ce3UCzq/4DAhAAAIQgMAVBPx5c0X+klM15GfqKCbaTKw0Hqu1Z9s79nSaiZ8x+jO1ZmNn9TO9oIUABCAAgb8E/F67wiT9N9vdQkvNdf4JlNLPqeFnq/mn6nxSnsjxk87OWSEAAQhA4DEEHvWs8TozJ/M4+SvxMzFXa3WOYt9xxPP5fObMMS4TuxLzjr8BZ4IABCBwJQG/167UST/9dgstNccHlBVst4jx66X4DAhAAAIQgMBpAo941niN2f49Vn42h/TF3mncta9TjHQ+5dNcVusjK73bUUzZd33xGRCAAAQgcJaA32dXMqfvzLuFlpqzDygl3nuQv5LXY5Sn2KtqeL1P8Z2r2H7K2TknBCAAAQg8joA/b66oupvf4+Vn+pzRZvKd0KinYt9x1M7nazPnfnTcO/4enAkCEIDAFQT8/rySP/0E3C201NzXBxTFeg/ytbdqladYDV/zde1jxwRgOGaEAgIQgAAE9gn482Y/2+8Mu/k9Xv7vKr9XZrS/o8+vqB/Z8xWem1HnKtaHr8c917kfY7JxJcdOrPeADwEIQAACdQJ+n60r+qs/nxId7W6hTurmlmq6QGuyvjfrK0exPnw97rnuDv5V/e3mPcnQc92BOT1AAAIQgMB9CNz9GeH9yR/Rk65mR7FX7N+ljyvOVnLG83md3p7r3I8xZZ4dO7HZGuggAAEIfDIBv8+ucEjf0XcLLTUX/gmUksP7kL+SO+ZSDuV0qz23vt/yXS+/pdW6dC0rXda28sT1TL4Y05vHfD1tbS/G1+a1ONYgAAEIQOBzCMRnwx1PPttj1PfmV563V7e3N+qpF1v2njG8p1jf97L9xZhsnGrHeK1jIQABCEBgn4DfY1eypZ9Uu4WWmrvwA4qfp/gacd33pCm2potrrpcfNXEuXc1GbWZeyxPXMnlcE+Nrc9cXPztiXGaezY0OAhCAAATej4A/J+54Ou+v+KMR9aP5KN/K/qhmb79VrxdT22vlOb3utWu5fb/4mRFjsnHKvRuvPFgIQAACEPhNwO+xv3fHK7knQfhgME57RtE6nK/Ln62ouGJ9+Lp836/50kVb02otajXXfs1KU2xt+L78mk5r0rjVnqzvyddez0or29NqT1q32pP1Pfe1j4UABCAAgc8icPdngfdX/NGI+ux8lHd2X3VLnPxo416rRowrcx+1fa257rSvGrK1/NqTrWl8TTq3vp/xPVZ+Jg4NBCAAAQiMCei+WuzKSEftFlpqLvx/4VEO70W+9rJWccX68HX5vt/ypXXb0mrdtfK1F632ZeO+z6Uptjdc19NGXU+rejFG6z3rMVmdYnp69iAAAQhA4H0J6DlQ7N2G97bTX8zTml95/lgzW2smLmo1z9aa1Sl/sa3hmp5O8VGfiVGs7IkcyoWFAAQgAIGfBPwe+3MnN2s/MUL8bqGQLjXt1fQ9+amk4b9JiTHK5TZqanPXy6/p4pq0snFfc+0XOxoZrWtGOaN2pC/9xZgTPXuOmD/Tk8fjQwACEIDAexDw58HdTuS9nXpOxZxxfhWD1TqzcVFf5lcMr9PL77pML1GfiYn1T+SIOZlDAAIQgMBfAn6PXWGSfirtFlpqrvFPoJRc3o/8bA3pi43D9+RHTW0urduaLq65vvit4bqWRusjre/3atbyeaz2a9Z1oxozWtWKMZprHwsBCEAAAp9BQPf/Yu80vK/TvcXcPr+KgdfInmclpvQf47L1smeP+XtxM9pTvceap8/fOy97EIAABN6dgN9jV86aftvYLbTU3OQHlMwDZnQO35ef6V1atytxrRjPW/zRkL6m055sTeNr0rn1/Zrv2uL3hmt7Ot/zGPddgw8BCEAAAu9PQM+Au51Ufcle0Z9yR/uIWpkaq33FuDI/OTz/KK9rM31EfSYm9nAiR8zJHAIQgAAE/hLwe+wKk/QTabfQUnOdDygln/ckf1RHumJrw/fl13RxTVq3UVObu774rRF1PW3JIX0tn/Zka5q4Jq1s3I9z6WTjvs+lOWE9Lz4EIAABCLw/AT077nRS9SR7ZW+q4faKep6/+KMxq4/5YnymZsxRm3ve2n5cc32mh6jPxIxqruSIOZlDAAIQgMBfAn6fXmEyfgJ+Zd0ttNTcwgeU3kMmcwbXyM/0Lq3blbhWjOd1v6VvrXus/JZ2Z125ZXu5pDlhe3XYgwAEIACB9yPgz467nO7RPXk9+adZKK/sKL90siN93Fec26hZmc/mc33xRyPqMzEx54kcMSdzCEAAAhD4S8DvsStMxk+Cr6y7hZaaG3xAKTm9L/m1WtqTrWlm8sV45XUbNbW564vfGlEX5624uL4aF/OM5tk6PV3c83msr724zhwCEIAABN6bgO7/xd5heD+P6inWvKJurDFiPauP+WL8iTN5zlivNfeYTCsTefsAACAASURBVA9Rn4mJtU/kiDmZQwACEIDAXwJ+j11hkn7b2C201NziB5Tawyrbv+vkZ3qX1u1KXC/Gc7f8XnzZq8WNYlb2Y51WjqyuFc86BCAAAQh8NgF/jjybhHopfRT/UUN13Z6u7bkzZ5vVx35jfKZmzOHzmM/3en4rTusxVutuo2Y091j5oxj2IQABCEAgR0D31WJXRjpqt9BSc3xAqWLz36LlVwP5gNLCwjoEIAABCLwoAT0Hn92++pB9dD+qK3u6vvLK9vJL47anr+15rPyaLrumHFfY2EOsEfdH8934UX72IQABCHwyAb/HrnB4+Q8o5dAOQb7D0Fqxo+Fa+aOYsi+t25W4TEyrnteWH/Np3W3UnJh7/uK3RlbXimcdAhCAAAQ+m4CeI8+koB5kn9GLasue7kF5ZXv5pYm2FxP3YmyZ74xavlNrsa9a3qhpzXdiWzlZhwAEIACBbwJ+n/1ezXvpp9FuoXxL38psTdfJ/87y8+OGr9d8xbut6eKa6+VHTW0urWxN01tTXMvG2Jouak7MY51WzqgrcwYEIAABCEAgS0DPkaz+tE71ZU/nz+ZTfdlsXFanvLKjOOncjmJ83+Pk+/6srxxX2NhLrUbUtOY7sa2crEMAAhCAwDcBv89+r+a99N9WdwvlW/pWztR0rXxlinOt16y0bmu6uOZ6+VFTm0srW9Nk1hRfsx4/2ndtz1eelkb7slld0a+MUZ2VnMRAAAIQgMD9CTzz/q/ass+kpR6KvWJ4/kyNqM/EeN+78Z6r+LV8rTXF9vZ9T3pZ35OvvZGV3u0ohn0IQAACEMgT2L2/pp+yu4XyR/pWztR0bcv/ztz2arFt9ffOqbjvjD+9kj8zRn3U9rO5Vd9zaC1a14zyR+1IH2uVuXLU9liDAAQgAIH3JaD7/8qzY5fKM2vH3q/uxfPLjz34XBq3vj/yPU7+KOaKfdWWzdaQXvbquGx+dBCAAAQ+nYDuy8WujHTUbqGl5hL/ElnP6z1G33U9P8aVeWaciqvVUu7aXm1Netmo0brbqOnNFZfRrGhLzMxQjdm4mRpoIQABCEDgngSe9QzYqXv6eeW9nM6tX322RtTP9hXj1cej7Wofj457NBfqQQACEHhVAn5/XjlD+m+qu4WWmnvzDyjOVH6Nk/aKzQzX12LivuazuXt65ZSd0WZiPJ/0xTIgAAEIQOCzCOgZ8MhTq+bKc2cntnXGK3LGWl4jc+6o1zzmbc2ll23prl5XfdmZeoqRHcVKJzvSsw8BCEAAAvMEdI8tdmWko3YLLTU3+QGl1PA+5c/UVozbTLzr5Y/ipHNbixntxxjXFz+OuO/zqPW562p5V7UlLub2ueeNflYX45hDAAIQgMDrE3jGM2Cn5k5s69fynMW/aqzUiTEz/cXYq841yrvTx2zsrH7UO/sQgAAEIPCbgN9rf++OV9JP2t1C41Z+K1Zqeoz835nbK4px21b/3PEY+T8Vf2faq9mMvqbxtZjX9+RHTW1etLX1stYbtZievlfHc410oxrsQwACEIDAexGIz4irT7daz+Pk13rVnmxN42vSyfreaV81ZGP+2rrWoo2xcT6rj/En5zu9zMTOaE+ej1wQgAAEPo2A329Xzt7/m7Bl3C1kqdLuSk2PKf7siPEzOWqxim/t9dbVe02jvZqN+pqmrEVddt7Kp/VaHu31bC0uu9bLyx4EIAABCLwnAX9GXH1Cr3XCj/32cma0UXN6Xuuv1IjrsW7cL/PR8JiR9up97yXTe+wnG5/VxfzMIQABCEBgjoDfb+ci/6rHT7GvrLuFZpvzevJLjuKPhvSzWo/r+a36vZi4p7PEdZ+rjq9FX5pWPt+v+TFfZr6bpxavtUz9qFEsFgIQgAAEPo+AnglXnVz5T9pWr6s1WvlOrmd7q9Vsxbq2pvH9Z/ineop5/Cxxr8wZEIAABCBwHQG/765USd+ldwuNmvP8Gb+Xz+OzOo/J+K28M7E9refv6Xp7nqPn93L43okcynciVy8HexCAAAQg8P4E9Ewp9qrhNU75rV5n87fyXLU+6q9XdxQb93u5HrUXeyrzlVHL01pbyU8MBCAAAQjkCfj9Nx/1rUw/CXYLfZese54/49ez/F1VfE9T9qRbsb3cvXweF3W+577ryrrPW77HZ/xWnrKeGb342t5uzkw8GghAAAIQeG8C/ny56qRe44Q/6jNbY5Tnqv1Wf9l6rXhfz+Z6hM77Kv7qiHlq89XcxEEAAhCAQJ6A33/zUd/K9JNgt9B3STwIQAACEIAABCCwT4B3k32GZHgOgXLtMiAAAQhA4PEEdt8d0nfv3UKPR0NFCEAAAhCAAATenYDeT979nJwPAhCAAAQgAIF9AnpvKHZlpKN2C600RwwEIAABCEAAAhBoEeDdpEWGdQhAAAIQgAAEagR23x34gFKjyhoEIAABCEAAArcnsPsSdPsD0iAEIAABCEAAAkcJ7L478AHl6M9BMghAAAIQgAAEHkVg9yXoUX1SBwIQgAAEIACBexDYfXfgA8o9fke6gAAEIAABCEBgksDuS9BkOeQQgAAEIAABCLw4gd13Bz6gvPgFQPsQgAAEIACBTyXgL0HFZ0AAAhCAAAQgAIEeAX936Olae+m3jd1CrQZYhwAEIAABCEAAAisEeDdZoUYMBCAAAQhA4HMJ7L478AHlc68dTg4BCEAAAhB4aQK7L0EvfXiahwAEIAABCEBgmsDuuwMfUKaREwABCEAAAhCAwF0I6EXoLv3QBwQgAAEIQAAC9yWg94ZiV0Y6arfQSnPEQAACEIAABCAAgR4BvZ/0NOxBAAIQgAAEIACBQkDvDW/zASVzkIxm9/JwsNEXeLeqJ63mbrVXbBmy0vi++7X92prni/E+l68eWnOtZ+1KvhiTrYXuz48//PCAxyP/LK3UKjEafr1qTVZ7o7nr5N/RlnPEvuLZMvuKqeVTvDSaR+ux0rbWfH/Gr9VsxUsb98t6dihHK0b7tXy9vZq+rLXqtPRxfVRzJf8oRvuy6inOfd335MtKV7OukV+s/BjTWo+6u8/f5Rx350x/EIAABE4Q0HNp9d6dfkvxQvj8ZY1rgGuAa4BrgGuAa4BrgGuAa4BrgGuAa4BrgGvgztdA/Ojivca9zJwPKF//zYiDxOcmwDXANcA1wDXANcA1wDXANcA1wDXANcA1wDXwvtdA5oNJ1Cx/QFGid7ugyrmuPNNs/l39bLzO3oprrZe42t5/i/YfvfyKN/kPVzV8Uflqe67L+DGX54x7mXy7GtWs5fHeavuZNeU/kStTD815Aq3fUOsnKyqnbCa3tLIlRr6s8tTm2pOVJtqYt6WPOuWRvmWlq1mP8f3WetFouL63HnVx3sonXdl3f6SX9s62dQat+5nvfA56e9+Xcn5bfluuAa4BroH7XwP+3pD1v9/kBhF+AQykbEMAAhCAAAQgAIGHEND7yUOKUQQCEIAABCAAgSME9Px2GxP39opW+4orcx9x32Oi1uN6/s8KHaWKrxbqpGYLAhCAAAQgAAEILBHQ+8lSMEEQgAAEIAABCHwUAb03rH7X4APKR10uHBYCEIAABCDwPgR2X4LehwQngQAEIAABCEAgQ2D33YEPKBnKaCAAAQhAAAIQuB2B3Zeg2x2IhiAAAQhAAAIQuJTA7rsDH1Au/XlIDgEIQAACEIDAVQR2X4Ku6ou8EIAABCAAAQjck8DuuwMfUO75u9IVBCAAAQhAAAIDArsvQYP0bEMAAhCAAAQg8GYEdt8d+IDyZhcEx4EABCAAAQh8CoHdl6BP4cQ5IQABCEAAAhD4S2D33YEPKFxJEIAABCAAAQi8JIHdl6CXPDRNQwACEIAABCCwTGD33YEPKMvoCYQABCAAAQhA4JkEdl+Cntk7tSEAAQhAAAIQeDyB3XcHPqA8/jejIgQgAAEIQAACBwjsvgQdaIEUEIAABCAAAQi8EIHdd4dbfUDxw+z4L/T70SoEIAABCEAAAosE/F1hMcVymNeu+cuJG4G1GmXt1GjlP1njVK/vkMd5nz6P53b/ZB3PK/9k/pJLed2erkE+CEDg8wjs3lPST97dQqOfxvOf9ke12YcABCAAAQhA4PUI+PvCI7v3uiN/t69Rfu3v1FGOkd2pQew3gcj5e2fPi3lb850qrZy+vpO/xHqulr9bg3gIQOBzCfh9ZYXCbT6glOb9MFf4K4CIgQAEIAABCEDgngT8XeERHXq9WX+lv2fWKP3W6q+cg5g6S/E9wUe5ii3D5zV/pWYtT2/t6hor+YmBAAQg4PetFRq3+oDiB/CDue+amu/aml+LYQ0CEIAABCAAgdcj4M/5q7v3Wiv+bH+xRis+q6vFZ2OzulqNT1+L7GrzXUYxZy1fRlOL01qML/M4oibuj+aZ+IxmVId9CEDgswn4fWSFxO+7XyPLbqFG2u6y15TfDbBN6WvWZLgQgAAEIAABCLwoAX/GX30EryW/VlN7NVvT19ZibE2jtRmtYoqdiZvReo1VX/VW4+8Sp3OM7G6/MX8rX1YX42NcmbdG1LZ0cT3GXVEj1mQOAQh8JgG/36wQaN8BQ7bdQiFdauo15acCv0SKiXYmB1oIQAACEIAABO5JwJ/vV3bodYo/GlGv+ShO+9LLar1lpZNt6XxdWlnfq/nSydY0p9YeUeNUr6M88Syaux3l6O17nuKPxqy+5JuNcf2oH+17TPF7I2pH+l4u9iAAgc8j4PeQldP371CWcbeQpUq7XlN+OvhLqLhoZ/OghwAEIAABCEDgXgT82X5lZyt1PMb9UZ+uLX5mzMbM6ksPMSbbW6b/qFGtuP4uc51PdvVcinc7yuXa4o/GrL7km42Z1a/UGJ2TfQhA4HMI+D1n5dTjO+dX1t1CS839+TN9E451vG/3o445BCAAAQhAAAKvReARz/WdGh4rf0RYOtmRvuxL67YX57riZ8dqXDZ/0XmNmbhX0voZZ/jHM67kiTGj+lEfe6jNZ2Nm9aVmjClzBgQgAIEMAb9/ZPRRk77b7BaKhTNzryk/E+caxUXrGnwIQAACEIAABF6PgD/br+p+p4bHyu/1KY3bnt73PKb4rRF1PW3MEWPj/om51ziR7445/Iwz/ONZVvPMxM1o1V+MKfPWmNHGHDE27jOHAAQgUCPg947a/mitfUcLkbuFQrrU1GvKTwWaSHHRmgQXAhCAAAQgAIEXJODP9qvaV42V/Ip128vjOvk9ve9JL+t77mvfre/3fI+R39PP7imn7Gz8q+h1PtmVvhXrNpvHY4rfGlHX08YcMTbua57VSe82xs7053nwIQCBzyLg946Vk7fvmiHbbqGQLjX1mvJTgSZSXLQmwYUABCAAAQhA4AUJ+LP9ru17j8XvjRltzBNjW7Wyupi/zHdia/l87crcXucOfjzrSk8xR5lnRzY2q6vVjbE1TVnL6mrxMbbMGRCAAARGBPzeMdLW9tN3mt1CteKjNa8pfxQT9xUXbdRl5jFHnGdy1DSeJ+77nvtRl517jpqfzYMOAhCAAAQg8GwC/hx7di+t+t5j8Vsj6nraWo5sfNTVcvXWduM9d8yVmXt8yx/lacU9Yr3W20rdnTzZ2Kib6TPGlnkcGU2M8fluvOfChwAEPoeA3ztWTv37btbIsluokba77DXldwMqm4qLtiJtLsXY0byZ6GujF++xPV3ZmxmjXHF/JjdaCEAAAhCAwDMI+LPrGfUzNb3H4rdG1PW0tRyZ+IymltvXTuQo+Wp5MmveS/Qz8a6J8Y+Ye335K3UVKzuTQzFuY7zvyY+a3lwxbqPe9+RHTW+uGLc9PXsQgAAECoHde0b7SR747hYK6VJTryk/FWgixUVrkq4b47LzVtJRvOJGOu1L37PSRlti4prPeznZgwAEIAABCDybwN2fWd5f8Xsjakf6mCsTn9HEvHF+Iody1nKN1hQb7SiutR/zPGIee5mtGePLfGZk4jOaXs1MfEbTq1H2TuQY1WAfAhB4LwJ+31g5WfqOu1toqbnN/zfG3rP72V48pvi9EbUtfU23s9brqezF3C191Gne0rMOAQhAAAIQeDYBPauKvePw/kY9Ru1IH8+bic9oYt44P5Ej5tQ85tb6yM7ERW2ZP3rEHmbrx/iVM8QcsYe4P1sjE5/RxL7i/ESOmJM5BCDw3gT8vrFy0vRTY7fQUnOLH1C81+hn+1iJizFl3hs1va/VYn1ffk2nNWmKHQ3Xyh/FsA8BCEAAAhB4FgE9qzLPuGf0ONOfa+XP9qw42RivdbdRM5p7rPxRTHZf+WQzcdLKjmKki3YUd2o/1i3z2bGboxYf+6hpZvrMxGc0o5oncoxqsA8BCLwXAb9vrJwsfdfeLbTUXOUDivcx48/Ur+XNxM/G1fRlrTdqMS29a1uauO4x8qOGOQQgAAEIQOAOBPScKvZuw3vL9Bf1mZh45phjtH+ixkqO2Jfmo/6lc3si5uQZvLeaH/tdqb2boxYf+6hpaudprWXiM5pW/rJeiy9rDAhAAAI9An7v6Olae+m7zG6hVgO9da+56/fqxL1YK+735jG2zFtjRus5Ypzvue86X+/5HiO/p2cPAhCAAAQg8CwCek4Ve7cx25vr5c+cSTFuY7zvyY+a0Vxxbkcx2X3PWfzRmNV7vhibqefxq/6Jurs5MvEZTY9BJj6j6dUoeydyjGqwDwEIvBcBv2+snGz8dPrKultoqbnKP4HSyuP99fxWvK/HeN8b+TG2zFtjRus5YpzvyY8a9dFaL3G1Pa0pLxYCEIAABCBwFwJ6RhV7p+F9ZXuLMdk4nTsTn9EoX8ueyJHN3dJpPfai9YyNsWX+iHGi7m6OTHxG0+OVic9oejXK3okcoxrsQwAC70XA7xsrJ0s/LXYLLTU38QEl5vd+a37Uaz6jVYzbWnxZq42atqaLazEu7pd51OzOazVYgwAEIAABCDyTgD/bntlHrL3Sl8fIj3l7c8W4jXrfkx81o7ni3I5isvues/ijMav3fDFWc9dc4auO29k6Hit/Jodi3MZ435MfNb25YtxGve/Jj5rRXHFuRzHsQwACn01g934xfjp98d0ttPIzeU35s3kUV7O1XFldLVZr2RxZnfLKxjitu42anbnnxYcABCAAAQjchYA/2169Jz+L/JkzKcZtjPc9+VEzmivO7Sgmu+85i98bUTvSx1y1+NkcMWdmXqubiXPNbo5MfEbjPUU/E5/RxLxxfiJHzMkcAhB4bwJ+31g5af/pZBl3C1mqtOs15aeDTajYaE3yz42aMp8d2RxZXawf4+J+mWc0tTjWIAABCEAAAq9CwJ91d+jZ+yn+zIixV8Tv1ijniTlmzjjSzuSO2jKfHSdyPKPmbt+Z+Iymd/ZMfEYzW6PkZEAAAhDoEfB7T0/X2kvfZXYLtRrorXtN+T19a0+xNRtjMpoYE+fZHFndKH/cL/OYu6ZhDQIQgAAEIPDKBPxZ9+xzqJfSR/Fnh+LdzuTwOPkxXutuo6Y39zj5Pf3snnLK9uKlcdvT1/Y8Vn5Nd3JNddzO5vdY+TM5FOM2xvue/KjpzRXjNup9T37U9OaKcdvTswcBCECgENi9Z6Sf8LuFVn4uryl/JU8EpVzFxuF78qNmNFec21qM78uv6eKatLJxv8y1J1vTsAYBCEAAAhB4ZQJ6xhX7zKE+Sg+rvSiH25kzeZz8GK91t1HTm3uc/J5+dk85ZXvx0rjt6Wt7Hiu/pju5pjpuV/J7fPFnRoytxWc0vZqZ+Ixmt0Yvnj0IQOAzCfi9Z4VA+o67W2ipuY1/iWys5/27n9FFzWju+eXXYrTntqaLa64vfm1kNLU41iAAAQhAAAKvQsCfdc/q2XtoPZOzve3kirGtXqIu21vRxdhWjZmcro35fS/6UbvSy4kcsa/R/FTNmGdU1/djbJnXRtTVNK21TGzUtPrI1piNb+VlHQIQeG8Cfu9ZOWn9jlnJtFuoknK45DXlD4MaAsVHG+Vxv8xnRzZHVhfrx7i4X+ZRU+YrYzVupRYxEIAABCAAgRkC/qybiTul9fonnpcx30zOGNs6Y9Tt1JiJbfXj67E334t+1K70ciJH7Gs0P1VzJ082NqurnTnG1jRlLaurxcfYMmdAAAIQGBHwe8dIW9tP32l2C9WKj9a8pvxRTGtf8dFGfdwv89kRc7Tioy5bK8bV8kdNNrfn8hy+jg8BCEAAAhC4A4FnPqdUu3BYecbW+Cmn25qutuYxvX6irqeNdWJs3N+dz+SPWs1nelCM7Ezsqla13K7k8nj52TzSy7bitO+2pfV118v3ffe1L+t7I18xsiM9+xCAAAQKAd0zil0Z6ajdQkvN3eR/wjML11n1YqOup3V+Mc735EeN5trP2JWYTF40EIAABCAAgRME9Jwq9tHjqtqeN3uuGDOKi/osu9W4q/LHfkbnjn3E+Lh/xTzWnO3Ze4q5fK/lx5hR/ahv5fX1mZioHfWjOqtxisdCAAKfS8DvHysU0m8cu4WWmnvCB5TSp59V/kz/ipFtxWrfbUvr664vfmtEXU8bc3hs3GMOAQhAAAIQuAOBZz2rduoqtsVP+7Itna9LK+t7NV862ZomrkkrG/d358rrdpTTtfJHMdqX3q32rrReT/5qPcXLZvJIKzuKkU52pC/70sqOYqSTHelXamRyooEABD6DgO41xa6MdNRuoaXmDn1A8d7db/XkGvdbel93ffF7I2pHeuWKcVqPNuo0j7o4l0427jOHAAQgAAEI3IGAnlPFPmrs1MzGui5ztll9YTUbM6uf/T1i/jIfjVpMJm7l/KNesvu1nrOxNV3MV9P42tX6UuuONZwBPgQg8NkE/B61QmL8dPrKultoqbknfUApvfp55WfOIK1sL0Yatz299lxf/NaIutrcY0f7rsWHAAQgAAEIPJuAP7ce0ctqPY+T3+tXGtlTWs+j3LK+F31pZOP+qbnyy8a8tXWtuY1xtbnri/+oEevu1o75eueY0XqembgZ7akapSYDAhCAQJaA36eyMa5L33F2C3nRrO815Wdji04xNTvKU4spa70RY3rashf1o/zKF+O0XrNROzuv5WQNAhCAAAQgcAcC/ky7uh+vtetneo01ajFRU+YzI8a3YrO6Vnx2Pdbx8/ie5/N1910TfdcV/5Ej1j5RP+ZsnSeri/HZuKgr8+yIsb24GW0vD3sQgMBnEvB7yAqB9J1tt9Bsc16v55e8vf3aXraXWqzWYg5fL/5oSF+zvdiaflSvFTNa7/XBHgQgAAEIQODZBPw5dlUvXuOUn+m1VsvjRvuubfmjHKP9Vt7V9Vq92lrMX9NozbVac+v7j/C9tvs7tT1P8ePw/bJX08SYOPcctfi4X9PEnHEec8T9Mo+alTq1vKxBAAKfQ8DvIyun/n2XbWTZLdRI+2/Z81/l/ys24az00ks/k095ZmKKtjZmctTiWYMABCAAAQjcjYA/267qzWuc8rO9ztTL5oy6R9SINXvzUT+t2FFcbb+V6+R6re5obaX+KKfvr+QvMZ5j5N+5xmpvxEEAAu9BwO9fKyeq/227kmm3UCXljyXPf9L/UWRjku1pVCKbp+g0ZmI8TvFuR7lciw8BCEAAAhC4MwF/pl3Vp9c44c/2Oao5m6+lf1SdVn2tt/rQ/si24uP6KM+p/Vg3M1+tfWVu9fQuNXQeLAQg8HkE/D62cvrvv6UPoncLDdL/t11qvNJ4tX5fiS29QgACEIAABEYEHvFuMurhHfff5f3Gr49nn2lUf7Q/c5094tzPqHGS0QxPtBCAwHsR8PvXysnSXyx2C600RwwEIAABCEAAAhBoEeDdpEWGdQhAAAIQgAAEagR23x34gFKjyhoEIAABCEAAArcnsPsSdPsD0iAEIAABCEAAAkcJ7L478AHl6M9BMghAAAIQgAAEHkVg9yXoUX1SBwIQgAAEIACBexDYfXfgA8o9fke6gAAEIAABCEBgksDuS9BkOeQQgAAEIAABCLw4gd13Bz6gvPgFQPsQgAAEIACBTyWw+xL0qdw4NwQgAAEIQOBTCey+O/AB5VOvHM4NAQhAAAIQeHECuy9BL3582ocABCAAAQhAYJLA7rsDH1AmgSOHAAQgAAEIQOAeBHZfgu5xCrqAAAQgAAEIQOBRBHbfHfiA8qhfijoQgAAEIAABCBwlsPsSdLQZkkEAAhCAAAQgcHsCu+8OfEC5/U9MgxCAAAQgAAEI1AjsvgTVcrIGAQhAAAIQgMD7Eth9d+ADyvteG5wMAhCAAAQg8NYEdl+C3hoOh4MABCAAAQhA4BeB3XcHPqD8QsrCKxEofwCyw/+wzMQpf4z3HL7X07su+iUurq3OW7l8vdan75+uvZrv1eKyDLM6P3/rN/NcNU3J0RuqoTw9re+N8roW//kE/Hd+fjc/O1i9lnQm2Z9Zz81W+zvXAZkgAAEIQAACEDhBQO8Mq8/2/lu1deiF8P8c+4suLGHJNcA1wDXANcA1wDXANcA1wDXANcA1wDXwCteAPhHEXlvrUae59D1b02pNtsTLd9vKm9G0Yv+r1dv0PS+Ezx9urgGuAa4BrgGuAa4BrgGuAa4BrgGuAa4BrgGugVe9Bvx7R9bnn0D5wwX/qhd8tu/yh0Ha+AdD66vW88Uc2tN6nGtdVvvFas2t72f8EqtcLb3yt/Yz660cXr+nkS5TK6tp1SvxvXqK62nUg2s8LvotfW1da8V6Hq1rrTePeyWmDMXu2pO5dnshnudXuQbiNRnnXCdcJ1wDXANcA1wDXANcA7Vr4L+XiMn/SH9AmcyLHAIQgAAEIAABCEAAAhCAAAQgAIEXJhA/POgovq41Wd9zf7RftGV4TMuv6WL+mkb5pJ21fECZJYYeAhCAAAQgAAEIQAACEIAABCAAgY8jwAeUj/vJOTAEIAABCEAAAhCAAAQgAAEIQAACswT4gDJLDD0EIAABCEAAAhCAAAQgAAEIQAACH0eADygf95NzYAhAAAIQgAAEIAABCEAAAhCAAARmCfABZZYYeghAicz0NgAAAMJJREFUAAIQgAAEIAABCEAAAhCAAAQ+jgAfUD7uJ+fAEIAABCAAAQhAAAIQgAAEIAABCMwS4APKLDH0EIAABCAAAQhAAAIQgAAEIAABCHwcAT6gfNxPzoEhAAEIQAACEIAABCAAAQhAAAIQmCXAB5RZYughAAEIQAACEIAABCAAAQhAAAIQ+DgCfED5uJ+cA0MAAhCAAAQgAAEIQAACEIAABCAwS4APKLPE0EMAAhCAAAQgAAEIQAACEIAABCDwcQT+DzidOIEIqgOlAAAAAElFTkSuQmCC"}}},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:30.794385Z","iopub.status.busy":"2020-09-16T09:09:30.792972Z","iopub.status.idle":"2020-09-16T09:09:30.795425Z","shell.execute_reply":"2020-09-16T09:09:30.795941Z"},"papermill":{"duration":0.059671,"end_time":"2020-09-16T09:09:30.796051","exception":false,"start_time":"2020-09-16T09:09:30.736380","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"def crop_image(img: np.ndarray):\n    edge_pixel_value = img[0, 0]\n    mask = img != edge_pixel_value\n    return img[np.ix_(mask.any(1),mask.any(0))]\n\ndef resize_image(img: np.ndarray,reshape=(512,512)):\n    img = cv2.resize(img,(512,512))\n    return img\n\ndef preprocess_img(img,resize_type):\n    if resize_type == 'resize':\n        img = [resize_image(im) for im in img]\n    if resize_type == 'crop':\n        img = [crop_image(im) for im in img]\n        \n    return np.array(img, dtype=np.int64)","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:30.899640Z","iopub.status.busy":"2020-09-16T09:09:30.898917Z","iopub.status.idle":"2020-09-16T09:09:30.902764Z","shell.execute_reply":"2020-09-16T09:09:30.902288Z"},"papermill":{"duration":0.06064,"end_time":"2020-09-16T09:09:30.902860","exception":false,"start_time":"2020-09-16T09:09:30.842220","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"class Test_Generate(Dataset):\n    def __init__(self,imgs_dicom,resize_type='no'):\n        self.imgs_dicom = imgs_dicom\n        self.resize_type = resize_type\n        \n    def __getitem__(self,index):\n        \n        slice_img = self.imgs_dicom[index].pixel_array\n        slice_img = (slice_img-slice_img.min())/(slice_img.max()-slice_img.min())\n        slice_img = (slice_img*255).astype(np.uint8)\n        if self.resize_type == 'crop':\n            slice_img = crop_image(slice_img)\n        elif self.resize_type == 'resize':\n            slice_img = cv2.resize(slice_img,(512,512))\n            \n        slice_img = slice_img[None,:,:]\n        slice_img = (slice_img/255).astype(np.float32)\n        return slice_img\n        \n    def __len__(self):\n        return len(self.imgs_dicom)","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:09:31.103500Z","iopub.status.busy":"2020-09-16T09:09:31.102595Z","iopub.status.idle":"2020-09-16T09:09:36.380909Z","shell.execute_reply":"2020-09-16T09:09:36.380153Z"},"papermill":{"duration":5.335577,"end_time":"2020-09-16T09:09:36.381027","exception":false,"start_time":"2020-09-16T09:09:31.045450","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"#the model has been trained from [this](https://www.kaggle.com/hfutybx/unet-densenet121-lung-of-segmentation)\ndevice =  torch.device('cuda:0')\nmodel = smp.Unet('densenet121', classes=1, in_channels=1,activation='sigmoid',encoder_weights=None).to(device)\nmodel.load_state_dict(torch.load(f'/kaggle/input/2020osic/best_lung_Unet_densenet121.pth'))\nbatch = 8\n\ndef Unet_mask(model: nn.Module, dataloader: DataLoader):\n    model.eval()\n    outs = []\n    for idx, sample in enumerate(dataloader):\n        image = sample\n        image = image.to(device)\n        with torch.no_grad():\n            out = model(image)\n        out = out.cpu().data.numpy()\n        out = np.where(out>0.5,1,0)\n        out = np.squeeze(out,axis=1)\n        outs.append(out)\n\n    outs = np.concatenate(outs)\n    return outs","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:10:22.961852Z","iopub.status.busy":"2020-09-16T09:10:22.961067Z","iopub.status.idle":"2020-09-16T09:10:22.964951Z","shell.execute_reply":"2020-09-16T09:10:22.964481Z"},"papermill":{"duration":0.064921,"end_time":"2020-09-16T09:10:22.965048","exception":false,"start_time":"2020-09-16T09:10:22.900127","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"def caculate_lung_volume(patient_scans,patient_masks):\n    \"\"\"\n    caculate volume of lung from mask\n    Parameters: list dicom scans,list patient CT Mask\n    Returns: volume cm³　(float)\n    \"\"\"\n    lung_volume = 0\n    for i in range(len(patient_masks)):\n        \n        pixel_spacing = patient_scans[i].PixelSpacing\n        slice_thickness = patient_scans[i].SliceThickness\n        lung_volume += np.count_nonzero(patient_masks[i])*pixel_spacing[0]*pixel_spacing[1]*slice_thickness\n        \n    return lung_volume*0.001","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:10:23.200485Z","iopub.status.busy":"2020-09-16T09:10:23.198811Z","iopub.status.idle":"2020-09-16T09:10:23.201292Z","shell.execute_reply":"2020-09-16T09:10:23.201782Z"},"papermill":{"duration":0.067317,"end_time":"2020-09-16T09:10:23.201898","exception":false,"start_time":"2020-09-16T09:10:23.134581","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"def caculate_histgram_statistical(patient_images,patient_masks,thresh = [-600,0]):\n    \"\"\"\n    caculate hisgram kurthosis of lung hounsfield\n    Parameters: list patient CT image 512*512,thresh divide lung\n    Returns: histgram statistical characteristic(Mean,Skew,Kurthosis)\n    \"\"\"\n    statistical_characteristic = dict(Mean=0,Skew=0,Kurthosis=0)\n    num_slices = len(patient_images)\n    \n    #patient_images = patient_images[int(num_slices*0.1):int(num_slices*0.9)]\n    #patient_masks = patient_masks[int(num_slices*0.1):int(num_slices*0.9)]\n    patient_images = patient_masks*patient_images\n    patient_images_nonzero = patient_images[np.nonzero(patient_images)]\n    \n    s_pixel = patient_images_nonzero.flatten()\n    s_pixel = s_pixel[np.where((s_pixel>thresh[0])&(s_pixel<thresh[1]))]\n    \n    statistical_characteristic['Mean'] = np.mean(s_pixel)\n    statistical_characteristic['Skew'] = skew(s_pixel)\n    statistical_characteristic['Kurthosis'] = kurtosis(s_pixel)\n    \n    return statistical_characteristic","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def extract_ct_features(dicom_path: str, dicom_df: pd.DataFrame)->pd.DataFrame:\n    lung_stat_pd = pd.DataFrame(columns=['Patient','Volume','Mean','Skew','Kurthosis'])\n\n    for i in tqdm(range(len(dicom_df))):\n        path = os.path.join(dicom_path, dicom_df.iloc[i].Patient)\n        lung_stat_pd.loc[i,'Patient'] = dicom_df.iloc[i].Patient\n        patient_scans = load_scan(path)\n\n        ds = Test_Generate(patient_scans, dicom_df.iloc[i].resize_type)\n        loader = DataLoader(ds, batch_size=batch, shuffle=False, num_workers=4)\n        masks = Unet_mask(model, loader)\n\n\n        patient_images = transform_to_hu(patient_scans)\n        patient_images = preprocess_img(patient_images, dicom_df.loc[i,'resize_type'])\n\n        lung_stat_pd.loc[i,'Volume'] = caculate_lung_volume(patient_scans,masks)                           \n        #patient_images = resize_image(patient_images) if dicom_pd.iloc[i].resize_type=='resize' else patient_images\n        #patient_images = resize_image(patient_masks) if dicom_pd.iloc[i].resize_type=='resize' else patient_images\n\n        statistical_characteristic = caculate_histgram_statistical(patient_images,masks)\n        lung_stat_pd.loc[i,'Mean'] = statistical_characteristic['Mean']\n        lung_stat_pd.loc[i,'Skew'] = statistical_characteristic['Skew']\n        lung_stat_pd.loc[i,'Kurthosis'] = statistical_characteristic['Kurthosis']\n        \n    return lung_stat_pd","execution_count":null,"outputs":[]},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:26:47.011607Z","iopub.status.busy":"2020-09-16T09:26:47.009764Z","iopub.status.idle":"2020-09-16T09:26:47.039307Z","shell.execute_reply":"2020-09-16T09:26:47.039801Z"},"papermill":{"duration":0.093492,"end_time":"2020-09-16T09:26:47.039938","exception":false,"start_time":"2020-09-16T09:26:46.946446","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"train_feature = extract_ct_features(train_dicom_path, train_dicom_df)\ntrain_dicom_feature = pd.merge(train_dicom_df, train_feature, on=['Patient'])\ntrain_dicom_feature = train_dicom_feature.drop(['list_dicom', 'height','width','resize_type','n_dicom'], axis=1)\n#train_dicom_feature.head()","execution_count":null,"outputs":[]},{"metadata":{"papermill":{"duration":0.059543,"end_time":"2020-09-16T09:26:47.160248","exception":false,"start_time":"2020-09-16T09:26:47.100705","status":"completed"},"tags":[]},"cell_type":"markdown","source":"Volume is 3000~4000ml is normal"},{"metadata":{"execution":{"iopub.execute_input":"2020-09-16T09:26:47.429027Z","iopub.status.busy":"2020-09-16T09:26:47.428446Z","iopub.status.idle":"2020-09-16T09:26:47.690586Z","shell.execute_reply":"2020-09-16T09:26:47.689518Z"},"papermill":{"duration":0.323332,"end_time":"2020-09-16T09:26:47.690712","exception":false,"start_time":"2020-09-16T09:26:47.367380","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"train_dicom_feature.to_csv(f'{WORKING}/CT_feature-train.csv',index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":false},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#merge/not merge with tabular data\nif 0:\n    test_df = pd.read_csv(f'{INPUT}/test.csv')\n    temp_df = pd.DataFrame(columns=test_df.columns)\n    for i in range(len(test_dicom_df)):\n        patient_df = test_df[test_df.Patient==test_dicom_df.iloc[i].Patient]\n        zeroweek = patient_df['Weeks'].min()\n        #if sum(patient_pd.Weeks==zeroweek)>1:\n        #    print(pd.unique(patient_pd.Patient))\n        temp_df = temp_df.append(patient_df[patient_df.Weeks==zeroweek].iloc[0])\n    test_dicom_df = pd.merge(test_dicom_df, temp_df, on=['Patient'])\n    test_dicom_df.head()","execution_count":null,"outputs":[]},{"metadata":{"papermill":{"duration":0.056385,"end_time":"2020-09-16T09:26:48.040108","exception":false,"start_time":"2020-09-16T09:26:47.983723","status":"completed"},"tags":[],"trusted":true},"cell_type":"code","source":"test_feature = extract_ct_features(test_dicom_path, test_dicom_df)\ntest_dicom_feature = pd.merge(test_dicom_df, test_feature, on=['Patient'])\ntest_dicom_feature = test_dicom_feature.drop(['list_dicom', 'height','width','resize_type','n_dicom'], axis=1)\n#test_dicom_feature.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_dicom_feature.to_csv(f'{WORKING}/CT_feature-test.csv',index=False)","execution_count":null,"outputs":[]},{"metadata":{"papermill":{"duration":0.057107,"end_time":"2020-09-16T09:26:47.804753","exception":false,"start_time":"2020-09-16T09:26:47.747646","status":"completed"},"tags":[]},"cell_type":"markdown","source":"# References:\n1. https://www.kaggle.com/hfutybx/osic-feature-extract-from-ct\n2. https://www.kaggle.com/allunia/pulmonary-fibrosis-dicom-preprocessing\n3. https://www.kaggle.com/aadhavvignesh/lung-segmentation-by-marker-controlled-watershed\n4. https://www.kaggle.com/currypurin/osic-image-shape-eda-and-preprocess\n5. https://kaggle.com/kugane/lazy-lung-cropping/notebook"},{"metadata":{"trusted":false},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}