{ "cells": [ { "cell_type": "code", "execution_count": 65, "metadata": {}, "outputs": [], "source": [ "import sys\n", "import os\n", "import numpy\n", "import matplotlib.pyplot\n", "import chardet\n", "import json\n", "import importlib\n", "import cModel\n", "importlib.reload(cModel)\n", "import runSolver \n", "importlib.reload(runSolver)\n", "import propagateErrorLN\n", "importlib.reload(propagateErrorLN)\n", "#you should get nixSuite via git clone https://git0.fmf.uni-lj.si/studen/nixSuite.git\"\n", "#if you don't put it to $HOME/software/src/, you should update the path\"\n", "#load a solution\n", "def downloadKey(key, doDownload=True):\n", " fh=os.path.expanduser('~')\n", " setupFileSrc=os.path.join(fh,'software','src','PBPK','setup','setupFast.json')\n", " setup=runSolver.parseSetup(setupFileSrc)\n", " ref=runSolver.getRunRef(setup,key)\n", " localDir=os.path.join(fh,'temp',ref)\n", " setup['localDir']=localDir\n", " setup['startFromRef']=ref\n", " if doDownload:\n", " runSolver.loadSolutionFromRef(setup,True)\n", " return {'localDir':localDir,'ref':ref,'key':key}\n", "\n", "def connectDB(host):\n", " nixSuite=os.path.join(os.path.expanduser('~'),'software','src','nixSuite')\n", " sys.path.append(os.path.join(nixSuite,'wrapper'))\n", " import nixWrapper\n", " nixWrapper.loadLibrary('labkeyInterface',1)\n", " import labkeyInterface\n", " import labkeyDatabaseBrowser\n", " import labkeyFileBrowser\n", " importlib.reload(labkeyDatabaseBrowser)\n", " \n", " #check connectivity. This checks the configuration in $HOME/.labkey/network.json, \n", " #where paths to certificates are stored\n", " net=labkeyInterface.labkeyInterface()\n", " fconfig=os.path.join(os.path.expanduser('~'),'.labkey','{}.json'.format(host))\n", " net.init(fconfig)\n", " net.getCSRF()\n", " db=labkeyDatabaseBrowser.labkeyDB(net)\n", " return db\n", " \n", "def getWeeklyPortion(x):\n", " if x==None:\n", " return 0\n", " if x==1:\n", " #never\n", " return 0\n", " if x==2:\n", " #less than 1 /month\n", " return 6/52.\n", " if x==3:\n", " #1-3 per month\n", " return 2/4.345\n", " if x==4:\n", " #once per week\n", " return 1\n", " if x==5:\n", " #twice per week\n", " return 2\n", " if x==6:\n", " #three times per week\n", " return 3\n", " if x==7:\n", " #four times per week\n", " return 4\n", " if x==8:\n", " #five times per week\n", " return 5\n", " if x==9:\n", " #six times per week\n", " return 6\n", " if x==10:\n", " #at least once per day\n", " return 10.5\n", " return 0\n", "\n", "def getDataScaled(variable,q0):\n", " db=connectDB('merlin')\n", " project='PBPK/PHIME'\n", " schema='study'\n", " query='data'\n", "\n", " ds=db.selectRows(project,schema,query,[])\n", " thg=[]\n", " \n", " #HG is hair growth in g/week, to match concentration measurement of ng/g\n", " #normally given as 2g/month,\n", " #from\n", " #strand weight: w=0.3 mg /15 cm/6 inch length/strand\n", " w=0.3e-3/6\n", " #growth of G=0.5 inch/month\n", " G=0.5\n", " #total count of C=0.8-1.2e5 strands on a head\n", " C=1e5\n", " #G*w*C\n", " HG=G*w*C/4.35\n", " print('HG {} g/week'.format(HG))\n", "\n", " for r in ds['rows']:\n", " if r[variable]==None:\n", " continue\n", " #q is intake [ng/week]\n", " q=0\n", " #multiply by grams of average portion\n", " q=q+getWeeklyPortion(r['freshfish'])*150\n", " q=q+getWeeklyPortion(r['tinnedfish'])*80\n", " q=q+getWeeklyPortion(r['frozenfish'])*150\n", " #ignore weekly uptake below 50g\n", " if q<50: \n", " continue\n", " #this is intake in ng/week\n", " q*=q0\n", " #a mass of 1000 g/week is equivalent to 3M ng/week or 300 ng/min, tested for 48.6 ng/min\n", " #portion of Hg trapped in hair alpha=cH*HG/q\n", " thg.append(r[variable]*HG/q)\n", " return thg\n", " " ] }, { "cell_type": "code", "execution_count": 71, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "remoteSourcesURL https://git0.fmf.uni-lj.si/studen/nixSuite/raw/master/remoteResources/resources.json\n", "{'labkeyInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/labkeyInterface/archive/master.zip', 'branch': 'master', 'modules': []}, 'irAEMM': {'url': 'https://git0.fmf.uni-lj.si/studen/iraemm/archive/master.zip', 'branch': 'master', 'modules': ['iraemmBrowser']}, 'SlicerLabkeyExtension': {'url': 'https://git0.fmf.uni-lj.si/studen/SlicerLabkeyExtension/archive/SlicerExtensionIndex.zip', 'branch': 'SlicerExtensionIndex', 'modules': ['labkeyBrowser']}, 'limfomiPET': {'url': 'https://git0.fmf.uni-lj.si/studen/limfomiPET/archive/master.zip', 'branch': 'master', 'modules': ['imageBrowser']}, 'parseConfig': {'url': 'https://git0.fmf.uni-lj.si/studen/parseConfig/archive/master.zip', 'branch': 'master', 'modules': []}, 'orthancInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/orthancInterface/archive/master.zip', 'branch': 'master', 'modules': []}}\n", "User: andrej studen CSRF: b33922e96cef7ec47352de4af81bbc88\n", "HG 0.574712643678161 g/week\n", "remoteSourcesURL https://git0.fmf.uni-lj.si/studen/nixSuite/raw/master/remoteResources/resources.json\n", "{'labkeyInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/labkeyInterface/archive/master.zip', 'branch': 'master', 'modules': []}, 'irAEMM': {'url': 'https://git0.fmf.uni-lj.si/studen/iraemm/archive/master.zip', 'branch': 'master', 'modules': ['iraemmBrowser']}, 'SlicerLabkeyExtension': {'url': 'https://git0.fmf.uni-lj.si/studen/SlicerLabkeyExtension/archive/SlicerExtensionIndex.zip', 'branch': 'SlicerExtensionIndex', 'modules': ['labkeyBrowser']}, 'limfomiPET': {'url': 'https://git0.fmf.uni-lj.si/studen/limfomiPET/archive/master.zip', 'branch': 'master', 'modules': ['imageBrowser']}, 'parseConfig': {'url': 'https://git0.fmf.uni-lj.si/studen/parseConfig/archive/master.zip', 'branch': 'master', 'modules': []}, 'orthancInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/orthancInterface/archive/master.zip', 'branch': 'master', 'modules': []}}\n", "User: andrej studen CSRF: 6c2277c833f32f81249f198f67e02b07\n", "Found 1 rows for 181\n", "remoteSourcesURL https://git0.fmf.uni-lj.si/studen/nixSuite/raw/master/remoteResources/resources.json\n", "{'labkeyInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/labkeyInterface/archive/master.zip', 'branch': 'master', 'modules': []}, 'irAEMM': {'url': 'https://git0.fmf.uni-lj.si/studen/iraemm/archive/master.zip', 'branch': 'master', 'modules': ['iraemmBrowser']}, 'SlicerLabkeyExtension': {'url': 'https://git0.fmf.uni-lj.si/studen/SlicerLabkeyExtension/archive/SlicerExtensionIndex.zip', 'branch': 'SlicerExtensionIndex', 'modules': ['labkeyBrowser']}, 'limfomiPET': {'url': 'https://git0.fmf.uni-lj.si/studen/limfomiPET/archive/master.zip', 'branch': 'master', 'modules': ['imageBrowser']}, 'parseConfig': {'url': 'https://git0.fmf.uni-lj.si/studen/parseConfig/archive/master.zip', 'branch': 'master', 'modules': []}, 'orthancInterface': {'url': 'https://git0.fmf.uni-lj.si/studen/orthancInterface/archive/master.zip', 'branch': 'master', 'modules': []}}\n", "User: andrej studen CSRF: 3c9248cbfa5a372ae753d5ce68f63d78\n", "Found 1 rows for 179\n", "t0 2116800.0 tmax 2628000.0\n", "Parsing [1/24]\n", "Parsing [2/24]\n", "Parsing [3/24]\n", "Parsing [4/24]\n", "Parsing [5/24]\n", "Parsing [6/24]\n", "Parsing [7/24]\n", "Parsing [8/24]\n", "Parsing [9/24]\n", "Parsing [10/24]\n", "Parsing [11/24]\n", "Parsing [12/24]\n", "Parsing [13/24]\n", "Parsing [14/24]\n", "Parsing [15/24]\n", "Parsing [16/24]\n", "Parsing [17/24]\n", "Parsing [18/24]\n", "Parsing [19/24]\n", "Parsing [20/24]\n", "Parsing [21/24]\n", "Parsing [22/24]\n", "Parsing [23/24]\n", "Parsing [24/24]\n", "t0 2116800.0 tmax 2628000.0\n", "Parsing [1/24]\n", "Parsing [2/24]\n", "Parsing [3/24]\n", "Parsing [4/24]\n", "Parsing [5/24]\n", "Parsing [6/24]\n", "Parsing [7/24]\n", "Parsing [8/24]\n", "Parsing [9/24]\n", "Parsing [10/24]\n", "Parsing [11/24]\n", "Parsing [12/24]\n", "Parsing [13/24]\n", "Parsing [14/24]\n", "Parsing [15/24]\n", "Parsing [16/24]\n", "Parsing [17/24]\n", "Parsing [18/24]\n", "Parsing [19/24]\n", "Parsing [20/24]\n", "Parsing [21/24]\n", "Parsing [22/24]\n", "Parsing [23/24]\n", "Parsing [24/24]\n", "modelFile: /home/studen/temp/1673618348/model.json True\n", "parameterFile: /home/studen/temp/1673618348/parameters.json True\n", "bodyWeight -0.0019/0.26\n", "kHScaled 0.0079/0.25\n", "kKU:IPaper 1.4e-13/0.3\n", "kBBr:IPaper 1.1e-11/0.3\n", "kBScaled -0.00042/0.3\n", "kidneyPC -4.3e-05/0.3\n", "kBH:IPaper 1.8e-10/0.3\n", "kFScaled -0.00092/0.36\n", "dBLPaper -0.0054/0.3\n", "kBF:IPaper -5.1e-12/0.3\n", "fatPC -0.00011/0.3\n", "kKB:IPaper 7.6e-12/0.3\n", "slowlyPerfusedPC 0.0076/0.3\n", "brainPC -0.00016/0.3\n", "kRBCScaled -5e-09/0.3\n", "kRScaled 0.00078/0.3\n", "kLF:IPaper -5.5e-12/0.3\n", "kBU:IPaper 4.9e-13/0.3\n", "rbcPC -0.00077/0.3\n", "kBRScaled 7.5e-07/0.3\n", "kBrB:IPaper -7.8e-12/0.3\n", "gutPC -0.0001/0.7\n", "liverPC -0.00076/0.3\n", "kBK:IPaper -1.4e-12/0.3\n", "richlyPerfusedPC -0.00027/0.3\n", "brainBloodPC -0.00018/0.3\n", "kLB:IPaper -9e-08/0.3\n", "sigmaAcc=9.23322748349282e-06, dydpLN=[-0.0074262771098825295, 0.03169150752353229, 4.635148799007175e-13, 3.6272112812964996e-11, -0.0013867145751247662, -0.0001423623685739955, 6.013092022414007e-10, -0.002552045129243579, -0.01790811109460983, -1.7015112724085004e-11, -0.00036435788602090795, 2.545309234678283e-11, 0.025368169141099968, -0.0005313022462264693, -1.661001000413734e-08, 0.0025871224334682846, -1.8197643752355288e-11, 1.6331310527909675e-12, -0.0025624902151878243, 2.492026157979815e-06, -2.5854621014725752e-11, -0.00014905161176985318, -0.0025434814398012547, -4.668837990805158e-12, -0.0008901394255844984, -0.0005908213527405435, -2.9851220107701666e-07] cvLN=[0.26, 0.25, 0.3, 0.3, 0.3, 0.3, 0.3, 0.36, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.7, 0.3, 0.3, 0.3, 0.3, 0.3]\n", "qa=[-7.42627711e-03 3.16915075e-02 4.63514880e-13 3.62721128e-11\n", " -1.38671458e-03 -1.42362369e-04 6.01309202e-10 -2.55204513e-03\n", " -1.79081111e-02 -1.70151127e-11 -3.64357886e-04 2.54530923e-11\n", " 2.53681691e-02 -5.31302246e-04 -1.66100100e-08 2.58712243e-03\n", " -1.81976438e-11 1.63313105e-12 -2.56249022e-03 2.49202616e-06\n", " -2.58546210e-11 -1.49051612e-04 -2.54348144e-03 -4.66883799e-12\n", " -8.90139426e-04 -5.90821353e-04 -2.98512201e-07]\n", "A=0.059649291789389575,0.03704746964271335\n", "cvPrime=[0.26 0.25 0.3 0.3 0.3 0.3 0.3 0.36 0.3 0.3 0.3 0.3 0.3 0.3\n", " 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.7 0.3 0.3 0.3 0.3 0.3 ]\n", "sigmaS=[0.25575992 0.24622068 0.29356038 0.29356038 0.29356038 0.29356038\n", " 0.29356038 0.34908966 0.29356038 0.29356038 0.29356038 0.29356038\n", " 0.29356038 0.29356038 0.29356038 0.29356038 0.29356038 0.29356038\n", " 0.29356038 0.29356038 0.29356038 0.63148723 0.29356038 0.29356038\n", " 0.29356038 0.29356038 0.29356038]\n", "muS=[ -1.63988195 -0.66274584 -25.62375375 -21.26378776 -3.32835157\n", " -5.60466848 -18.45572768 -2.7362368 -0.77003516 -21.54445281\n", " -4.66490759 -21.61799972 -0.89807592 -4.28771312 -14.66079393\n", " -3.1810249 -21.47726262 -24.3643379 -2.71430938 -10.12623031\n", " -21.12606539 -5.71505076 -2.72175511 -22.83764452 -3.77166607\n", " -4.18153049 -11.77198865]\n", "Reading for sigma=0.24622067706923975\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.25575992365543143\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.34908965575742934\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.6314872286573718\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "Reading for sigma=0.293560379208524\n", "B=0.010327726481017901\n", "201/201\n" ] }, { "data": { "text/plain": [ "(array([ 6., 16., 22., 27., 20., 30., 26., 17., 17., 24., 21., 28., 17.,\n", " 23., 27., 26., 14., 16., 16., 15., 16., 17., 9., 13., 6., 9.,\n", " 12., 6., 11., 6., 8., 12., 5., 8., 4., 3., 2., 4., 3.,\n", " 2., 1., 0., 3., 3., 2., 2., 4., 0., 1., 2., 4., 4.,\n", " 3., 1., 1., 4., 3., 2., 2., 1., 0., 0., 1., 0., 1.,\n", " 0., 2., 2., 0., 0., 0., 0., 0., 2., 0., 0., 1., 0.,\n", " 0., 1., 0., 1., 0., 0., 0., 1., 0., 0., 0., 0., 1.,\n", " 0., 0., 0., 0., 0., 2., 0., 0., 0.]),\n", " array([0. , 0.00164648, 0.00329295, 0.00493943, 0.00658591,\n", " 0.00823239, 0.00987886, 0.01152534, 0.01317182, 0.0148183 ,\n", " 0.01646477, 0.01811125, 0.01975773, 0.02140421, 0.02305068,\n", " 0.02469716, 0.02634364, 0.02799012, 0.02963659, 0.03128307,\n", " 0.03292955, 0.03457603, 0.0362225 , 0.03786898, 0.03951546,\n", " 0.04116194, 0.04280841, 0.04445489, 0.04610137, 0.04774785,\n", " 0.04939432, 0.0510408 , 0.05268728, 0.05433376, 0.05598023,\n", " 0.05762671, 0.05927319, 0.06091966, 0.06256614, 0.06421262,\n", " 0.0658591 , 0.06750557, 0.06915205, 0.07079853, 0.07244501,\n", " 0.07409148, 0.07573796, 0.07738444, 0.07903092, 0.08067739,\n", " 0.08232387, 0.08397035, 0.08561683, 0.0872633 , 0.08890978,\n", " 0.09055626, 0.09220274, 0.09384921, 0.09549569, 0.09714217,\n", " 0.09878865, 0.10043512, 0.1020816 , 0.10372808, 0.10537456,\n", " 0.10702103, 0.10866751, 0.11031399, 0.11196047, 0.11360694,\n", " 0.11525342, 0.1168999 , 0.11854638, 0.12019285, 0.12183933,\n", " 0.12348581, 0.12513228, 0.12677876, 0.12842524, 0.13007172,\n", " 0.13171819, 0.13336467, 0.13501115, 0.13665763, 0.1383041 ,\n", " 0.13995058, 0.14159706, 0.14324354, 0.14489001, 0.14653649,\n", " 0.14818297, 0.14982945, 0.15147592, 0.1531224 , 0.15476888,\n", " 0.15641536, 0.15806183, 0.15970831, 0.16135479, 0.16300127,\n", " 0.16464774]),\n", " )" ] }, "execution_count": 71, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAh8AAAGeCAYAAAA0WWMxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABIrklEQVR4nO3de3xT9f0/8Nc5SZPe0pTeW3qhXMutqAi1XgClcnEqTn7zuonOeRvOr6LTdd+p4C44/U79ziH63VTcJqJuovMyHBdBRUColDuFlkILvUGhSS80aXI+vz/aBFIKNG2Sc5K8no9HHpLk5OSdY6Cvfj6f8z6SEEKAiIiIKEBktQsgIiKi8MLwQURERAHF8EFEREQBxfBBREREAcXwQURERAHF8EFEREQBxfBBREREAcXwQURERAHF8EFEREQBpVe7gO4URUFNTQ1MJhMkSVK7HCIiIuoFIQSam5uRkZEBWT7P2IbwwiuvvCLGjh0rTCaTMJlM4pJLLhGfffaZ+/nJkycLAB63++67z5u3ENXV1WfsgzfeeOONN954C45bdXX1eX/WezXykZmZiWeffRbDhg2DEAJvvfUWZs2aha1bt2L06NEAgHvuuQfPPPOM+zXR0dHevAVMJhMAoLq6GnFxcV69loiIiNRhtVqRlZXl/jl+Ll6Fj+uuu87j/m9/+1ssXrwYGzdudIeP6OhopKWlebNbD66plri4OIYPIiKiINObJRN9XnDqdDqxbNkytLa2orCw0P3422+/jaSkJIwZMwbFxcVoa2s7535sNhusVqvHjYiIiEKX1wtOd+zYgcLCQrS3tyM2NhbLly/HqFGjAAC33XYbcnJykJGRge3bt+OJJ55AWVkZPvjgg7Pub+HChViwYEHfPwEREREFFUkIIbx5gd1uR1VVFSwWC/7xj3/gL3/5C9atW+cOIKdbs2YNpk6divLycgwZMqTH/dlsNthsNvd915yRxWLhtAsREVGQsFqtMJvNvfr57XX46K6oqAhDhgzBa6+9dsZzra2tiI2NxYoVKzB9+vRe7c+b4omIiEgbvPn53e8mY4qieIxcnK60tBQAkJ6e3t+3ISIiohDh1ZqP4uJizJw5E9nZ2WhubsbSpUuxdu1afP7556ioqMDSpUtxzTXXIDExEdu3b8cjjzyCSZMmIT8/31/1ExERUZDxKnw0NDTgjjvuQG1tLcxmM/Lz8/H555/j6quvRnV1NVatWoWXXnoJra2tyMrKwuzZs/GrX/3KX7UTERFREOr3mg9f45oPIiKi4BPQNR9ERERE3mD4ICIiooBi+CAiIqKAYvggIiKigGL4oH4pq2vGS6v2ofJYq9qlEBFRkPD62i5EALDziAUvr9mPz3fVAwCWfHMQb945ARdmD1C5MiIi0jqGD/LK7hornv98L74oOwoAkCQg1RSJOms7bvvzJrz6o/GYPDxZ5SqJiEjLOO1CvdbYYsPNr23AF2VHIUvA9y8ciP88PAmrH52MScOTcbLDibuXbMZHpUfULpWIiDSMIx/Ua6+srUCzzYG8NBNe/eF4DEqKcT/3lzsuxmPvb8O/ttXgv5aV4nirHXddlqtitUREpFUc+aBeqbWcxN82HgIA/PKakR7BAwAMehkv3XwB7rx0EABgwce7safWGugyiYgoCDB8UK/8cXU57A4FBbkJuGJYUo/byLKEp68bhRmj0wAAf/7qQCBLJCKiIMHwQed18Fgr3ttSDQD4+fQRkCTprNtKkoQHpgwBAPyrtAa1lpMBqZGIiIIHwwed14ur9sGpCFw5IhkXD0o47/bjsuJRkJsAhyKwZP1B/xdIRERBheGDzmlvnRX/2lYDAHh02ohev+7eSYMBAEs3VaG5vcMvtRERUXBi+KBz+sN/9kEI4Htj0zFmoLnXr7tyRAqGJMeg2ebAu5ur/VghEREFG4YPOqvS6ias3F0PWQIeuXq4V6+VZQn3XNE5+vHG15XocCr+KJGIiIIQwwed1bubqwAAN1w4EENTYr1+/Q0XDkRSrBE1lnZ8tqPW1+UREVGQYvigHimKwKo9DQA6O5n2RWSEDndemgMA+L8vD0AI4bP6iIgoeDF8UI+2H7HgaLMNsUY9CnIT+7yf2wtyEBWhw64aK76paPRhhUREFKwYPqhHq3Z3Xq128ohkGPR9/5oMiDHgposzAXSe+UJERMTwQT1ataczfFw9MrXf+7qha9pm3b6jsDu48JSIKNwxfNAZqo+3YW9dM3SyhCkjkvu9v3GZ8UiKNaLF5sC3lcd9UCEREQUzhg86g2vUY8KgAYiPNvR7f7Is4aq8ZI99ExFR+GL4oDO4AkKRD6ZcXKZ27Wv13nqe9UJEFOYYPsiD5WQHNh3onBq5epTvwscVw5Jg0MuoPn4S+xtafLZfIiIKPgwf5GHdvqNwKALDUmKRkxjjs/1GG/S4dEjnKbsrd3PqhYgonDF8kAfXKbZFPhz1cHFPvXDdBxFRWGP4ILcOp4Ivyjq7mvpyvYfL1LwUAMDW6iYca7H5fP9ERBQcGD7IbXPlcTS3O5AUa8AFWfE+339GfBRGZ8RBCOCLvQ0+3z8REQUHhg9yW9k1HXJVXgp0suSX9zg19cLwQUQUrhg+CAAghPDLKbbdFY3snHr5av9R2BxOv70PERFpF8MHAQCONJ1E9fGT0MsSLh+W5Lf3GZNhRorJiFa7ExsPsNspEVE4YvggAMB3VU0AgJHpcYg26P32PrIsYWrX6AfPeiEiCk8MHwQA2Fp1AgBwUXa8399rat6pdR/sdkpEFH4YPgjAqZGPi3IG+P29LhuaBKNexpGmk6g4ym6nREThhuGD0N7hxO4aCwDgomz/h48ogw4Xdo2wfFt5wu/vR0RE2sLwQdh5xIIOp0BSrBGZA6IC8p4TByUAALYc5KJTIqJww/BB+O609R6S5J/+Ht1d3BU+vmX4ICIKOwwfhO8ONQEIzHoPl4tyBkCWgMMnTqLWcjJg70tEROrzKnwsXrwY+fn5iIuLQ1xcHAoLC/Hvf//b/Xx7ezvmzp2LxMRExMbGYvbs2aiv5+mUWiaEOG3kI3DhI9aox+gMMwBg80Gu+yAiCidehY/MzEw8++yzKCkpwZYtW3DVVVdh1qxZ2LVrFwDgkUcewccff4z3338f69atQ01NDW688Ua/FE6+caTpJBqabdDLEsYONAf0vS8e1Bl2Nldy6oWIKJx41U3quuuu87j/29/+FosXL8bGjRuRmZmJ119/HUuXLsVVV10FAHjzzTcxcuRIbNy4EZdcconvqiaf2Xpac7Eogy6g7z1xUALeXH8Qm7nug4gorPR5zYfT6cSyZcvQ2tqKwsJClJSUoKOjA0VFRe5t8vLykJ2djQ0bNpx1PzabDVar1eNGgfNdAJuLdedadFpW3wzLyY6Avz8REanD6z7aO3bsQGFhIdrb2xEbG4vly5dj1KhRKC0thcFgQHx8vMf2qampqKurO+v+Fi5ciAULFnhduObM72HKYr4l8HV4KZDNxbpLNhmRmxSDymOt+O7QCVyZlxLwGoiIKPC8HvkYMWIESktLsWnTJjzwwAOYM2cOdu/e3ecCiouLYbFY3Lfq6uo+74u8E+jmYj2Z0LXug6fcEhGFD69HPgwGA4YOHQoAGD9+PDZv3oz//d//xc033wy73Y6mpiaP0Y/6+nqkpaWddX9GoxFGo9H7yqnf1Ggu1t3FgxLw3pbDXHRKRBRG+t3nQ1EU2Gw2jB8/HhEREVi9erX7ubKyMlRVVaGwsLC/b0N+oEZzse5cnU63H7agvcOpSg1ERBRYXo18FBcXY+bMmcjOzkZzczOWLl2KtWvX4vPPP4fZbMbdd9+NefPmISEhAXFxcfjZz36GwsJCnumiUa7mYheqNOUCADmJ0Ug2GXG02Ybthy2YmJugWi1ERBQYXoWPhoYG3HHHHaitrYXZbEZ+fj4+//xzXH311QCAF198EbIsY/bs2bDZbJg+fTpeeeUVvxRO/ePZXCxetTokScKEQQPw2Y46bD54nOGDiCgMeBU+Xn/99XM+HxkZiUWLFmHRokX9Kor8r8bS7m4ulp8Zr2otEwYluMMHERGFPl7bJUx9d6hz1EON5mLdTeha91Fy6AScilC1FiIi8j+GjzC1rboJAHChilMuLiPT4xBr1KO53YGyuma1yyEiIj9j+AhTe7t+yI/OiFO5EkAnS+4mZ5x6ISIKfQwfYUgIgT21nW3s89LUDx8AMKErfGw5xCvcEhGFOoaPMHS02YbGVjtkCRiealK7HADABV3TPzsON6laBxER+R/DRxja0zXlMigpRvXFpi5jB3ZeG+dgYxssbbzIHBFRKGP4CEOuKZeR6dqYcgGA+GgDchKjAQA7jmj/gnxERNR3DB9haG9X+BilofABnBr92H6kSd1CiIjIrxg+wtCe2s5pl7w0baz3cBnX1exsezVHPoiIQhnDR5ixOZyoONoCQFvTLgAwNrNz5IPTLkREoY3hI8yUN7TAoQjEReqRbo5UuxwPYwaaIUnAkaaTONZiU7scIiLyE4aPMOOachmZHgdJklSuxlOsUY8hybEAgB2HOfpBRBSqGD7CzF4NnulyunzXolOGDyKikOXVVW1JI+abe3isdz+s99S5woe2Fpu6jM0044OtR7CDZ7wQEYUsjnyEkc626qemXbQov+uMl22HLRCCV7glIgpFDB9h5GizDcc11la9u1HpcdDJEo4221Bv5aJTIqJQxPARRnZ3rffITYpBZIQ22qp3F2XQYVhK56LT7bzOCxFRSGL4CCN7u67pkqfRKRcXd7MxLjolIgpJDB9hZI9G26p352o2tp3NxoiIQhLDRxjZq9G26t3lu8LH4SYuOiUiCkEMH2FCy23VuxuRZoJBJ6OprQOHT5xUuxwiIvIxho8w4Wqrbo6K0Fxb9e6Meh3yuvqQcN0HEVHoYfgIE6dfyVZrbdV7MnbgqakXIiIKLexw2hv96CiqFa7FpiOrlgLzZ5x6QqOfY1xmPN7eVMWRDyKiEMSRjzCx19VWXapSuZLecZ3xsvOIBYrCRadERKGE4SMMeLRVlw+pXE3vDEuJRWSEjGabAwcbW9Uuh4iIfIjhIwwca7HjeKsdEhQMk46oXU6v6HUy8tI6z8pxdWYlIqLQwPARBsobOk+xzZKOIkqyq1xN743K6Awfu2oYPoiIQgnDRxgo7+rvMTRIRj1cXJ1YdzN8EBGFFIaPMFBe37neY6hUo3Il3nGNfHDahYgotDB8hIFgHfno7EkCHG22oaG5Xe1yiIjIRxg+woBrzcdQObjCR7RBj9ykGACnmqQREVHwY/gIcdb2DtRbbQCCb+QDAEZndPb74LoPIqLQwQ6n/tS9M6oK3URdox4pJiPiOnq4SJsGajyXUelx+HhbDXbVaKsuIiLqO458hDhX+BiWGqtyJX3DRadERKGH4SPEVbjWeyQHafjoOt228lgr2uwOlashIiJfYPgIce7FpinBGT6STUYkm4wQAthbx0WnREShgOEjxLlOsx0SpOEDYLMxIqJQw/ARwto7nKg63gYAGJZiUrmavuO6DyKi0OJV+Fi4cCEmTJgAk8mElJQU3HDDDSgrK/PYZsqUKZAkyeN2//33+7Ro6p0DR1shBGCOikBSrEHtcvqMIx9ERKHFq/Cxbt06zJ07Fxs3bsTKlSvR0dGBadOmobXV85Ln99xzD2pra9235557zqdFU++4O5umxEKSJJWr6bvRXSMfe+uscCpC5WqIiKi/vOrzsWLFCo/7S5YsQUpKCkpKSjBp0iT349HR0UhLS/NNhdRn5UF+potLTmIMog06tNmdqDzWgqFBPIVERET9XPNhsXQ2fkpISPB4/O2330ZSUhLGjBmD4uJitLW1nXUfNpsNVqvV40a+Ud7QeXZIsPb4cNHJEvLSOgPHLk69EBEFvT6HD0VR8PDDD+Oyyy7DmDFj3I/fdttt+Pvf/44vvvgCxcXF+Nvf/oYf/vCHZ93PwoULYTab3besrKy+lkTduEY+gvlMFxcuOiUiCh19bq8+d+5c7Ny5E19//bXH4/fee6/7z2PHjkV6ejqmTp2KiooKDBky5Iz9FBcXY968ee77VquVAcQHHE4Flcc61+IE+7QLAIxK5zVeiIhCRZ/Cx4MPPohPPvkEX375JTIzM8+5bUFBAQCgvLy8x/BhNBphNBr7UgadQ9XxNnQ4BaIidBgYH6V2Of3mHvmosUIIEdQLaImIwp1X0y5CCDz44INYvnw51qxZg9zc3PO+prS0FACQnp7epwKpb/a7p1xiIMvB/4N6RKoJsgQ0ttrR0GxTuxwiIuoHr0Y+5s6di6VLl+Kjjz6CyWRCXV0dAMBsNiMqKgoVFRVYunQprrnmGiQmJmL79u145JFHMGnSJOTn5/vlA1DPQuVMF5cogw6Dk2NR3tCC3TVWpMZFql0SERH1kVcjH4sXL4bFYsGUKVOQnp7uvr377rsAAIPBgFWrVmHatGnIy8vDo48+itmzZ+Pjjz/2S/F0dhVBfk2XnozmolMiopDg1ciHEOdu8JSVlYV169b1qyDyjdMbjIWKUelx+Ki0hotOiYiCHK/tEoIURZx2NdvQaciV19VmfU8dwwcRUTBj+AhBtdZ2tNmd0MsSchKj1S7HZ0amdwapg8dacdLuVLkaIiLqqz73+SA/mW/u4TGLV7twjXoMSopBhC508mVyrBGJMQY0ttqxr74Z47Li1S6JiIj6IHR+MpFbqJ3p4iJJEkamn7rIHBERBSeGjxBUEYKLTV1c13jZU9usciVERNRXDB8h6EBX+BicHKNyJb7nGvnYw9NtiYiCFsNHCDpwtPOaLoNDbNoFAPLSXSMf1vOe+k1ERNrE8BFimts73O3Hc5NCb+RjaEos9LIEa7sDtZZ2tcshIqI+YPgIMQePtQEAkmINMEdFqFyN7xn1OgzpGtHholMiouDE8BFiDhzrWu+RFHpTLi6npl646JSIKBgxfISYCvd6j9CbcnHholMiouDG8BFiXGe6hOJ6D5dTp9syfBARBSN2OA0xlcdOO9Olp26pIWBU18hH5bFWtHc4ERmhU7kiIiLyBkc+QogQ4rTwEbojH8kmIxJiDFAEsL++Re1yiIjISwwfIaSu64JyOllCdkLoXFCuu84265x6ISIKVgwfIaSya7FpdkJ0SF1Qrid5aV2LTnm6LRFR0Antn1BhpsI15RLCi01dXItO9/J0WyKioMPwEUJC+Zou3blPt61jm3UiomDD8BFCQvmaLt0NTYmFTpbQ1NaBeqtN7XKIiMgLDB8hxHWmSyj3+HCJjNBhSNcIDxedEhEFF4aPEGFzOHH4ROd1XcJh2gXgolMiomDF8BEiDjW2QRGAyahHcqxR7XICwrXug4tOiYiCCzuc9iQIO4OevthUkiSVqwmMPPb6ICIKShz5CBEHwmi9h8vIrmmXA11t1omIKDgwfISIcDrTxSU1zogB0RFwKgLlDWyzTkQULBg+QkQ49fhwkSTp1KJTTr0QEQUNho8QccDd3TR8Rj6A0xad1nHRKRFRsGD4CAEnRCya2joAAIOSQveCcj3holMiouDD8BECDoh0AECGORLRhvA6gWlU+qlpF7ZZJyIKDgwfIaBCyQAQXotNXYamxEKWgBNtHWhoZpt1IqJgwPARAlwjH+G02NQlMkLnDl2ceiEiCg4MHyGgUqQBCK8eH6dzX+GWnU6JiIJCeC0QCFbn6bh6QITvtAsA5KWZ8PE2YC+v8UJEFBQ48hHknELCIZEKABgctiMfnWe88BovRETBgeEjyB0RybAjAga9jIz4KLXLUYVr2qXiaAtsDrZZJyLSOoaPIHfAtd4jMQY6OTwuKNddWlwkzFERcLDNOhFRUGD4CHKu9R7hutgU6GyzPtLdbIxTL0REWsfwEeRcZ7qE42m2p3Nd42UvT7clItI8r8LHwoULMWHCBJhMJqSkpOCGG25AWVmZxzbt7e2YO3cuEhMTERsbi9mzZ6O+vt6nRdMplV09PsJ55AM4teh0D894ISLSPK/Cx7p16zB37lxs3LgRK1euREdHB6ZNm4bW1lb3No888gg+/vhjvP/++1i3bh1qampw4403+rxw6nRACd8GY6c7vdcH26wTEWmbV30+VqxY4XF/yZIlSElJQUlJCSZNmgSLxYLXX38dS5cuxVVXXQUAePPNNzFy5Ehs3LgRl1xyie8qJ5wUBtQgCQCQG2ZXs+1ueKoJsgQcb7XjaLMNKXGRapdERERn0a81HxaLBQCQkJAAACgpKUFHRweKiorc2+Tl5SE7OxsbNmzocR82mw1Wq9XjRr1zsKu/RzyakRBjULkadUVG6NxTT3vquOiUiEjL+tzhVFEUPPzww7jsssswZswYAEBdXR0MBgPi4+M9tk1NTUVdXV2P+1m4cCEWLFjQ1zLUc56uo4HgXu8h9Xxsw01eehwqjrZib60Vk4cnq10OERGdRZ9HPubOnYudO3di2bJl/SqguLgYFovFfauuru7X/sLJAXf4qFW5Em0YmeY63ZajZ0REWtankY8HH3wQn3zyCb788ktkZma6H09LS4PdbkdTU5PH6Ed9fT3S0tJ63JfRaITRaOxLGWHPvdhUZvgATi063ctpFyIiTfNq5EMIgQcffBDLly/HmjVrkJub6/H8+PHjERERgdWrV7sfKysrQ1VVFQoLC31TMbm5pl0Gc+QDQOe0CwCUN7DNOhGRlnk18jF37lwsXboUH330EUwmk3sdh9lsRlRUFMxmM+6++27MmzcPCQkJiIuLw89+9jMUFhbyTBcfE4LTLt1lmCMRF6mHtd2BioZWjMqIU7skIiLqgVcjH4sXL4bFYsGUKVOQnp7uvr377rvubV588UVce+21mD17NiZNmoS0tDR88MEHPi883J2ACRZ0nl47SGITN6CzzXqeu98H130QEWmVVyMfvWneFBkZiUWLFmHRokV9LorOz9VWPQPHECXZVa5GO0ammfBt5XHsZadTIiLN4rVdglSF0nlBOS429eRadLqbIx9ERJrF8BGkXCMfXO/hiW3WiYi0j+EjSLHBWM9GpJ1qs97QbFO7HCIi6gHDR5ByjXwMlmpUrkRbIiN0GJLcuRB3dw2nXoiItIjhIwg5hXRa+ODIR3euU2y57oOISJsYPoJQjUiEHQZEwIGB0lG1y9GcUa5Fpxz5ICLSJIaPIORa75Ej1UEncVFldxz5ICLSNoaPIHSAbdXPyXXGy8HGVrTaHCpXQ0RE3TF8BKFKtlU/p6RYI1JMRgjBi8wREWkRw0cQOjXywcWmZ8OpFyIi7WL4CEIHlK6RD3Y3PSsuOiUi0i6GjyDTLiJQg0QAXPNxLhz5ICLSLoaPIHNIpEJAhgmtSAR/sJ6Na+SjrM4Kp8IzgoiItIThI8icvt5DklQuRsNyEmMQFaFDe4eCymOtapdDRESnYfgIMgd4pkuv6GQJeekmAJx60aopU6bg4Ycf7vPrlyxZgvj4eJ/VQ0SBw/ARZFyLTQdzsel5cdFpaLv55puxb98+tcsgoj7Qq10AeeeAyADAC8r1BhedhraoqChERUWd9Xm73Q6DwRDAioiotzjyEUSEACq6wscQho/z4siH9imKgscffxwJCQlIS0vD/Pnz3c+98MILGDt2LGJiYpCVlYWf/vSnaGlpcT/ffdpl/vz5uOCCC/CXv/wFubm5iIyMDOAnISJvMHwEkUbEwYoYSFCQywZj5zUizQRJAo612NDQ3K52OdSDt956CzExMdi0aROee+45PPPMM1i5ciUAQJZl/PGPf8SuXbvw1ltvYc2aNXj88cfPub/y8nL885//xAcffIDS0tIAfAIi6gtOuwQR16jHQOkYIqUOlavRvmiDHrlJMThwtBV7apuRYuJvwlqTn5+Pp59+GgAwbNgw/OlPf8Lq1atx9dVXeyxGHTRoEH7zm9/g/vvvxyuvvHLW/dntdvz1r39FcnKyv0snon7gyEcQcS825ZkuvcapF23Lz8/3uJ+eno6GhgYAwKpVqzB16lQMHDgQJpMJP/rRj9DY2Ii2traz7i8nJ4fBgygIMHwEEa738B4XnWpbRESEx31JkqAoCg4ePIhrr70W+fn5+Oc//4mSkhIsWrQIQOfoxtnExMT4tV4i8g1OuwSRUw3GOPLRW6dGPiwqV0LeKCkpgaIo+MMf/gBZ7vwd6b333lO5KiLyFYaPQJpv7tfLzzny0c99B0T3Guf7PxC4Rj4qj7XipN2JKIPO7+9J/Td06FB0dHTg5ZdfxnXXXYf169fj1VdfVbssIvIRTrsECZvQo1qkAACGypx26a0UUySSYg1QBLCnjlMvwWLcuHF44YUX8Pvf/x5jxozB22+/jYULF6pdFhH5iCSE0NRVt6xWK8xmMywWC+Li4tQpQoOjCPuUgZhmfx6xaMMO40/8d10Xf45GqDDyAQBz3vgW6/Ydxa9njcaPCgcF5D2JiMKNNz+/OfIRJA64p1xqeUE5L40d2Bl6dh7hyAcRkRYwfASJCrZV77MxAzsT+I4jXHRKRKQFDB9BoqKrx8cQrvfw2piukY999c2wOZwqV0NERAwfQeLUyAdPs/XWwPgoDIiOgEMRKKtrVrscIqKwx/ARBIQ41eODDca8J0mSe/SD6z6IiNTH8BEEjsKM5q4LyuVI9WqXE5Rc4YPrPoiI1MfwEQRcox5Z0lFeUK6PxmR0ho9d7HRKRKQ6djjVYE+P7ioUldd79OYYde/ZobHj6jrddm9tM+wOBQY9czcRkVr4L3AQOMALyvVbVkIU4iL1sDsV7KvnolMiIjUxfAQB9vjov9MXnXLqhYhIXQwfQcB9pgt7fPTLWC46JSLSBIYPjWsXEagWyQDY46O/RvN0WyIiTWD40LhDIhUCMkxoRTL4G3t/uEY+9tRa4XAqKldDRBS+vA4fX375Ja677jpkZGRAkiR8+OGHHs/feeedkCTJ4zZjxgxf1Rt2Tu9sygvK9U9OQjRijXrYHArKj7aoXQ4RUdjyOny0trZi3LhxWLRo0Vm3mTFjBmpra923d955p19FhjN2NvUdWZYwOqPrInOHOYpERKQWr/t8zJw5EzNnzjznNkajEWlpaX0uik5x9fgYInO9hy+MGWjGpsrj2FVjxQ/ULoaIKEz5Zc3H2rVrkZKSghEjRuCBBx5AY2PjWbe12WywWq0eNzqFPT58i2e8EBGpz+cdTmfMmIEbb7wRubm5qKiowC9/+UvMnDkTGzZsgE6nO2P7hQsXYsGCBb4uIyQIAVR0TbtovseHxjqano2r18fuGiucioBO5kIaIqJA83n4uOWWW9x/Hjt2LPLz8zFkyBCsXbsWU6dOPWP74uJizJs3z33farUiKyvL12UFpQbEowXRkHlBOZ/JTYpBtEGHNrsTB462YFiqSe2SiIjCjt9PtR08eDCSkpJQXl7e4/NGoxFxcXEeN+q0X8kEAAyS6mCUHCpXExp0py063clOp0REqvB7+Dh8+DAaGxuRnp7u77cKOftEZ/gYJh1WuZLQMrrrCrc7DnN9ERGRGryedmlpafEYxaisrERpaSkSEhKQkJCABQsWYPbs2UhLS0NFRQUef/xxDB06FNOnT/dp4eFgvxgIABgmHVG5ktDiWnS6/XCTuoUQEYUpr8PHli1bcOWVV7rvu9ZrzJkzB4sXL8b27dvx1ltvoampCRkZGZg2bRp+/etfw2g0+q7qMOGadhkmM3z40riseACd0y4dTgUROjb6JSIKJK/Dx5QpUyCEOOvzn3/+eb8Kok5CnD7ywWkXXxqcFANTpB7N7Q6U1TW7z4AhIqLA4K98GnUUZlgQCxkKLyjnY7Is4YKu0Y9tnHohIgo4hg+Nck255Ej1iJQ6VK4m9IzLjAcAlFY1qVoHEVE4YvjQKNeUy1AuNvUL18hHaXWTqnUQEYUjnzcZI9/Y7z7NNsDhI0g6lfaXa9Fp+dEWNLd3wBQZoW5BRERhhCMfGrVf6Rz5GC5zsak/JJuMGBgfBSF4hVsiokBj+NAgIU41GBvKM1385oLseADAVk69EBEFFMOHBh1DHJpgggQFQ7V+QbkgdiHXfRARqYLhQ4NcZ7pkSw0808WPxp0WPs7Vu4aIiHyL4UODytlWPSDGZJihkyUcbbah1tKudjlERGGD4UODeEG5wIgy6JCXZgLAqRciokBi+NAg15kuvKaL/7mmXrYxfBARBQzDhwa5pl2Gc+TD71zNxnjGCxFR4DB8aEyjMKERZkhQMIRnuvid64yXHYctcDgVdYshIgoT7HCqMa626pnSMURJdpWrUUH3Dqvz/dsAbHByLGKNerTYHNhX34JRGXF+fT8iIuLIh+a4TrPllEtg6GQJ+ZmdgYdXuCUiCgyGD43Zz86mAee+yByvcEtEFBAMHxrjmnYZzjNdAmYcO50SEQUUw4fGuKZd2OMjcFyLTvc1NKPF5lC3GCKiMMDwoSHHhQnH0Ln+gGe6BE5KXCQyzJEQAtjO0Q8iIr9j+NCQcpEBAMiUGhAj2VSuJryMH5QAANh88ITKlRARhT6GDw3Z555y4XqPQJswaAAAYMuh4ypXQkQU+hg+NGS/YPhQy4SukY/vDp1gszEiIj9j+NCQPUo2ACBPrlK5kvAzPNUEU6QerXYn9tQ2q10OEVFIY4dTjRAC2CM6w8co6ZDK1QSZ7l1RAa87o+pkCeNzBmBt2VFsPngcY7sajwW64yoRUTjgyIdGHBZJaEYMDOjgmS4qcU29cN0HEZF/MXxoxG4xCAAwVDqCCMmpbjFhasJpZ7wIIVSuhogodDF8aIRrymUkp1xUk59phkEn42izDYca29Quh4goZDF8aIRrselILjZVTWSEzn2Ruc0HOfVCROQvDB8asUfkAOBiU7Vd7Fr3wWZjRER+w/ChAc0iClUiFQBHPtTmaja2mYtOiYj8huFDA/Z2rfdIRyMGSC0qVxPeLs7pHPk4cLQVx1rY4p6IyB8YPjTg1HoPTrmozRwdgRGpJgCceiEi8heGDw1wrfcYKXHKRQsudl3nhYtOiYj8gh1ONWA3Rz78z4suqBMGJeDtTVXYfIgjH0RE/sCRD5U5hYQykQWAZ7poxYTcznUfu45Y0CaMKldDRBR6GD5UVinS0Q4jotCOHKle7XIIwMD4KGSYI+FQBEqVIWqXQ0QUchg+VObqbDpCqoZOYktvrXD1+9gsRqhcCRFR6GH4UNkepWuxKft7aIpr6mWTMlLlSoiIQg/Dh8pcIx9c76EthYMTAQBblOFoFxEqV0NEFFq8Dh9ffvklrrvuOmRkZECSJHz44Ycezwsh8NRTTyE9PR1RUVEoKirC/v37fVVvyNndNfIxime6aMqQ5BikxUXCDgNKlOFql0NEFFK8Dh+tra0YN24cFi1a1OPzzz33HP74xz/i1VdfxaZNmxATE4Pp06ejvb2938WGmuPChHp0Du+PkKpVroZOJ0kSLh3aOfrxtTJG5WqIiEKL130+Zs6ciZkzZ/b4nBACL730En71q19h1qxZAIC//vWvSE1NxYcffohbbrmlf9WGGFdn0xypDrESw5nWXD40CR98dwTrlTEA3lW7HCKikOHTNR+VlZWoq6tDUVGR+zGz2YyCggJs2LChx9fYbDZYrVaPW7hgZ1Ntu2xoEgBgh8iFRcSoXA0RUejwafioq6sDAKSmpno8npqa6n6uu4ULF8JsNrtvWVlZvixJ09jZVNtS4yIxVDoMARkblFFql0NEFDJUP9uluLgYFovFfauuDp+1D7u7Rj54pot2XS7vBMB1H0REvuTT8JGWlgYAqK/37NRZX1/vfq47o9GIuLg4j1s4sAsdKsRAAOzxoWWXdYWP9QwfREQ+49PwkZubi7S0NKxevdr9mNVqxaZNm1BYWOjLtwp6e0U2OqBHPJoxEMfULofOokDeAxkKKkU6johEtcshIgoJXp/t0tLSgvLycvf9yspKlJaWIiEhAdnZ2Xj44Yfxm9/8BsOGDUNubi6efPJJZGRk4IYbbvBl3UFvW9c1Q8bJFZAklYuhs4qTTmKcVIGtYhjWO8fgJrULIiIKAV6Hjy1btuDKK6903583bx4AYM6cOViyZAkef/xxtLa24t5770VTUxMuv/xyrFixApGRkb6rOgRsVYYCAMZJFSpXQudzubwTW53DsF5h+CAi8gWvw8eUKVMgxNkvgCZJEp555hk888wz/Sos1G0TnSMfF8gMH1p3qbwTLzu/j/XKaAghIHGoioioX1Q/2yUcWUWUe7HpOIYPzbtI3o9I2HAM8dhX36J2OUREQY/hQwXbu9Z7ZEkNSJSaVa6GzscoOTBBLgMAfF3OxcFERP3F8KEC95SLVH6eLUkrXP0+1jN8EBH1m9drPqj/tp52pktYm28O7Ov64TJ5BwBg04FGdDgVROiY24mI+or/ggaYEEBp15kuXGwaPEZJVRiAZrTanfju0Am1yyEiCmoMHwFWg0QcQzx0cGKMVKl2OdRLsiQwWd4GAFizt0HlaoiIghvDR4C5movlSVWIlDpUroa8MVX3HQBg1Z7682xJRETnwvARYK7wwSmX4DNZ3ga9LKHiaCsqj7WqXQ4RUdBi+AgwdjYNXnHSSUzMTQAArOboBxFRnzF8BJBDyNghcgEAF8g8zTYYTR2ZCoBTL0RE/cHwEUD7xUCcRCRi0YYhUo3a5VAfFI1MAQBsPngCljau2SEi6guGjwByrfcYK1dCJ539+jikXTmJMRiWEgunIrB2H896ISLqC4aPACoVXf092Nk0qLmmXlbvYfggIuoLdjgNoFJ2NtW2XnZOLRqZglfXVWBtWUNnt9NfD+i2H4sfiiMiCh0c+QiQNmHEPpEFgKfZBrsLswcgIcYAa7sDWw6y2ykRkbcYPgJkh8iFAhlpaESaxB9YwUwnS7hyROfCU571QkTkPYaPANnGKZeQ4jrrZfWeegiuHSYi8grDR4BsVkYA4JRLqLhieDIMOhkHG9tQITLULoeIKKgwfASAU0jYpIwEAFwi71a5GvKFWKMeBYO7up0qF6lcDRFRcGH4CIA9IhtWxCAWbRjLK9mGjKtHdZ5y+7nzYpUrISIKLgwfAbBRGQ0AmCCXQS8pKldDvjJ9dBokCfhODEe1kqR2OUREQYPhIwA2KKMAAIWccgkpqXGRuCQ3EQDwsVKocjVERMGD4cPPHELGt0oeAIaPUDTrgs7Fpv9yXqpyJUREwYPhw892ixw0IxomtGKUdFDtcsjHZo5JRwQc2CtysF8ZqHY5RERBgeHDzzZ0rfcokPfyYnIhyBwdgcnyNgAc/SAi6i2GDz9zrffgKbah6zrdNwCAj5RL2XCMiKgXGD78qEPo3M3FuN4jdF0tf4cotKNKpGKbGKJ2OUREmsfw4Uc7xSC0IgpmtGCkVKV2OeQn0ZINV8slAIB/OXnWCxHR+TB8+JFryqVA3gOZ6z1C2vW6DQCAT5yFcCr8f01EdC4MH37kWmzKKZfQN0neBjNa0IAB2HSgUe1yiIg0jeHDT+xChy3KcAAMH+HAIDkxU/ctAOBf22pUroaISNsYPvxkhxiMk4hEAqwYLh1WuxwKgOvlzrNe/r2zDnYH2+gTEZ0Nw4efnDrFlus9wkWBvAcpOAHLyQ6s3lOvdjlERJqlV7uAgJtvDsjbuNZ7sL9HP/X1/1dvXufj74JOEviBbh0WOW/A25uqMHNsuk/3T0QUKjjy4Qc2oUeJMgwA13uEm1t0ayBJwNflx1B5rFXtcoiINInhww82KKPQDiNScAJDpSNql0MBlCUfw5UjUgAAb288pHI1RETaxPDhByuViwEAV+tKIEkqF0MBd3tBNgDgH98dRnuHU+VqiIi0h+HDxxQhYaVzPADganmLytWQGqaMSMHA+Cg0tXXg0+21apdDRKQ5DB8+tl0MRgMGIBZtXO8RpnSyhFsnZgEA3t7EqRciou58Hj7mz58PSZI8bnl5eb5+G836j7NzymWyvA1GyaFyNaSWmyZkQS9L+K6qCbtrrGqXQ0SkKX4Z+Rg9ejRqa2vdt6+//tofb6NJK5XOKZdpOk65hLMUUySmj0kDAPydox9ERB78Ej70ej3S0tLct6SkJH+8jeZUKmnYLzKhhwNT5G1ql0Mqcy08/WjrEbTYOApGROTil/Cxf/9+ZGRkYPDgwbj99ttRVXX2y8nbbDZYrVaPW7ByjXoUyrthltpUrobUVjg4EYOTY9Bqd2L5Vp5yTUTk4vPwUVBQgCVLlmDFihVYvHgxKisrccUVV6C5ubnH7RcuXAiz2ey+ZWVl+bqkgHGt97haLlG5EvKb+eYzb2chSRJuL8gBACxZXwlFYZt9IiLAD+Fj5syZ+MEPfoD8/HxMnz4dn332GZqamvDee+/1uH1xcTEsFov7Vl1d7euSAuKYiEOJ6OxqWqRj+KBON12cibhIPSqOtmLFrjq1yyEi0gS/n2obHx+P4cOHo7y8vMfnjUYj4uLiPG7BaLXzIgjIGCsdQIZ0XO1ySCNMkRG487JcAMDLa8ohBEc/iIj8Hj5aWlpQUVGB9PTQvsiWa73H1Rz1oG5+fNkgxBh02FNrxZq9DWqXQ0SkOp+Hj8ceewzr1q3DwYMH8c033+D73/8+dDodbr31Vl+/lWa0CSO+UsYCAKaxqyl1Ex9twI8KBwHg6AcREeCH8HH48GHceuutGDFiBG666SYkJiZi48aNSE5O9vVbacaXyljYYECW1IARUnCuWSH/+skVuYiMkFFa3YT15Y1ql0NEpCq9r3e4bNkyX+9S8z51XgKg81ouvJAc9SQp1ohbJ2bjzfUH8fKa/bh8WHj0viEi6gmv7dJPJ0QsPlcmAAC+r1uvcjWkZfdOGgyDTsamyuP4tpKLkokofDF89NNy5+WwIwKjpIMYI1WqXQ5pWLo5Cv/v4kwAwJ++6PnsLyKicMDw0Q9CAO86pwAAbtWt4ZQLndcDk4dAJ0v4ct9RlBzi6AcRhSefr/kIJ6ViCMpENoyw43rdN2qXQ1rRU9fT+RYAQFZCNH4wPhPLNldjwce78eFPL4MsSz2/rus1mnGOz9Xv/fhivz3tW2vHkIgAcOSjX951XgkA+J68iddyoV57dNoImIx6bD9swT++O6x2OUREAcfw0UctIhL/cl4KALhZ/4XK1VAwSTYZ8dDUzlb8z60oQ3N7h8oVEREFFsNHH33ivARtiMRgqQYTpb1ql0NBZs6lgzA4KQbHWmz40xouPiWi8MLw0UfLuqZcbtat5UJT8ppBL+PJa0cBAN5YX4kDR1tUroiIKHAYPvpgr5KFUjEMejhwo+5LtcuhIHVlXgqmjEhGh1Pgt5/uUbscIqKAYfjoA9eoR5H8HZIlq8rVUDB78tpR0MsSVu9twFpnvtrlEBEFBMOHl9pFBJY7LwcA3KJbo3I1FOyGJMfirssGAQCedPwYLSJS3YKIiAKA4cNLf3cWwYJYZEoNuELeoXY5FAIemjoMA+OjUC1SsMBxh9rlEBH5HZuMeaFNGPGq43oAwIO6D6GTeGl06qVzNL8yRUbgxZsvwM2vrcf7zim4St6KmbrNwdHUy1/8+dm19lmJwhBHPrzwd2cRjsGMLKkBs3VfqV0OhZCJuQl4QPcxAKC44yeoF/HqFkRE5EcMH73UKox41XEdAOBnuuWIkJwqV0Sh5mH9PzBWOoAmmPBYx/1QBM/hJqLQxPDRS391TsNxxCFHqsONHPUgPzBITrwY8QoiYcNXSj7eck5TuyQiIr9g+OiFFhGJ/3NcCwD4mX459JKickUUqobKNfhv/dsAgIWOW7FDyVW5IiIi32P46IW3nNNxAibkSrW4QV6vdjkU4n6oW4Wp8neww4C77Y+hRiSoXRIRkU8xfJxHs4jC/zm+BwB4SP8BRz3I7yQJeDFiEUZIVWjAAPzY/jiaRZTaZRER+QzDx3n82XENLIjFYKkG18vfqF0OhYk46STeMDyPZJzAXpGNBzt+BofgX1ciCg381+wcypRMLHbOAgA8qn+ffT0ooAZKjXjd8D+IhA3rlAvwtONOCMHvIBEFP4aPs3AIGT/vuA8d0KNILsE18ia1S6IwlC9X4n8jFkGCgredRXh13QG1SyIi6jd2OD2LPzu/h+1iCOLQit9GvA6JLRfIl3rTdbTLdN0W/LdYit84fojfr9gL+6pf4yHd8s7vpK+6dXpRj6b23Rda7/ja2+PlrxrZFVbbztEtOZhw5KMH5UoGXnTMBgA8qf8bUqUmdQuisHe37jPM078PAHjR8QP8znEbOANDRMGKIx/dOIWEn3fcBzsMmCKX4v/pvlS7JCJIEvCQfjlicRLPOO7An53XogXR+I0ioJM5LEdEwYUjH9286ZyJrWIYYtGG30X8hdMtpCk/1q/Ac/rXIEPBO86r8Mi7pbA52OqfiIILw8dptiu5eN5xEwDgv/VvI0M6rnJFRGe6Sb8Of4x4GXo48K9tNfjBqxtQfbxN7bKIiHqN4aPLYZGEH9t/DhsMuEr+DrfovlC7JKKzula3CW9EPI/46AhsP2zB9/74FVbtrle7LCKiXmH4AGAVUfix/ec4hnjkSYc6T23kdAtp3CTdDnz60BW4ICse1nYHfvLXLVj47z1wONmFl4i0LezDh13o8EDHI9gnspCK43jT8DxM0km1yyLqlYHxUXjvvkLcddkgAMBr6w7gxsXfYMfh4Dz9jojCQ1iHDyGA/3bcjfXKGESjHa8bnkc613lQkDHoZTx93Wi8cvtFMBn12H7YgusXfY0nP9wJS1uH2uUREZ0hbMOHEMCzjlvwvnMKZChYFPFHjJEPqV0WUZ9dMzYdqx+djBsuyIAQwN82HsJVf1iLfzivgFNwHpGItCMs+3zYhQ6Pd9yHD5XLAQAL9Etwpa5U3aKI+qJbt8MUAC/Nt+CmCVl46qNdKG9owWN4AK9IszBX/xGul79BhOTHU3NDpVNqX7pI9qYzaCh3D9X6MfPnsQ+RrqOBFHYjH1YRhTs7nsCHyuXQw4Hn9K/hR/pVapdF5FOXDknCZw9dgSdm5CEOrTggMvBoxwO4yv4HvO24CjYRlr93EJFGhFX4qGk6iR/Yn8Y3yhjE4CTeiHgeN+nXqV0WkV8Y9DIemDIE640P4Qn9O0iEBdUiBf/t+AkKbIswv+MO7FJy1C6TiMJQ2Pz6s6fWirve3Iw6kY0UnMCbhucwmms8KAyYpJN4QP8x7tR9jmXOK/F/jmtRi0Qscc7AEucMjJIO4v/pvsQ03RZkSsfULpeIwkDYhA+9LKHN7sAw6TCWGH6PgVKj2iURBVSUZMdd+s9xh+4/+EoZi/edk7FSuRi7xSA84xiEZxx3IE+qwlT5O0zVfYdxUgV0Eq9eR0S+FzbhY1iqCX//SQFy/nwXzFKr2uUQqUYnCUzRbccU3XacELH4yHkpPnMWYIsYgb0iG3ud2VjkvAEmtOEieR8myGUYL+3DBXIFoiS72uUTUQjwW/hYtGgRnn/+edTV1WHcuHF4+eWXMXHiRH+9Xa/kZ8YDDB5EbgOkFtyp/w/u1P8HJ0Qs1injsMp5EdYp+WhGDNYpF2CdcgEAQAcnBku1yJOqkCdXYaRUhaGNbciIj4ReF1bLx4ion/wSPt59913MmzcPr776KgoKCvDSSy9h+vTpKCsrQ0pKij/ekoj6aYDUght063GDbj2cQsJekY0tyghs7rrVIwH7RSb2i0x8rFza+aLnv4BelpA5IAo5iTHI6bgTaVIj0qXjSMNxpEknkCw1IQbtvGQBEbn5JXy88MILuOeee3DXXXcBAF599VV8+umneOONN/CLX/zCH29JRD6kkwRGS4cwWj6EOfgPhADqMQB7lSzsETnYq2Rjr8hCpW4Q7A4FBxvbcLCxDcC0HvdngB2JaEaCZEWC1Iw4tCJOakMcWmGWWmHCSURLNsSgHdFoR4zUDiPsiEQHIiU7ImGHER0woAMGOBhkiIKcz8OH3W5HSUkJiouL3Y/JsoyioiJs2LDhjO1tNhtsNpv7vtVq9XVJRNRPkgSk4QTSdCcwBdvdjytPNaHO2o6Dja041NiGqo9+g3qRgDoMQJ1IQJ1IQBsiYYcBtUhErUgEfLCG1QA7jHBADyf0cMAAB/SSE3o4oYMCHZzQQ4EMBToo0MMJWVKgg4Dc9bgMARkCEhRIgMdjgIAEQILovC3b2nlfkiABgAR0/QmShM4/ddwLqevDubKR9M9Tx0qSAHT8BN0PgLR8h+eH6/gxepOtpA93nH+jvui4q4f32tmL193peb8vr+npdR13nnk8PurFvnvxXlJf9tMDqeMOzwf+tcsn++3RWd5LiM7vlUBnB28B0fXfrvui874iBBQBpMQZ8cSMPP/VeR6ScFXsIzU1NRg4cCC++eYbFBYWuh9//PHHsW7dOmzatMlj+/nz52PBggVn7MdisSAuLs6XpRGRCtrsDjS22HG89dStub0DlpMOWNs7YDnZgTa7Ay02J9psDrTYHGizO9He0XVzKLA7eKVeIl8anByDNY9O8ek+rVYrzGZzr35+q362S3FxMebNm+e+b7VakZWVpWJFRORL0QY9ohP0yEqI7vM+FEXA7lRg6woidmfnfx1OBR1OgQ6nAoeiwOEUcCoCDsXzv4o49V9FCCgKTv1ZwP0boXDfF6f9Bnnab5Wn/UYJnPqtsvPPp7bp7vTXezzew2ftza+Doo/DR779VbP/AlqOjz58X/bS17fu7f9nCZJ7KrJzZE46NfrWNUonS4Asdz4qSxIGREf0rSgf8Xn4SEpKgk6nQ319vcfj9fX1SEtLO2N7o9EIo9Ho6zKIKITIsoRIWYfICJ3apRCRD/j8/DiDwYDx48dj9erV7scURcHq1as9pmGIiIgoPPll2mXevHmYM2cOLr74YkycOBEvvfQSWltb3We/EBERUfjyS/i4+eabcfToUTz11FOoq6vDBRdcgBUrViA1NdUfb0dERERBxOdnu/SXN6tliYiISBu8+fnNnshEREQUUAwfREREFFAMH0RERBRQDB9EREQUUAwfREREFFAMH0RERBRQDB9EREQUUAwfREREFFAMH0RERBRQfmmv3h+uhqtWq1XlSoiIiKi3XD+3e9M4XXPho7m5GQCQlZWlciVERETkrebmZpjN5nNuo7lruyiKgpqaGphMJkiS5NN9W61WZGVlobq6mteNAY9HT3hMPPF4nInHxBOPx5nC9ZgIIdDc3IyMjAzI8rlXdWhu5EOWZWRmZvr1PeLi4sLqC3E+PB5n4jHxxONxJh4TTzweZwrHY3K+EQ8XLjglIiKigGL4ICIiooAKq/BhNBrx9NNPw2g0ql2KJvB4nInHxBOPx5l4TDzxeJyJx+T8NLfglIiIiEJbWI18EBERkfoYPoiIiCigGD6IiIgooBg+iIiIKKAYPoiIiCiggip8LFq0CIMGDUJkZCQKCgrw7bffnnP7999/H3l5eYiMjMTYsWPx2WefeTwvhMBTTz2F9PR0REVFoaioCPv37/fY5vjx47j99tsRFxeH+Ph43H333WhpafH5Z+srXx6Tjo4OPPHEExg7dixiYmKQkZGBO+64AzU1NR77GDRoECRJ8rg9++yzfvl83vL1d+TOO+8847POmDHDYxstf0d8fTy6HwvX7fnnn3dvo+XvB+DdMdm1axdmz57t/kwvvfRSn/bZ3t6OuXPnIjExEbGxsZg9ezbq6+t9+bH6zNfHY+HChZgwYQJMJhNSUlJwww03oKyszGObKVOmnPEduf/++3390frM18dk/vz5Z3zevLw8j220/B3xCxEkli1bJgwGg3jjjTfErl27xD333CPi4+NFfX19j9uvX79e6HQ68dxzz4ndu3eLX/3qVyIiIkLs2LHDvc2zzz4rzGaz+PDDD8W2bdvE9ddfL3Jzc8XJkyfd28yYMUOMGzdObNy4UXz11Vdi6NCh4tZbb/X75+0NXx+TpqYmUVRUJN59912xd+9esWHDBjFx4kQxfvx4j/3k5OSIZ555RtTW1rpvLS0tfv+85+OP78icOXPEjBkzPD7r8ePHPfaj1e+IP47H6cehtrZWvPHGG0KSJFFRUeHeRqvfDyG8PybffvuteOyxx8Q777wj0tLSxIsvvtinfd5///0iKytLrF69WmzZskVccskl4tJLL/XXx+w1fxyP6dOnizfffFPs3LlTlJaWimuuuUZkZ2d7fAcmT54s7rnnHo/viMVi8dfH9Io/jsnTTz8tRo8e7fF5jx496rGNVr8j/hI04WPixIli7ty57vtOp1NkZGSIhQsX9rj9TTfdJL73ve95PFZQUCDuu+8+IYQQiqKItLQ08fzzz7ufb2pqEkajUbzzzjtCCCF2794tAIjNmze7t/n3v/8tJEkSR44c8dln6ytfH5OefPvttwKAOHTokPuxnJycHv+Cqc0fx2POnDli1qxZZ31PLX9HAvH9mDVrlrjqqqs8HtPq90MI74/J6c72uc63z6amJhERESHef/999zZ79uwRAMSGDRv68Wn6zx/Ho7uGhgYBQKxbt8792OTJk8V//dd/9aVkv/PHMXn66afFuHHjzvo6LX9H/CUopl3sdjtKSkpQVFTkfkyWZRQVFWHDhg09vmbDhg0e2wPA9OnT3dtXVlairq7OYxuz2YyCggL3Nhs2bEB8fDwuvvhi9zZFRUWQZRmbNm3y2efrC38ck55YLBZIkoT4+HiPx5999lkkJibiwgsvxPPPPw+Hw9H3D+MD/jwea9euRUpKCkaMGIEHHngAjY2NHvvQ4nckEN+P+vp6fPrpp7j77rvPeE5r3w+gb8fEF/ssKSlBR0eHxzZ5eXnIzs7u8/v6gj+OR08sFgsAICEhwePxt99+G0lJSRgzZgyKi4vR1tbms/fsK38ek/379yMjIwODBw/G7bffjqqqKvdzWv2O+JPmrmrbk2PHjsHpdCI1NdXj8dTUVOzdu7fH19TV1fW4fV1dnft512Pn2iYlJcXjeb1ej4SEBPc2avHHMemuvb0dTzzxBG699VaPKzM+9NBDuOiii5CQkIBvvvkGxcXFqK2txQsvvNDPT9V3/joeM2bMwI033ojc3FxUVFTgl7/8JWbOnIkNGzZAp9Np9jsSiO/HW2+9BZPJhBtvvNHjcS1+P4C+HRNf7LOurg4Gg+GMAH+uYxsI/jge3SmKgocffhiXXXYZxowZ4378tttuQ05ODjIyMrB9+3Y88cQTKCsrwwcffOCT9+0rfx2TgoICLFmyBCNGjEBtbS0WLFiAK664Ajt37oTJZNLsd8SfgiJ8UOB1dHTgpptughACixcv9nhu3rx57j/n5+fDYDDgvvvuw8KFC0PuWga33HKL+89jx45Ffn4+hgwZgrVr12Lq1KkqVqa+N954A7fffjsiIyM9Hg+n7wed29y5c7Fz5058/fXXHo/fe++97j+PHTsW6enpmDp1KioqKjBkyJBAl+l3M2fOdP85Pz8fBQUFyMnJwXvvvdfjyGE4CIppl6SkJOh0ujNW/tbX1yMtLa3H16SlpZ1ze9d/z7dNQ0ODx/MOhwPHjx8/6/sGij+OiYsreBw6dAgrV670GPXoSUFBARwOBw4ePOj9B/ERfx6P0w0ePBhJSUkoLy9370OL3xF/H4+vvvoKZWVl+MlPfnLeWrTw/QD6dkx8sc+0tDTY7XY0NTX57H19wR/H43QPPvggPvnkE3zxxRfIzMw857YFBQUA4P57pRZ/HxOX+Ph4DB8+3OPfES1+R/wpKMKHwWDA+PHjsXr1avdjiqJg9erVKCws7PE1hYWFHtsDwMqVK93b5+bmIi0tzWMbq9WKTZs2ubcpLCxEU1MTSkpK3NusWbMGiqK4/7KoxR/HBDgVPPbv349Vq1YhMTHxvLWUlpZCluUzph8CyV/Ho7vDhw+jsbER6enp7n1o8Tvi7+Px+uuvY/z48Rg3btx5a9HC9wPo2zHxxT7Hjx+PiIgIj23KyspQVVXV5/f1BX8cD6CzhcGDDz6I5cuXY82aNcjNzT3va0pLSwHA/fdKLf46Jt21tLSgoqLC/Xm1+h3xK7VXvPbWsmXLhNFoFEuWLBG7d+8W9957r4iPjxd1dXVCCCF+9KMfiV/84hfu7devXy/0er34n//5H7Fnzx7x9NNP93iqbXx8vPjoo4/E9u3bxaxZs3o81fbCCy8UmzZtEl9//bUYNmyYJk6jFML3x8Rut4vrr79eZGZmitLSUo/Twmw2mxBCiG+++Ua8+OKLorS0VFRUVIi///3vIjk5Wdxxxx2BPwDd+Pp4NDc3i8cee0xs2LBBVFZWilWrVomLLrpIDBs2TLS3t7v3o9XviD/+zgghhMViEdHR0WLx4sVnvKeWvx9CeH9MbDab2Lp1q9i6datIT08Xjz32mNi6davYv39/r/cpROdplNnZ2WLNmjViy5YtorCwUBQWFgbug5+FP47HAw88IMxms1i7dq3HvyFtbW1CCCHKy8vFM888I7Zs2SIqKyvFRx99JAYPHiwmTZoU2A9/Fv44Jo8++qhYu3atqKysFOvXrxdFRUUiKSlJNDQ0uLfR6nfEX4ImfAghxMsvvyyys7OFwWAQEydOFBs3bnQ/N3nyZDFnzhyP7d977z0xfPhwYTAYxOjRo8Wnn37q8byiKOLJJ58Uqampwmg0iqlTp4qysjKPbRobG8Wtt94qYmNjRVxcnLjrrrtEc3Oz3z6jt3x5TCorKwWAHm9ffPGFEEKIkpISUVBQIMxms4iMjBQjR44Uv/vd7zx+GKvJl8ejra1NTJs2TSQnJ4uIiAiRk5Mj7rnnHo8fKkJo+zvi678zQgjx2muviaioKNHU1HTGc1r/fgjh3TE529+JyZMn93qfQghx8uRJ8dOf/lQMGDBAREdHi+9///uitrbWnx+z13x9PM72b8ibb74phBCiqqpKTJo0SSQkJAij0SiGDh0qfv7zn2umz4cQvj8mN998s0hPTxcGg0EMHDhQ3HzzzaK8vNzjPbX8HfEHSQghAjTIQkRERBQcaz6IiIgodDB8EBERUUAxfBAREVFAMXwQERFRQDF8EBERUUAxfBAREVFAMXwQERFRQDF8EBERUUAxfBAREVFAMXwQERFRQDF8EBERUUD9f/4zsjovgtHwAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#multiply by average content, EFSA Summary, Tuna fish 3.04 mg/kg or 3040 ng/g\n", "#best fit to the data for q0=100 ng/g of fish, a level of 30 below tunna content, \n", "#a level of 2 below fish muscle\n", "q0=0.05e3\n", "thg=getDataScaled('H_THg',q0)\n", "\n", "#5year data\n", "humanHG=181\n", "humanHGplusI=179\n", "\n", "data={}\n", "data['humanHG']=downloadKey(humanHG,0)\n", "data['humanHGplusI']=downloadKey(humanHGplusI,0)\n", "for x in data:\n", " data[x]['solution']=runSolver.loadSolutionFromDir(data[x]['localDir'],True)\n", "volume='hair'\n", "model=cModel.model()\n", "Q=data['humanHGplusI']['solution']\n", "#sys.parse(os.path.join(fh,'software','src','Integra','models','cDiazepam.json'))\n", "setupFile=Q['setup']\n", "modelFile=Q['model']\n", "parameterFile=Q['parameters']\n", "print('modelFile: {} {}'.format(modelFile,os.path.isfile(modelFile)))\n", "print('parameterFile: {} {}'.format(parameterFile,os.path.isfile(parameterFile)))\n", "\n", "model.parse(modelFile,parameterFile)\n", "pars=model.parSetup['parameters']\n", "lutSE=Q['lutSE']\n", "lut=Q['lut']\n", "y0=Q['sol'][-1][lut[volume]]\n", "dydp=Q['sOut'][-1][lut[volume]]\n", "sigmaAcc=0\n", "cvLN=[]\n", "dydpLN=[]\n", "for p in lutSE:\n", " par=pars[p]\n", " try:\n", " par['dist']\n", " except KeyError:\n", " continue\n", " j=lutSE[p]\n", " if par['dist']=='normal':\n", " sigmaAcc+=dydp[j]*dydp[j]\n", " continue\n", " cv=par['cv']\n", " dydpLN.append(dydp[j]/cv)\n", " cvLN.append(cv)\n", " print('{:40} {:.2g}/{:.2g}'.format(p,dydp[j],cv))\n", "print('sigmaAcc={}, dydpLN={} cvLN={}'.format(sigmaAcc,dydpLN,cvLN))\n", "fx=numpy.linspace(0,5*y0,101)\n", "cvLN=numpy.array(cvLN)\n", "y=propagateErrorLN.calculateDistribution(fx,y0,numpy.ones(cvLN.shape),cvLN,dydpLN)\n", "#convolve with zero mean gaussian with a sigma=sqrt(sigmaAcc)\n", "fx1=numpy.concatenate((-numpy.flip(fx[1:]),fx))\n", "g=numpy.exp(-0.5*fx1*fx1/sigmaAcc)\n", "y=numpy.concatenate((numpy.zeros((len(fx)-1)),y))\n", "z=numpy.convolve(y,g,'same')\n", "print('{}/{}'.format(len(g),len(y)))\n", "z=z[-len(fx):]\n", "\n", "#step=(fx[-1]-fx[0])/(len(fx)-1)\n", "z=z*len(thg)/numpy.sum(z)\n", "matplotlib.pyplot.plot(fx,z)\n", "matplotlib.pyplot.text(0.5*fx[-1],0.5*numpy.max(z),volume) \n", "#fx=numpy.linspace(0,50,100,endpoint=True)\n", "matplotlib.pyplot.hist(thg,fx)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.6" } }, "nbformat": 4, "nbformat_minor": 4 }