{
"cells": [
{
"cell_type": "markdown",
"id": "500b07b7-5f43-40c0-ba80-bc6cd759f9f4",
"metadata": {},
"source": [
"# Fl\u00e4chenkartierung von Objekten in Kacheln\n",
"\n",
"In diesem Notebook werden wir Zellkerne in Kacheln segmentieren und ihre Fl\u00e4che messen. Anschlie\u00dfend speichern wir die resultierende Fl\u00e4chenkarte wieder als Kacheln in einer Zarr-Datei. Diese Strategie kann verwendet werden, um Daten zu verarbeiten, die als Ganzes nicht in den Arbeitsspeicher des Computers passen."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "e6a9300d-1f11-4a3b-94bb-a136ba69f09d",
"metadata": {},
"outputs": [],
"source": [
"import zarr\n",
"import dask.array as da\n",
"import numpy as np\n",
"from skimage.io import imread\n",
"import pyclesperanto_prototype as cle\n",
"from pyclesperanto_prototype import imshow\n",
"from numcodecs import Blosc"
]
},
{
"cell_type": "markdown",
"id": "8959f8d4-a6d6-4a2d-b4b7-9378d2ceec01",
"metadata": {},
"source": [
"Zu Demonstrationszwecken verwenden wir einen Datensatz, der von Theresa Suckert, OncoRay, Universit\u00e4tsklinikum Carl Gustav Carus, TU Dresden, zur Verf\u00fcgung gestellt wurde. Der Datensatz ist unter der [Lizenz: CC-BY 4.0](https://creativecommons.org/licenses/by/4.0/) lizenziert. Wir verwenden hier eine zugeschnittene Version, die als 8-Bit-Bild neu gespeichert wurde, um sie mit dem Notebook bereitstellen zu k\u00f6nnen. Das vollst\u00e4ndige 16-Bit-Bild im CZI-Dateiformat finden Sie [online](https://zenodo.org/record/4276076#.YX1F-55BxaQ)."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "cc2eeeb8-eb5e-49fc-8569-cdff5e143e5e",
"metadata": {},
"outputs": [],
"source": [
"image = imread('../../data/P1_H_C3H_M004_17-cropped.tif')[1]\n",
"\n",
"# Zu Testzwecken beschneiden wir das Bild noch mehr.\n",
"# Kommentieren Sie die folgende Zeile aus, um es auf den gesamten 5000x2000 Pixeln auszuf\u00fchren\n",
"image = image[1000:1500, 1000:1500]\n",
"\n",
"# Komprimieren UND das Numpy-Array in ein Zarr-Array umwandeln\n",
"compressor = Blosc(cname='zstd', clevel=3, shuffle=Blosc.BITSHUFFLE)\n",
"\n",
"# Bild in Zarr-Array konvertieren\n",
"chunk_size = (100, 100)\n",
"zarray = zarr.array(image, chunks=chunk_size, compressor=compressor)\n",
"\n",
"# Zarr auf Festplatte speichern\n",
"zarr_filename = '../../data/P1_H_C3H_M004_17-cropped.zarr'\n",
"zarr.convenience.save(zarr_filename, zarray)"
]
},
{
"cell_type": "markdown",
"id": "d76246fe-7358-4e0c-8112-1f1fd0af4108",
"metadata": {},
"source": [
"## Objektfl\u00e4chenkarten in Kacheln\n",
"Dask bietet integrierte Unterst\u00fctzung f\u00fcr das Zarr-Dateiformat. Wir k\u00f6nnen Dask-Arrays direkt aus einer Zarr-Datei erstellen."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "2132d10e-1ec5-43eb-9c3c-a4d9358919cc",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"
\n",
"
\n",
"
\n",
" \n",
"
\n",
"
\n",
"
Array
\n",
"
Chunk
\n",
"
\n",
" \n",
" \n",
" \n",
"
\n",
"
Bytes
\n",
"
244.14 kiB
\n",
"
9.77 kiB
\n",
"
\n",
" \n",
"
\n",
"
Shape
\n",
"
(500, 500)
\n",
"
(100, 100)
\n",
"
\n",
"
\n",
"
Count
\n",
"
26 Tasks
\n",
"
25 Chunks
\n",
"
\n",
"
\n",
"
Type
\n",
"
uint8
\n",
"
numpy.ndarray
\n",
"
\n",
" \n",
"
\n",
"
\n",
"
\n",
" \n",
"
\n",
"
\n",
"
"
],
"text/plain": [
"dask.array"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"zarr_image = da.from_zarr(zarr_filename)\n",
"zarr_image"
]
},
{
"cell_type": "markdown",
"id": "c2721aa7-947e-4855-9325-c3e2b4746226",
"metadata": {},
"source": [
"Wir k\u00f6nnen die Bildverarbeitung direkt auf diesen gekachelten Datensatz anwenden."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "cba0b2e9-c7ac-43dc-b5b3-b0aadb43b425",
"metadata": {},
"outputs": [],
"source": [
"def area_map(image):\n",
" \"\"\"\n",
" Objekte in einem Bin\u00e4rbild markieren und ein Pixelz\u00e4hlkarten-Bild erstellen.\n",
" \"\"\"\n",
" print(\"Verarbeite Bild der Gr\u00f6\u00dfe\", image.shape)\n",
" \n",
" labels = cle.voronoi_otsu_labeling(image, spot_sigma=3.5)\n",
" result = cle.pixel_count_map(labels)\n",
" \n",
" print(result.shape)\n",
" \n",
" return np.asarray(result)"
]
},
{
"cell_type": "markdown",
"id": "40c716df-6dbf-4b94-a1e6-a090d142a395",
"metadata": {},
"source": [
"## Testen der gekachelten Bildverarbeitung\n",
"Wir sollten unseren Fl\u00e4chenkartierungsalgorithmus an einer einzelnen Kachel testen. In einem realen Szenario wird der Bildverarbeitungs-Workflow tats\u00e4chlich auf einzelnen Kacheln entwickelt, z.B. in einem Notebook wie diesem. Sobald wir sicher sind, dass der Algorithmus funktioniert, k\u00f6nnen wir ihn auf alle Kacheln anwenden."
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "697eb25f-e546-48e2-a744-a63e8db189bb",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAPsAAAD7CAYAAACscuKmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAkZUlEQVR4nO2dW4xd53Xf/2uf29yHHFKkaFI25Uq1rDpwHaiNEgeFYCVt6hoR0MKGU7hQXQd6SRPHDRBL7UPQhwJ+MIL4oS1AOA2MxGhsOEYlGEEuYKKHAK1gWnbiSLRkSaYkkiNyRM4M5z7nsvqw1tp7nzNnyBnOmTPnzP7/AOJwX883++y9/9+3vnURVQUh5PCTHHQDCCH9gQ87IQWBDzshBYEPOyEFgQ87IQWBDzshBWFPD7uI/JKIvCIir4nI071qFCGk98jdzrOLSAnAqwB+EcBlAN8F8Cuq+nLvmkcI6RXlPRz7TwG8pqpvAICI/DGAJwBs+7BXZURHZRwQsRWJbN0p3j3xEop9fTn/cupydNsx2mrlVsmW4w8F0nEVDtvf1yOkVLL/xD2Xv25xzdLPOMgXm83cvvvWxJ6wjhVs6kbXR2MvD/tpAG/nli8D+JnOnUTkKQBPAcCIjOPRkY8DlYptq1WzHeNC+4XVesP2qZTblnWznp27tM0oxH/Y1spKuioZGbF16+s7+NOGB6lU25a1vnlALRlQ/KEuTU3Z4uiorS/nbv2G31v1etty7NNaWEx31dg2oLyg57fdtpeHvdvbY8t7T1XPATgHAFPJMdVmC1IyxdW17MGTkZq3yJvU9H38M14C6X4AdG3N1tV8XdL+8Lc9CPFm7+gpDAU5FQqFips2mZoEALRuLQHoeNiH8W/tNf63N+OB9c94+QOZAKT3i/q9N+AP9m7Zi4HuMoD7cstnAFzdW3MIIfvFXpT9uwAeFJH7AVwB8GkA//a2R6ia8vibE9LlXdOp4KFkE+N2CldzIDe8z3XtAUCqlfZzAGgtL9u6sm0bqu5uTplDbULrWzdutq3f7jjSTrfh3FDdE3fBXT/sqtoQkf8I4M8BlAD8L1V9qWctI4T0lL0oO1T1TwH8aY/aQgjZR/b0sN8taVc0Z0CLrngy7tbS6L6HIaXsBrZ8Nz7OU27/Mzq79UBm2Br6rpob3WJYcruuelyXw2ZoGiSGaTaE7rKEFIQDUfag2zRRy6fGk6kJW87NcQId8+wdb9V0vRvmUvU7TOzC6EZF3yN+T4ZRF8juy+hx6hGb+oTflzpu65P5pfSY5tVrtq1h+6S9zD7/PlR2QgrCgSp7N0LtW4u3bDmm4vztGtNqANBaXbV18aYM77vletsx+fMScidKx2bsPzNHAADNYxPpts1Re2Sk4c46o3bvrdxr91ppw9ZPvZ7paBL37C3vabqiN2/O23Kfpkip7IQUhMFRdn+7bbEgpy6f5ogTag5sVfTOc1HNyU6Iey6ZNv/51ntPAQAWP2jj8fJ6prwtnxQauWn3Z7Jp92Vt0e7B6i33s88Feemo2ZYE3kPw+7MUvvfuGAXs7zieyk5IQRgcZXe2vNk6xzP55dTdttl9X0J2QKeiz/+ULa/PmDpXF7P7amTBlDwUXUu2z8jcxrbn16or+JgpfHPEl6vH7NgXs55p890be/lTbguVnZCCMHDKvh3dxjIck5O7JRkfzxbc6r74kI3RG+60OXHFFLe8lil7Zdm9Nj30Oqn7tpbbiSpuRypnOqpud6pPtD9uJe8dpOHX+wyVnZCCwIedkIIwNN14QnpJpDsDgMY90X13g9yydckTt5vVbmbGt2TTHbd8ak3qbqirWVc82XCHmdGc81e1XVPLq3aO0nKXYeg+ZheishNSEKjspJjkjGLSMHWuLttnGOTCGJfPAdiYsOmzVOHdEFdaNRft5njV12fHlFc8fHu90XZMsuTu3tXuAV29hspOSEGgspNCouvZOFx82ixcXlvuKBNj7fyYe+14uW3f2g07T6j12kkLr47xPwBMvmk9Btn0abu6b/MAGc31MpKxMfvOXBr0XkFlJ6QgDL6ys+IJ2Qfy1vi4o+rjprCps4tb46WZq0Lk60rrHpjlTjTwz0bN7tfVE5mOTr/qY393xGkcM4eeVtW/bz3nMPZqFujVa6jshBSEgVH2CDOM6i5p2Sd/G0apKN3I5iZj3BXpfqj6ZKeklYaQjaVrC+2JSlsV08K8so9d8+QqPo6PQJjYZ3zWtk9czuX69zn5xnELcd2csrF6c8TOMXXherpvYx/vYSo7IQWBDzshBeFguvFudItpBgBIZo4CALTmueY2Io+cG0Aanouulivs6FVctZmtAwB4Nc7DVrGV9I781Ftp0aa5Kl4YNKbRovvdHMtcXyNyrbxm92Oz5k41GxHn3mz7zNP0rn953baNzFkNhMblK3v+e3YClZ2QgtB/ZRdJs3emtbIB6KYb3qZsWmL1H1gWjzCSjF22zJzJ3EJ2qji+5caWjjrbbRVnGPtOcuTvB/VMxmneoxmv5e7LeQNdKHrJ1RkjHgAT7rPe22yOZjpaWqv7p92fpVvW42y9+oY3gNllCSE9pK/KLkmCZGIiHXfrSE55j9nbdPWMTU/cel+4JXr2zps+lh/Nxuey5uMuV/Io95wGOeSnV1zlqfCkk6ZXHUp8HB/OqzphNqXGZHaflpfa759yTA17ME192u7bVPkBJEum5LJqn423L/ey+TuGyk5IQejvmD0RSLUCHfMkX0n2ronQwFbFRkoVTyAQObvLC2a5TNUcgG50ZPSMaq7xts1Vj0lrxO1jcgAy3KSzN+9aHndZMdfV2mo2q6MxOxT3lt9HkTm2sui9g7msRmHr2hwAoHnAs0NUdkIKQn+VvVwGjs+k1nNZzpz+yz7eHl+1cdKo19QqXzdLaSh6Xs2l4m/XCBV0F1uZt7dqvuJrWhXW621x7E62Iw0v9U/JVRKO6kNJrd23I61O5PdnYwAr6FLZCSkI/VV2EWitDFkzVW1TaR8LlW6akpdiPO+ec7GvjGdedzpu8+wxXoLPh0bVzGRlLdt3uffJAEgx6FazYBi9M++o7CJyn4j8tYhcFJGXROTzvn5GRP5SRH7sn0f3v7mEkLtlJ934BoDfUtUPAngUwK+JyMMAngZwXlUfBHDelwkhA8odu/GqOgtg1v+/JCIXAZwG8ASAx3y3rwF4HsAXb3uyVguyugHx7nXb5Jd311HuKIUTy0e8+N5E5mK7dspca5fusz9j5Yytn3nJznzkxSxOWEZ8us8NdFtKQxNyyNmVgU5EzgL4CIAXAJz0F0G8EE5sc8xTInJBRC5sNte67UII6QM7NtCJyASAPwHwm6p6Szpzw22Dqp4DcA4ApmsnVTbq6VSZ5JxqUoUNQ1zHtFp8NiazKY/l91jz5z9sU3nTZ2yKZHnRAm2OfC+XLaTDASemUAgpCjtSdhGpwB70r6vqt331NRE55dtPAbi+3fGEkIPnjsouJuG/D+Ciqv5ubtNzAJ4E8CX/fPaO3yYCrVUgdVfxVmvrPh1urDriSu6hg1FLCwAqq7bv2Ns2rl971xT99A9cxRduZefx3HWpC616jS6O2UlB2Ek3/qMA/h2AH4rID3zdf4Y95N8Ukc8BeAvAJ/elhYSQnrATa/zfIIvj7+TxXX2bCFCrQj1TbN5dFh5emNa9iqCWUGdX4PJS5lRzxBMGTFyu+TY7RyS40Hwvod6eOZSKTooG3WUJKQh9dZdt1UpYPTsFaZjijr6dU97IAR/jeE9Tpas+Xefr22p0+T6V2fbvSVU7P2PQEajAUFdSNKjshBSEviq7JoL6WILxqz4+X1zONsace4SrpmN4t56HEufG2uEFF152GvP3Htqa7wWkKavEv0fdBnD6PXbaK1f38JcRMvhQ2QkpCHzYCSkIfe3Gl1brmP7+9XQ6rW1qLLrknvcrDGeRbSYKO6b7AWiesKjalfstC436q2viksWuJ5c6LHcAJKrI+Ewcu++kKFDZCSkI/c1Us1mHXp7NyjCXtn/XREaaUPQIjInsNACwctYUfeEBM9BVF6w3MPmaT9PlyjtvmXojpGBQ2QkpCP2delNtz92lua/3KbGYakvH6mOm5OrZYRvHxtNDGqNeDdbH35NXbDotuWE2gVbORXYYc4YR0kuo7IQUhIOpz+50C0bRhiu6j7FjXC+eTqo0mtXdGrtmzZ980xxmKm+8AwBozL6zTy0mZHihshNSEA5U2bsS8+uu+lFhM5JUJpuZhb36ertbbINVXgjZFio7IQWBDzshBWHwuvHbEIUYmzduHnBLCBlOqOyEFAQ+7IQUBD7shBQEPuyEFAQ+7IQUBD7shBQEPuyEFAQ+7IQUBD7shBSEofGgI+SwIrkkqhEAloyM2AoP9U7rHjh3k4yFyk5IQeDDTkhBYDeekH7jpcxKM1b3IC1TDmwpgxYFTdHStmVZWUkPiXwOd+raU9kJKQhUdkL6hFRMwUvHZ2xFZFIeqWU7ec7F5oStWztl2ZXrY7a+tmCGutErS9l5r1nYtzQawNa0jilUdkIKwo6VXURKAC4AuKKqnxCRGQDfAHAWwCUAn1LV+f1oJCHDTKro955oW6+VrY9fKPqtB6w+wvV/4ufwofv0K2VfnkiPGZu3zMtSLkMasm07dqPsnwdwMbf8NIDzqvoggPO+TAgZUHak7CJyBsC/AvDfAPwnX/0EgMf8/18D8DyAL/a2eYQMPzFG1zF3lHEHGdmw7MitY6NbjhH3oRmbNT2Oqke1RZP4ynKXwXmplFr6u7FTZf89AL8NoJVbd1JVZwHAP090OQ4i8pSIXBCRC3WwqCIhB8UdlV1EPgHguqp+T0Qe2+0XqOo5AOcAYEpm9A67E3IoSN1dkas83KHouma1EGRzLN23VLd9Jt80hZ563fS1OWqW+2TDFL20uJadf8kt88nttXsn3fiPAvhlEfk4gBEAUyLyRwCuicgpVZ0VkVMAru/gXISQA+KO3XhVfUZVz6jqWQCfBvBXqvoZAM8BeNJ3exLAs/vWSkLIntmLU82XAHxTRD4H4C0An+xNkwgZfmRiPLfgRjMvbQaPbJOaTbNFeXIAkKZ120uxbyxfW/B9N9vO0fadpeS2BrpdPeyq+jzM6g5VvQHg8d0cTwg5OOguS0gvCWWNOPQcOlZr20cabrBrtnL7uGEvtq1acMvtFD2lXAa2F3a6yxJSFKjshPQQKXto6mQ2Zm8eM9fWVsXUvtL0suSxg+ZmpDdcwSPENfFeQGSz8ek1refG+ZOTtm5pKXfSrVDZCSkIVHZCeojkw1WdzWkLhFk/Zo/b9JqpcinG46HiyI3Ra3ZM4+QRAECy6WP4ZXemWco5s4Z1f3QUWO5NIAwhZIihshPSQ8STT6BW3bJt9F13db1h7q1adot9KdPc5vFpAMDSAzbOXzll+4zOmZJP/8hTWm1kY/ZIVaWNBsfshBAqOyE9JZI/Jsur6brRS+1z5mnSSN8XuXG+Ttk8+/IZU/SVM7ZvddEUPVn1Y3Lz7RrnSxLOsxNC+LATUhjYjSekh0T5pubVa+m6ZMqMbdrhSitjWzPUJOt2/MzL5lwzdcn2nXh9wY6JoJlc1z/cblEuAbK9flPZCSkIVHZCekgoe1uxRldjiek4z1ijq56pJq/SK2bEG5l1p5nIPhNusuGAk5uuQ5y/mc8atxUqOyEFgcpOyD6guRLLiauxjNoYXWu+vGbTaLqSTdMJPB9dVIvxHK/q43ydtKm50txi9l3LK9l3dpR2zkNlJ6QgUNkJ2Q9yYautW1axJYlxvCei6BriGtVa37XiSomP51vTFjLb6uKGCz+vrm+0n6sDKjshBYHKTsg+o3VX8sVbAACZnrINMb7Oj7PDDTbWudJLwxNP3nKX20Z2TIzZWysrUN3eIk9lJ6Qg8GEnpCAcym58lMiVnONBa339oJpDCIDcPbjhkXFjPs2W68a3VszRJvGSUbpq03KSj1kH0PIhQdt57wCVnZCCMHTK3uaGGIXyIkd3BAG4kUJzdo9Qe214IMFtpigI2Vf83mutmGGtrYqLb2suugtsZKv16bsw9t0NVHZCCsLQKHuUwM3X0Iophy2U/G3YzXVQ3e3wdpU1COkn3XqZvm4vSt4JlZ2QgjDwyr5F0fNvwYoreKu7I0Fb6OC4Ha9uuYxcYeqWUSo9OexQ2QkpCAOv7DGHmITFPT8O76yUGUkBPIi/bXw/bamBYr4ynae8TRofMmBE9dPO2Zcc6Ri3syY6obITUhQGXtljXj1N7ZOri5WStL+zki71tnD9hp0n3vSdqkAlGDgSt7MkM0dtRYddRde87lm+h9cxA9PKJ3YACv37UtkJKQh82AkpCDvqxovIEQBfBfAhWIKN/wDgFQDfAHAWwCUAn1LV+V43MM3WWWnPzAkgy7/drWuPbZxuvMsfQTLplFuBu3eDQAzXkmMz2TrPq67+W0liQ62YNo0MLVLJfn/1LDCxLal5ppcbN217gadYd6rsXwHwZ6r6EIAPA7gI4GkA51X1QQDnfZkQMqDcUdlFZArAPwPw7wFAVTcBbIrIEwAe892+BuB5AF/seQtjuqWLekul3LbPFvKGu3q96y5FftMPAtFjS6YnAQB66ni6rTFmv3myZr9R4mWK0187qqPkf+e4F6IHGD1Dz+wKN+oV8XffibK/H8AcgD8Qke+LyFdFZBzASVWdBQD/PNHtYBF5SkQuiMiFOjZ61nBCyO7YyZi9DOCnAfy6qr4gIl/BLrrsqnoOwDkAmJKZux4Yh1trknOUiXFZOt6Ot3mM5duSV7hLbUfVjDT0tYcBB2QHuAJHkgacNEVfO539vpsT9juOzJsKV+f80LLbW/x3TssYI1P9VLnjXvD7B6Hw69kxRfntd6LslwFcVtUXfPlbsIf/moicAgD/vL4/TSSE9II7KruqviMib4vIB1T1FQCPA3jZ/z0J4Ev++ey+tNBVO01akU9eseFv5BiPh1W2S82rUPmWv+HD5VJvU0GD7B9hJZdJG6tvzliKpsZopj9l740ldXdxHnUHq5L1xjan7feuLma51MOyI57sIa2RFuu7/N5FSWiyUw+6XwfwdRGpAngDwGdhvYJvisjnALwF4JP700RCSC/Y0cOuqj8A8EiXTY/3tDW3IRQ5uY1rZKfFXfMKH/m0IynAAVljU38BZCG46XyyV/VM8b/nMCbLTH+bso+76/Zb1uaz36W0butKSx6W7L/9yv0W1LR8ypaPvZT7LeOeSHuCPmb3T/UZnPz8Tcmvf+vmgu1zSMfw9KAjpCDwYSekIAx81FtKZ0ZOtGea7XpIrpsvnbHvfSIit6J7KeNj2cZw9/RupFR9GjBcPhvuVJJvuxucWktL+9XkvhBGMb1lf0fJDajla7mhV2Qg8hJIrWNWNqnll2Ni1q5pdTbLoQ4vfNg84pmJyu16tnHUhk7l9ezeqL3hjjdRJpndeELIMDM8yt6F3RjZ+m2QK+UCOgBAPJdem2unG41ak6Nt25JF772EO2i+B+M9hFJp2hYXFnvY6v4TQS2yYOqsXXpgqSFzzRR34k2rklJaiaCX7JjmpO3bHLFr1qraNV0/astrx215/Fr2O1Tf9im9jcPt4UllJ6QgDLWyDyJpnTkff4cq66SN1ZsTWRYdrfi004YHeiy1l+PVUVe0ek7tfPwa4/rEnVKGbgwfNhivZRaqI6Mj6S5p6GqULfaeTnmuY4o1FyQVit4Yt2u2MWWfy/fZN9Ru2LlG3s3OIWuHW9EDKjshBYHK3gO6OcqkVnO3IC89YAp84x9lKl1yQTn+Q1OZmidnKLkVWrqFcG5utH9PMGzKHoQ7dCh6zmVV4zqEvSXG5hHEEtdlNLsW1Tmzd9QnzKaxOdke/lzesPOXVjNlT2sJHHLXaSo7IQWByt4D8iqbHHcrvKtS0xMwrE+7Vfh9uTncTR9H3vQUSnU7T7LuSRqWV7d+WahZpNdy28DQhup6qGukEMvn+k9dX6vV9mM6rkE3KsvWG4gx/PpxU/TqgofWrueV3XsKhzwQhspOSEGgsveCXK25CGYJtS8v2PL0JbvUzdz4UpqmJGNzdnxpzb3sfD5Zpya2nF9W3WIfterq7fn0h07ZQ00j5Hgja7/UYkbDt0Uykqjus+Q9n82cSo/Z9W2M2jGldb/GV03Rpy/Z+eXqXHpMa7N7yrLDBpWdkILAh52QgsBu/B5IA3FyhiKJ/0e5Kv8cWbHu/ImVo+m+YTyK7nxjzGOvT9h0XQwBJDdNlGZe6XD/TQ1bK11y5Q8B3XIMimezCacZHXVjZAxzImBoLcsDIGJTbmNvmhvxyJwfs+nBLtcsf3xrMZuqHLqhz11CZSekIFDZ90Aoe2pIQs4RJNxkI2+5u3wm65kir500R5Jm1aefkvbP8rwZoDSXJTfLme7LoX5DHsQRDi1az65PlNduTnsmmZpd0+q8qXL8zfmMRDJvih554kvzFmATBR5b6TTb1jyFhx0qOyEFgcq+B1LV7pIXr9M9Q2tbK9qU10xd5j9gPYO6D1fv+dv2o5szE9l5PHimHCGaHsSRr3c2lHRLThLTl7NeQjBpd8BJFT2Xe1A7K//E1F6sd0VnRRhCyKGFyr4TonpJVASN3PNhje+mEmnwho/rYyy6kmWKlaaNK8sdBvSwzmvVjt04loV9Lp+2dfdEOOySHzzkY/ZuNL3yqizauDu93t6Til5At0Ck1DnHx/2HMUPvbqGyE1IQqOy3obPCKI74/HeMFcOy260+fCSZiNDTCFzJzcnXLpvl+ORNV6OKp6WKQBifT67eypRrejPCPtutybp++JQ9iPF1Os6OSq0RRNPIjdmXtxmzEyo7IUWByu6kKh5VRQHg5D0AgPX3HQEALJ02BRcX1SOvmedWVCwBgGTJvbkiVXKM1UsddemQWdBLnckYGu1JFCqXb6T/L3coePQUmodwzL4tVOu7gspOSEHgw05IQWA33omSztF1B4Clh48BAN75Gds29tACAKDesKmf9f9rBrsTL+bKDPu0UOmmTRclU2bcS5098l3QmErqiKdOAzy65USLuO8Vd6UNZxp2bckdoLITUhCo7E7UYIugCwBYus+dN86Y0e2pf/g3AICHalcBAL+68FkAQHUhyz4zfcmng3warXRz2bf4eSWX7XTdQzUjRDMccLx6TOvW1oyx4ZobziPDXhGG9A8qOyEFgcoeRJXVeuasUl61cXDlDVPaL6//CwDA2FFXYjWVXn5vptbr95jijs3aVN7Ri7a+NO9urV3G1hGOmU69RZhsqcu72NtJRSe7hcpOSEHYkbKLyBcA/CoscvOHAD4LYAzANwCcBXAJwKdUdX5fWtkH0mqiOeWtLbrKvxX53T0UddI+K6O27+bZzKnmYw+/DAC4vHoEAHD1D+8HABz/np0ryWU1TZW8o8586h7qbrmaS7tUxNBM0hvuqOwichrAbwB4RFU/BKAE4NMAngZwXlUfBHDelwkhA8pOx+xlAKMiUocp+lUAzwB4zLd/DcDzAL7Y4/b1jVDMZDNTztJm+/i65HXCpBVjdPs8cjyzmv+P0/8PAPB3m6b2/+aBLwAAaouWgGL0SGa5Ly+ZNb505V0AudRJw1q3jQw0d1R2Vb0C4MsA3gIwC2BRVf8CwElVnfV9ZgGc6Ha8iDwlIhdE5EIdBfLfJmTA2Ek3/iiAJwDcD+A9AMZF5DM7/QJVPaeqj6jqIxXU7nwAIWRf2Ek3/hcA/ERV5wBARL4N4OcAXBORU6o6KyKnAFzfx3buO5GlVday3sfKvebAsnqvdddH3vV93Ys1Si7PvZXlgv/50r8GAFxfsG57ZdmO3ZjyeHbNYt8rVV83c8bO/+aCbXiF3XjSe3Yy9fYWgEdFZExEBMDjAC4CeA7Ak77PkwCe3Z8mEkJ6wR2VXVVfEJFvAXgRQAPA9wGcAzAB4Jsi8jnYC+GT+9nQfadL0MnkZTPW1SdMjZOG5373OBUtmWqPvZ1dxrkbJwEAVVf02k07pryx1ZmmOeLKHl/dzYmGkB6xI2u8qv4OgN/pWL0BU3lCyBBAd1knso/q5dl03XiEnrZOAwDqE6a86zP2WV7xqbhm5i47dcmV3HPCRwmXxKfxFt+fXfJRL9U8+ZY75dxY6MWfQkhX2G8kpCBQ2TvIV/RsXrMJhrEXbezevP9eAICoBa7U5mNMn1WEiXF9dd7O06qWfB+71M0sBXyay65yyb6ncW2oJzTIgENlJ6QgUNlvQ7jQhsInqxaQMvlDD2o5YrXAKyeyefZk1Sffy54+ymuzVW7aMeM/yVnlr7wDAGgwXJX0ASo7IQWByr4LOgNU0oqjV65m6yK1VEd99qg/fleJIfOprJhY8rZE/v+87YUYVHZCCgIfdkIKArvxPWZLEcLd0FEaOjLIRq54IMuoo+Hwo+0FHru15bAjuUw/kbcvGTfjKXP1ZVDZCSkIVPYDJBQpsstKzYxLUfAxqsik2WeBVP1T1XfDX1pAMlcSGp677rArvJS33saHuYT13UJlJ6QgUNn7TDI+nv3f68DBVVrDEcen12TN3sX57LKh+pGvLsaoMmp+uLqRTTlFQg6In+eQTkdFEBOArOeTs3MQg8pOSEGgsveJZNJUPFxs82jV1XrMrfDhOBMqVc5Z4716azLhPYSk/X0dNesAAA0fz8cMgbv7phb9w6j0fu0Ou53ibqCyE1IQqOx9Ipk5AgDQkSzDro6a9T3CYLXslV9XfQ49VDun3qlyt1zBatYrkA0/ppFLrxU9gmlLfqn3WbbvyFPfmH1nL38SGTKo7IQUBD7shBQEduP3mTDMacUudRjhAKA54U40rchl56WaNzuMS+EaC0B9mkmqfmxiRrzozqOa5aVvjdk+K+91Y57b/SYXzcjXGZlHDjdUdkIKApV9n4i46nCBDcWVXH76UPIwyEndt7nraxjw5FYu2CV1j/VjwuXWDXSticy1du1eM+atHbV3+uj89kEz5PBDZSekIFDZ94lOh5WoISe58Xey6FNjPmZPp8p8+kxiGi03Dsemn7fDHVQnTcVbI9lPWl6z46cv2Wdpzcfm8xb2qV2q4JDDC5WdkIJAZd9vXIF1aRlAFsgCALrSHoaZ3wYAOuZJ5vM14CZMwRvTti1xtU6WzRW2tJSdM1nt6F1s2L6tRc+lx3x2hYLKTkhBoLL3iU7Vzq/TllnJI5Q1TWbh4/vmscn0mBs/Zf+/+WFT5Xsu2PqZ77qi50Ncl1bavk8Xb9lnow5SPKjshBQEKvt+ExbvSJ2UC1eNABeJFEqu6OkcemzPedSpv54rt8xzrrThXnerXoV2MzdOj9DWmAFw+4GUI+3VIQxxJdtCZSekIPBhJ6QgsBu/30SQST5DrKOjHhQTn9HdXrBP9fJSSS6e/fjf2k82c9G65OUb7Ua4fKbV5s0FWxdTd9GNr7IbX0So7IQUBCr7PtNasem0tMqLh6YCgFbt8i98cAoAMHrDegGjF13ZvZCk5gpKljo+U8ebLmGqUup4lydbs9WS4kBlJ6QgiPbRZVJE5gCsAHi3b1+6d45jeNo7TG0Fhqu9w9LW96nqPd029PVhBwARuaCqj/T1S/fAMLV3mNoKDFd7h6mt28FuPCEFgQ87IQXhIB72cwfwnXthmNo7TG0Fhqu9w9TWrvR9zE4IORjYjSekIPBhJ6Qg9O1hF5FfEpFXROQ1EXm6X9+7U0TkPhH5axG5KCIvicjnff2MiPyliPzYP48edFsDESmJyPdF5Du+PMhtPSIi3xKRH/k1/tlBba+IfMHvgb8Xkf8tIiOD2tbd0JeHXURKAP47gH8J4GEAvyIiD/fju3dBA8BvqeoHATwK4Ne8jU8DOK+qDwI478uDwucBXMwtD3JbvwLgz1T1IQAfhrV74NorIqcB/AaAR1T1QzDP5E9jANu6a1R13/8B+FkAf55bfgbAM/347j20+VkAvwjgFQCnfN0pAK8cdNu8LWdgN93HAHzH1w1qW6cA/ARuEM6tH7j2AjgN4G0AM7DYke8A+OeD2Nbd/utXNz4uYHDZ1w0kInIWwEcAvADgpKrOAoB/njjApuX5PQC/DSBf5mVQ2/p+AHMA/sCHHV8VkXEMYHtV9QqALwN4C8AsgEVV/QsMYFt3S78edumybiDn/ERkAsCfAPhNVb110O3phoh8AsB1Vf3eQbdlh5QB/DSA/6mqH4HFRwxkN9jH4k8AuB/AewCMi8hnDrZVvaFfD/tlAPflls8AuNqn794xIlKBPehfV9Vv++prInLKt58CcP2g2pfjowB+WUQuAfhjAB8TkT/CYLYVsN//sqq+4Mvfgj38g9jeXwDwE1WdU9U6gG8D+DkMZlt3Rb8e9u8CeFBE7heRKszg8VyfvntHiIgA+H0AF1X1d3ObngPwpP//SdhY/kBR1WdU9YyqnoVdy79S1c9gANsKAKr6DoC3ReQDvupxAC9jMNv7FoBHRWTM74nHYcbEQWzr7uij4ePjAF4F8DqA/3LQxoou7ft52NDi7wD8wP99HMAxmCHsx/45c9Bt7Wj3Y8gMdAPbVgD/GMAFv77/B8DRQW0vgP8K4EcA/h7AHwKoDWpbd/OP7rKEFAR60BFSEPiwE1IQ+LATUhD4sBNSEPiwE1IQ+LATUhD4sBNSEP4/38X4rnROknUAAAAASUVORK5CYII=\n",
"text/plain": [
"