{"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_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"<div class=\"alert alert-block alert-success\">\nThis is the accompanying Kaggle Notebook to my Medium article \n<a href = https://towardsdatascience.com/5-ideas-to-create-new-features-from-polygons-f8f902f5ad8f>\"5 Ideas to Create New Features from Polygons\"</a></div>\n\n![polygon copy.png](attachment:3810739b-f5b3-4133-a6f2-28ccc72bfb16.png)\n\n*Image by author*\n\nPolygon data can be useful in various applications of data science.\nFor example, to determine the energy usage of a building, you could use the floor plan of the building which could be represented as a polygon.\n\nThese polygons can be represented in **well-known text (WKT) format**. \nThe WKT format is a markup language to represent geometric 2D and 3D objects, such as points, lines, polygons, and so on. \nIn the WKT format, a polygon is represented by the coordinates of each point of the polygon.\nHere is an example of a polygon description in WKT format: \n* `\"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10))\"` \n* `\"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"`\n\nWhile you could parse the polygon coordinates from the WKT string and write the functions to calculate features like the polygon’s area or perimeter yourself, the **Shapely package** [1] does all of this for you out of the box.\nYou can simply load a polygon's WKT string with the Shapely package into a Shapely polygon as follows:","metadata":{},"attachments":{"3810739b-f5b3-4133-a6f2-28ccc72bfb16.png":{"image/png":"iVBORw0KGgoAAAANSUhEUgAAAzQAAAE4CAIAAAACTSFyAAAAAXNSR0IArs4c6QAAAERlWElmTU0AKgAAAAgAAYdpAAQAAAABAAAAGgAAAAAAA6ABAAMAAAABAAEAAKACAAQAAAABAAADNKADAAQAAAABAAABOAAAAAC5lQ7FAAAAHGlET1QAAAACAAAAAAAAAJwAAAAoAAAAnAAAAJwAADEPbX3syQAAMNtJREFUeAHsXXncTlUeL9llK1tCJLIlUk0rsqWmmkiDipG8FGX9pLGV7Fkiu3gJbyLLiApJGjIYWSKvNYXI+iJeNFTznXmmO6dz73O353nuc8+9X3/Uuef+zrm/3/ec33m+71l+54pfffAvb968V3j1LyUlxQcWUwUiQASIABEgAkSACBgjcIVxtre5JGfe4s2vEQEiQASIABEgAv5FgOTMv21DzYgAESACRIAIEIEQIkByFsJGp8lEgAgQASJABIiAfxEgOfNv21AzIkAEiAARIAJEIIQIkJyFsNFpMhEgAkSACBABIuBfBEjO/Ns21IwIEAEiQASIABEIIQIkZyFsdJpMBIgAESACRIAI+BcBkjP/tg01IwJEQCEEJkyY8Lq//6Wmpm7evPns2bMiqr/88svhw4dXrlw5atSoWNQfMGDAoEGDYqkh2GVnzZolws40ETBHgOTMHB++JQJEgAjYQqBFixZexdJ2/J2bb755yZIl5mZcvnx53rx59913n83aK1Wq9OKLL86cOXPr1q0a4fvXv/61c+dOZDZt2jRXrlw2qwqD2B//+Edz/PmWCIgIkJyJaDBNBIgAEXCJwNSpU/1JMtq0aZOZmWnfqtWrV999992GtmTJkqVOnTpjx449ePCgZYUZGRlDhw4tVqyYYVVhyyQ5s+wwFBARIDkT0WCaCBABIuASgX379vmQcAwZMsSdPQsWLLjppps0i2688cb+/fvb4WTS58ALBw8enD9/fq2qcCZIzqSOwUdzBEjOzPHhWyJABIiAXQRuuOEGXzGPvn37Gqp+8eLFLVu2rFmz5rvvvjMUiGSeP3/+5Zdfbty48WeffYataSaSlq9OnDjRpUuXHDly+AofL5UhObPsJBQQESA5E9FgmggQASLgHoFWrVp5+Xtv/q3atWvrGdXevXtbtmwp7gYrWbLkX//6VxdTYhJM2G0GMidlSo/ggs8888yVV15pqHnOnDmxZtqzZ8+5c+fi4AKOKZw5c+bYsWPgke+++25KSkrhwoUNCyqRSXImdQY+miNAcmaOD98SASJABOwiMH36dP8QBZy+lPSeMmWKSMtEVbNmzYot/Nu2bZOKmDz+9NNPn376KYhdrVq1ChUqFKktX7589evXx7FQE6K2bt26u+66S/t6wYIFn3vuuaVLl164cMHkc3h16dKlOXPm3HHHHVpZhRIkZ+aNy7cSAiRnEiB8JAJEgAi4RGD//v3+oQuvvPKKaAb25lvqhv3+kyZNEkvp02BICxcubNasWd68eU0qBF3DF6NRNEzpTZw4EXxl0aJFmHLTf8U8Z/HixVWrVjX5uuGr7NmzgxRidXXatGkI23HttdcaiiUok+TMvE35VkKA5EwChI9EgAgQAfcIYON8gn7dnVaLSaz169dHLAETiraYqFVbvHjx9957z9zyQ4cOYRlUK2KZKFGiBA6x/vzzz+bVuniLOtPS0izRzp07d7169RCDbdWqVdLMHFZLQUYtTYiXAMmZi1YOcxGSszC3Pm0nAkQgzghghS5eP+ex13PVVVfVqFFDPHQZrU5sREPkCztYYP0RjCdaPYb5lStXxgyZncqdymDWDXE9pDmwbNmyYaUVhOyLL74wn5a75557DBVORCbJmdPGDbk8yVnIOwDNJwJEIJ4IYDonET/tiasTE2wmE2aIKKtHZ9myZSBATlVCuDV9VVIODgEggO0LL7wA2oSZPKycXnfdddjE9uabb2LSThLWHk+dOtWtW7fy5cuj4AcffPDjjz9qr8wTKOLUCtfyJGfmbcG3EgIkZxIgfCQCRIAIuEfg+++/d/377X3BihUr7t6929BanJfE2UmohGga+lOf77zzjn1tsXrYqVMn80C4OItgvtMfdPDZZ59FMDlDbV1kjh492r4JsUuSnLloozAXITkLc+vTdiJABOKPQLly5WL/LfegBmzGMpxkQkwyLM6K+7EQ/wLnACSk+vTpY0fJa665BjNtUln9I/hf8+bNLSvEpv6XXnrp6NGj+hrs5Bw5cgTBdUE3q1WrZvmt+AqQnNlpIMpoCJCcaVAwQQSIABGIAwJt27aN7+96Impr0qQJYmHorZ09e7YWF0P87p/+9CdcvinKY0s+QqmJMvr0Lbfc8s0334iltLR+xg5RzUqVKqWvRMxBKBAQR8zqafVYJg4cODBjxozWrVvb2Xsnfiu+aZIzy5aigIgAyZmIBtN2EcDiAuKGe/YPR/ftakY5IpBsBLCFK76/63Gv7fHHH9fPhGEW7amnnjL5FlYVJWgxDuTJkydakXvvvddwZg6BcB966CGU0kdiw2mDaLXhAqhXX30Vs3qSDtEeEYMNNM7yOGe0z8U9n+QsWksx3xABkjNDWJhpgcCmTZviPniZVDhs2DALhfiaCPgGgR9++MGkMyf9FWJ9SUElgNz27dvt7I4fOHCgBHO/fv0MLcJ5ybNnz0rCmGwbOXKkeNhz1qxZkowYnzZSM4KAtGvXzj4ti1SIoww4TGCoW1IySc6khuajOQIkZ+b48K0xAlhW8HKAIzkzbgbm+hUBbLT30kHsfwtLljgRKcG2fPlynNm0UwkuEli7dq1YHNv8pUgWqAfxO86dOyeKIX3y5MkGDRpIX8F5zG+//VaUnDdvnigDxT766CNRQEpjsRX6Gx44xSED/1y4TnImNRwfzREgOTPHh2+NESA5M8aFuUTgvwggpoPIMPyTfv/996UmQuwJ7LK3r2HZsmVxb7pYSe/evcXiRYsW1Ye9AE8qXbq0KKalcdpArA1sT7sf/eqrr5a4oCi5ceNGHA4oUqRIpCosYup30X388ceW0Xc1TRKaIDkT245pSwRIziwhooABAiRnBqAwiwj8hgA4UEJ/6d1VXrNmzd8U/N//V65cqTEh+3WOGDFCrGfPnj1iWewbE98ivWHDBpzZFGXENCLlHj9+XCwSCeEBGZysFPMjacSVxf1Lt912m1hJJI2IaPqAHQiBppf0PofkTN+UzDFBgOTMBBy+iooAyVlUaPiCCPz667Fjx3wyYSOykA8//FBsHIRk069IivLR0rjBSYp8pu1Xw00D4ieQ/uqrryzXTBG5Vyz1/PPP49PSjFpEAPHVzO+PAgGVdtThfk/zItHMjG8+yZnYxExbIkByZgkRBQwQIDkzAIVZREBAoEqVKvH9dY+xtgIFCkixMJ544gnXda5Zs0aw9denn34aVWFH2sGDB8V8xBWzjI6Bgn379hVL9e/fH9NpUrxZUC7wGzsKP/nkkxJ3xG3udgomVIbkTGxipi0RIDmzhIgCBgiQnBmAwiwiICDw4osvJvTH3mnlWCsUtPsV4cdimdubMGGCWFuPHj2gD2LriJlIP/jgg3b0BBsTCw4ZMkTSFm8RmM1OVREZXPckVoi1Tpw8sF88EZIkZ2KLMG2JAMmZJUQUMECA5MwAFGYRAQGB+fPnJ+I33nWdIDeCdr/ihKPrqlBQCoHRs2dPZOJCJPETNhEoVqwYFljFguC10lzakiVLHGmL6GtYWRbrfPTRRx3VEHdhkjOxOZi2RIDkzBIiChggQHJmAAqziICAAOJyxTI1FXdygJCwgnb/Sb7yyitYPXTxIdxPJe27jyxrSvEsEALNHIGcOXNij5q0EgrFsCI8efJkUVtczelUzzFjxog1ROb2nFYSR3mSM7E5mLZEgOTMEiIKGCBAcmYACrOIwO8RuPXWW+P46x5jVbj4SB/HddeuXV27drWzLSzydQTdwD1I0uFKhJYtXrw4BIYPH/57AH7dv3//22+/jZkwrHgiyBlOUyKRkpIyePBgTIYZ3h+AA6SoauzYsWJVkHdq/qBBg8QaXn/9dac1xFee5ExsDqYtESA5s4SIAgYIkJwZgMIsIvB7BDp37hzfH/gYa+vevfvvFfz/E266nDlzJubScPv4/fffjyC6uPioRIkSZcqUQUTZRo0aIZgZblEzpFPadVUI7v//Gl2lUH8kfi8CxYkVnDp1qkKFCvbNr1u37unTp8Ua2rRpY794IiRJzsTmYNoSAZIzS4goYIAAyZkBKMwiAr9HAGwmET/zruvEIiPi7/9ex1ifcN+AeEvS1KlTXdeIi8+1m9TBC6XbP0G2MAMnXv2kxwGRQUAuP/nkE0kHHFNNejQNkjOpUfhojgDJmTk+fGuMAMmZMS7MJQICApjvyZIli55DJDEH+gwdOlSKqSGo7Cx54MABKWJItmzZpO1iNmtE6FrM0onIGNaDqbXFixdjyRKzkq1atcIKaceOHXG6c86cObhMM9q3cLZUrDkpaZKzaK3DfEMESM4MYWGmBQIkZxYA8TUR+C8ChoHsk0IOxI9Wr1590aJFUjAwRy2GfWbTp0+PFsMWE2CYNZSmvgzrx/Y1bErTX3YObRGYDYuthqUcZa5atcp8vk1EJnFpkjNHrUZhkjP2ATcIkJy5QY1lwoeAT+4OMuQcuCUTITBwtxKYlv2W+eGHH3AQ0s7N7mBXDRs27NWrF2LA4gZPnADAgiNuZMLSJw5ytmjRonLlyubHOTGXFiM/S01NxUkIQ/M9ziQ5s9/HKAkESM7YDdwgQHLmBjWWCR8CuDHJYxLg4nO4XglBX7FQOGrUKGxKW79+fXp6+t69e7Fq+e233+L+pWXLlo0fP759+/bVqlUzp1Muvm5eJH/+/LhJ08Uk3+rVq7UdbOaf8OYtyVn4vD8mi0nOYoIvtIVJzkLb9DTcEQLY5O4ulpg3jEGVr1StWhVb0E6ePGkJ/o4dOxDRA+u2fjON5Myy7SggIkByJqLBtF0ESM7sIkW50CNwxx13+I0oxF2fokWLYo0S/5CIe+VahaC5t99+e9u2bXECADegY5IP/xABZOTIkTjIWatWrWh74LQakpggOQv9SOAMAJIzZ3hROoIAyRl7AhGwiQCiiyWREyTu0ziY2axZM9AjXJT0L+Hf0aNHkdm0aVMIxPh1xDbDWirOWka2rM2ePXvgwIG4iMkPG/ydmkZyZtNfKBZBgOSMPcENAiRnblBjmVAi4PReSKe/+kmRx3Z+bEcTKJlBct++fZFrnVxo+MQTT6xbt86g0v9mZWRk4FACYqG5qDlZRUjOQun97o0mOXOPXZhLkpyFufVpuyMEzp49mzVr1mRxgrh/FxdipqWlRaNN+vwZM2bkyJHDvhq4Bx2Xsuvr0ecgLC2WOO3XnFxJkjNHXkNhkjP2ATcIkJy5QY1lwoqAYRyv5HIFd1/Pmzfv2rVr9TzJPOcf//gHCtr5IsLCWU7ISd8aN26cEkcuSM7C6v0u7SY5cwlcyIuRnIW8A9B8Rwj06NHDDjXxuQyCaCB0rcSNbD4iJq1lDI5y5cpJ29dsVo4wHz6HDuqRnDlyGQqTnLEPuEGA5MwNaiwTVgQQfNX/7MFSQ9x9bpMtGYoh4K3JJzC1tm3bNsOCdjJxbsCkcj+8IjkLq/e7tJvkzCVwIS9GchbyDkDzHSGQmZkZ+9HF5DKMUqVKYfOcHZ4UTQbFTbbwv/XWW9EK2snHnZulS5dOLkTmXyc5c+QyFCY5Yx9wgwDJmRvUWCbECNx7773mP94+f/vmm2/aIUnmMggPa2gmmB/4q3lZy7eIUmtYuU8ySc5C7P1uTCc5c4May5CcsQ8QAUcIYE3QJyzBhRpZsmTBlZrm9AjhzV5//fUaNWpggRI79DFJhjgan376qVjq0KFDqEqvAAqKYpbprVu34o7OWbNm4VLzc+fOReQxM4e7nvSV+ySH5MyRv1CY5Ix9wA0CJGduUGOZECOwYsUKn7AEF2rgkgNzwoTwsLjm3LDmxo0b49olrTjOY+rFcH2nJmCSuHjxIu5Qx33tYg3ggh06dDh8+DAKNm/eXHzlqzTJWYi9343pJGduUGMZkjP2ASLgCIHz5887CvflK2LRpk0bE8709ttvm2sLQoaYZJEannvuOUkYrM6kcu3ViRMnHnzwQams9ojoaAjYgci0Wo7fEiRnjvyFwiRn7ANuECA5c4May4QbAVz+6DfGYFOfvn37aiRJSmzatMnOWYdnn302UrBPnz7SR3FdplSn/vHChQv169eXCkqPIHmpqalSpn8eSc7C7f2OrSc5cwwZCwABkjN2AyLgFIHXXnvNP1zBkSZDhw7VE6ZIDm66tFMVgpxF1i5xZ7kkX69evWiVa/mjR4+WShk+VqtWzTDfD5kkZ079JeTyJGch7wAuzSc5cwkci4UYgc8//9wPLMGFDsOGDdN4kphAzFj70fkR5wxl33jjDUmBhg0binXq0z/99FPx4sWlUso9kpyF2PXdmE5y5gY1liE5Yx8gAk4RwH52XEypHKuAwtHI2dKlS+2bE5khc0HOsJnM/ld8K0ly5tRfQi5PchbyDuDSfJIzl8CxWLgRqFOnjm/Zg4li0cjZ9OnTTUpJr2655RZ3M2fTpk2TqlLxkeQs3K7v2HqSM8eQscAvv/yyY8eOVq1aYQ8NNqNgszBOYFWvXt3y7jxHQyrmGIoWLXrNNdegWvw2EHYiEAAE+vXr58gLfCIcjZxNmTLFvoYVK1Z0R84Q/9b+V3wrSXIWAP/10gSSMy/RVv5boGWXL1/WbwqJ5Ozfv3/AgAFFihRxPT5mz579ySefTEtL27dvn/YVxJb8/vvv8V18XXkEaUC4EVi9erVr70hiQZKz2MEnOQu36zu2nuTMMWShLfDzzz9rhMkkkZGRgZiQTscyTI8hltJ3331nUjNeXbp0iRQttD0wAIZjb3vu3LmdekfS5UnOYm8CkrMA+K+XJpCceYm2qt8ynzAzpFNTp07FNJjNEQ2Xrnz00UeG9RhmYhZNVSipd+gRsIzXZdNrvBQjOYsdbZKz0Lu+MwBIzpzhFU5pzFcZkiTzzEWLFhneoycNc7iD7+uvvzavSv+W/CycXTEAVg8cOFByAf8/kpzF3kYkZwFwXi9NIDnzEm0lv2WyyUzPmaQcHBcwH9Ry5cq1fv16qZTNR/IzJftT6JVWMTAEyZn5OGbnLclZ6F3fGQAkZ87wCpu0zX1mJnTqgQceMBm5cEjepKzlK6gXthahvaojgHnoq6++2sQpfPgqGjmbOXOmfW1vvfVWeLSLOGc8ral6n6f+LhAgOXMBWoiKuFvQFEkVJsaiDd/333+/KOkuzfMBIeqOQTEVMfGjOYU/83HnkqF7rlixwr7CkZsABg0aJBWxvCFg8uTJUhEVHzlzFhT39cgOkjOPgFbxM+A9hiOy08z77rvPcDD9+OOPnVall+fkmYpdK+Q662ePDB3EP5m9evXSux5yTp8+bf/OA2y2Q5EePXpIdtWuXduwci3TEQWUKvfPI8lZyL3eqfkkZ04RC5G8/d1mJ06cSE1NbdmyZd26dXERcteuXXGNoDa2vvXWW/ohEucAEFZAk9myZQv+PsYf6JMmTVq3bp2Wb5nA3F6ImoSmBgIBk+lkvaf4Iefpp5+O5oktWrSwoyHObu/duxeVNG/eXJKvVKlStMoj+ZmZmTjQLZVS7pHkLBC+650RJGfeYa3cl2yuaY4fPx5x/PVjJSbMIscwN23apH/bunXryMj74YcfVq1aVRIoW7YstrOYD9naW65sKte1Qq4w/uzJly+f1Of9/Fi+fHnN3aQEKJcdW1555ZVIwZtuukmyFLwNgaalaqXH3r17S6UMHwsXLmyY74dMkrOQe71T80nOnCIWInlpfDR8bN++vcnAh793v/jii3PnzulvdhozZgwq1K9xiLU9++yzFy5cMPyumElyFqJOGRRT8VMtdnX/p9PT00WnE9OLFy82D2r48MMPRxx527ZthpYuW7ZMrFCfPnXqVJUqVQzLapnZsmXDJXLao98SJGdB8V2P7CA58whoFT+jHyKlnNGjR1uOgLgfE5cv6f+2XrhwIY6AWRbv1KmT9FH9I7edqdi7Qq7z8OHDLTu/rwS6deumdz0tB9dSlStXTq9w1qxZUVD7E6tLly56GeTgzzCtqmgJTNHhdk7D4sjE1jdM4dsJrBithkTnk5yF3OWdmk9y5hSxsMhbngY4cuSIzY0gbdu2LVSokDT24QJN87+2Nfm1a9dGG68j+SRnYemUAbLzyy+/1Hq4EgmE/zh48KCJJ54/f37GjBm4GxchM7B2ibPYWIvcuXOnVuTAgQN58uQxNDZHjhx4q0lGS5w8ebJz584QliqpU6fOxo0bozE/SThZjyRnAXJfL0whOfMCZRW/YUnOxo4da3OYw2WC1157rSTcpEkTKSfaIySjDdaRfJIzFTtYyHVGpy1QoEC0Pu/P/MaNG5t7ovnbxx9/3MQuHCwwL669BUWbM2cOQnL06dNn4sSJu3btwqsdO3boSZvJ57x/RXIWcpd3aj7JmVPEwiJvSc4wUtsf4PRzbDitabM4/mS/ePGiNjTrEyRnYemUwbLzscces+kC/hHr16+f3gHt5PTt29fSiilTptipSi+DiB6YrrOsP7kCJGfBct+EW0NylnCIFf2AJTm788477Q92+pDoV111lf3i5kseJGeK9rGQqz1y5Ej7LuAfychRHj1DMskxDKajtwj7HD777DOTegxfYUNbo0aN9LX5LYfkLOT+7tR8kjOniIVF3pKc6eNfmIyGenJmIqx/tXXrVsNxOZJJchaWThksOxHbT9/Vlchp164djmCbuKT2CmIpKSn2jcLS5PTp07XilomjR49iw5n9+pMoSXIWLPdNuDUkZwmHWNEPkJwp2nBUWxUE4GKGAQKTSCDsfxqRz2bPni3GkZaIFLYivPfee4ZHOC2/giNEx44dkyrUPy5ZsuTGG2+0rM0nAiRnqjimT/QkOfNJQ/hODZIz3zUJFQocAkqsx5mQm5IlSyLYzfz587dv347j2/iHxLx58zp27Gh/U6lh/Qi+g4tGcKZVz8kyMjIQoRqnQQ0L+jaT5Cxw7ptYg0jOEouvurWTnKnbdtRcFQTsRAr0LdvwRrEiRYpg4fKZZ55BLLQnnniiWrVqjrareqOkna+QnKnilT7Rk+TMJw3hOzVIznzXJFQocAhEi5hv58eeMmohQHIWOPdNrEEkZ4nFV93aSc7UbTtqrgoC8DI/XwepFvvxubYkZ6p4pU/0JDnzSUP4Tg1LcuboGjue1vRdA1MhfyBgPxqzz8mHZ+rhuhGscv7hD3/w7Itx+RDJmT8cThktSM6UaSrvFdVvxRVzGjRoYH/MypUrlySMS/ekHJNH3M4pflpKM5SG932DX4wXAuPGjTPp+XwlIYC/Cc+cORMBf8SIEdJbPz+SnMXLZUJSD8lZSBrajZkSB5Ieu3fvHstQaP8wF24XMDmxD61Izty0Lsv4A4H09PRY/ChsZXEyVGs3OH7x4sVVQYDkTGs4JuwgQHJmB6WQyly6dEkiZOLjP//5T5vDIg5Y6SXr1q2rzzTMeeqpp8Tv6tNYgQ1pC9HsQCBQtGhRw57PTD0CJ06cENu8QoUKehl/5pCciQ3HtCUCJGeWEIVX4PLly3omJObgZLudcRBBj/RiKKtf69SLXXnllRs3bhQ/qk+TnIW3jwbC8qZNm+p7PnP0CCDkrNjgmZmZCoXVIDkT245pSwRIziwhCq+A5ZkA3J1StmxZ/Rgq5vTq1atmzZpiTiSN2YLx48fr86Wc3r1769mYmIPpvfC2EC0PBAITJ06Uuj0fDRFo2bKl2OCrV682FPNnJsmZ2HZMWyJAcmYJUagFzFc2QZJwJfldd91lOBrij9pBgwYhxrfhW2QuWLDgjTfewNxYNIEuXbqY7zaDAtxwFuoOGgjjd+3aFc0FmC8ikJqaKjb44MGDxbc+T5OciW3HtCUCJGeWEIVaANRHnKYyTOMSvWnTpuE2lSxZskTGR5x1Rzjvr7/+GvL169ePNmjefvvt4F6rVq3S38Ry2223LV682PBzYia4I9c0Q91Bg2K8Qhvbo7mzB/kHDx4UG7xhw4YefDRenyA5E9uOaUsESM4sIQq7gOXOM40tnT9//ttvv8Vap5YzcuRI86EN19dEhFFw7ty5EyZMmDNnzs6dO7UazBOcNgt77wyK/Tj1Yu4pfHvrrbeKrY2/CfPkyaMQLCRnYvMxbYkAyZklRGEXwNSU5eKmIYVasmSJ5XbdbNmyrVy50rC4ZSaZWdi7ZoDsnzx5skI8Iymq9uvXT2zwTz75JClquP4oyZnYfExbIkByZgkRBX51wc/S0tJy5sxpZyDDGui6dessqZgkgPk8NgwRCAwCe/futeMsYZbZs2eP2Nzt27dXCw2SM7H5mLZEIHTkDIHpc/OfKwQwGp48eVIiSfrH06dPYyO/o3ETNA5kTl9VtBzOmVk6NgWUQ6BkyZKOvCZUwjjxLTYoRgDlgsORnIktyLQlAqEjZ6Ea0eJuLGa5sLiwf/9+Q9p0+PDhoUOHut7a3Lhx423bthnWrGXyBIClS1NAUQQQJyLuDhuYCrEhVWxWbJlQzjSSM7EFmbZEgORMOR9PvsIIfoHTlCkpKQMHDhw1ahTiZTz//PN33nmn5Q4zS9VRc506dcaMGbNhwwbMwIGT4TgnLtbcvHkz1jF5MNPSnymgLgJTp061dJBwCpQrV07axmAz/LWv4CI5U9c3k6I5yZmv/JfK/A4BrEFHoqAhhHpS3IMfJQKeIYADy7/r/Xz4DYF3331XbAXM3GNk+O2lMv8nORMbkWlLBEjOlPHtMCsaYHJ25syZESNG4BZ53FJl6a4UCDYCpUuXDrObG9p+zz33SFPmnTt3NpT0eSbJWbCdN+7WkZz53KOp3n8QCCo5wxa66tWrR9o4e/bsuI4m7h7OChVCAKGb6fAiAjgqlJ6eLrbgoUOHbB4DF+vxQ5rkTGxHpi0RIDnzg9tSBwsEgkrOpMsBb775Zuyxs3RaCgQVgenTp1t4QsheI/yb1NatW7dWFAOSM6kp+WiOAMmZop4eLrWDSs4++OADqSFff/11c4/l2wAjgJtqpf4Q5se+fftKbb127VqTq3h9jhXJmdSafDRHgOTM5x5N9f6DQHjIGZZsEI/U3Gn5NsAIlC1blj4PBta7d2+plXE7HKaW1QWH5ExqUD6aI0Bypq6ze615iRIlMDgm5S/X8JAzNGqDBg3MnZZvA4xAmzZtvHZsn32vYMGCf/vb3/RNrPqGPJIzfZsyxwQBkjOfjUx+VQdLDJG4/GvWrMmbN2+81MRlxo8++mjhwoXNKwwVOQMUs2fPNnFavgowAggbYe4Lar1FHP+lS5deuHBh06ZN4J358uUz0R8Xmb/88ssZGRn69h05cqRJQSVekZzpm5U5JgiQnCnh10lWsmLFimIfGjZsmLlCkO/QoUOTJk1wAjGaJGbgUlNTI9VmZmY2atQomiTycR+UqEBg0vo9ZxEQrrvuOoTYCIyZNMQ+AjiNaOIIyr1CDxdtx3mXFStWvPbaa/D3GjVqlC9fvlKlSrVr137hhRfwB8nZs2dFYS2N8LxJmbCPL9okZ1qDMmEHAZKz+DpgMGt75JFHxM70448/FihQIJqpbdu21cJ5b926NdoVeB07dhTrBD+7/vrro9WJa51E4cCko5Ez4PDSSy8Fxkwa4ggBUJZojqBcPq73cGS7Xnj48OEBYGZoOJIzfeMyxwQBkjPlhrskKIxdICBPYjfCKqehHjfccMPFixdFyc8++0wviW3vJ06cEMWQHjJkiF4SOVj6lCQD82hCznAXFsPSBqahHRnSrl07Q0dQMfP99993ZLsojAvcmjdvrqLVhjqTnImNy7QlAiRnhn7ETBkB3HcpdiaMm2BsstAVV+Dic1EsksbEmySJMVcvdvDgQcM/kd977z29cDByTMgZELvjjjsi+/yCYSytsIkAFvgkf1H3EdtJly1bZtNwTQxXAqSlpZlMpasICMmZ1r5M2EGA5ExFN0+CzqVKlZLio+LKc70ee/bs0Xc7/eTZggUL9GLIAR2R6qxQoUKACYo5OQMUY8eONQSKmQFG4MiRI5IXqP542223jR49Gn99WbbayZMnx48fj02rqpus15/kzLL1KSAiQHKmdyLmGCMwbtw4sevg6iFpDC1TpowoIKbLlSunVYpTAtF2/vbr108TiyQWLlwo1hOwtCU5y58//w8//BAwq2mOJQKSZ0lOoe5j5cqVU1JSJkyYgOk0bEjF33I7d+5cv379nDlzXn311Vq1apkcIVLX6ojmJGeW3Z4CIgIkZ6q7vHf64wihRKq++OKLLFmyaBrgKIDYt8T0wIEDNbF69eqJr8T05s2bNTEk6tatK74NXtqSnAGEZs2aBc9wWmSOQPv27UVHYDoACJCcmfd5vpUQIDkLgNd7Z0LPnj2lDoRA3trnoy1WoggCBGTNmjUiKW1fkyrEkYKIGIIeGS6SSvJKP9ohZ0Bj+fLlSptJ5Z0iMHfuXM2tmAgGAiRnTr0g5PIkZ8FwfI+swCnL3bt3iz6DDWGIIovPX3311dKJTlEMaQSShRgomrROh+VRUbJbt24RY6ZNmybmBzJtk5xhUVg6AxtINGiUhsCxY8cMD8d45Of8TAIQIDnTujcTdhAgOUuAFwa6ypo1a0o79BH+G6cvO3fuLHa448ePz58/X8xJT0/Pli2bdE4ThwwQx0gU+/rrr4EfNqCImUFN2yRnAARxO4MKAu0yRKBKlSqBHkhCZxzJmWE/Z2Y0BEjOQjdGxG7w4MGD9f0Jp9/FTBy5Ao0Tc5AGXTt8+LCYOW/ePITcFHOQXrVqlZQjHRSV3qr7aJ+c5ciRQ5qzVNdqam4HAUQhjt1VWYN/ECA5s9PtKaMhQHLmH+dVRhPER/3kk0+0PmSYQORY2LN69WrDt1omCBzEPv/8cy3HMDFp0iTDfNUz7ZMzoFS/fn3V7aX+9hHADk5lRgQqagMBkjP7nZ+SQIDkzIZXUUSHAO4+X7t2bTQXQljwSAnzMKqLFy+OiD300EPRqkI+IiSF7eJzHd7/y5g1a5YJUHwVJAQQ8YvbzqI5gor5JGdBck8PbCE5U9HNfaEzTlMa3s3yzTffFClSRFOxT58+hv14//79hQoV0sRA1AzFJk6ciGgdJGcRoIoVK4a7GQyBYmbwENiyZctf/vIXnMLR3IQJdREgOQuehybUIpIzdZ3dF5rjqObKlSsjJy5xWvOdd97R33SOABzSpjHMupUsWVI0ADepf/nll2Jfx8xBmzZtIjIkZxpWHTp0EFFiOvAI4BbaAQMGSP6i9QcmVEGA5CzwrhpfA0nOVHFtX+uJABnXXnutySpMiRIlEFcTBzP79+9fu3ZtQ2NQyXPPPYepMlwM1aJFC8Tm0MRIzjQoMI+4YcOG+I4CrM3/CFy+fBlh9O+9916tJzChFgIkZ/73Ml9pSHKmloOHVNvu3bv7ym3ipYyjAwFa29eoUUOKZhIvfViP/xH46quvWrVqlStXLq0/MKEEAiRn/ncuX2lIcqaEX4ddyalTp/rKbeKljDtyht6AQxLx0oH1qIgAFv1xJVqpUqXCPjSoYz/JmYqOlkSdSc7Uce4Qa0pyJjV+vnz5pIhxSRxE+OlkIYC1ThzKuf/++6XuwUcfIkByliw3UfS7JGc+9GKqJCNAciYjcsUVf/7znxUddKh23BHYunVr69atudapdxP/5JCcxb3bB7tCkjP/OC81iYoAyZkhNMuWLQv28ETrHCGQkZGB2ztuuOEGw97CzOQiQHLmqDNTmOQsuQ7Lr9tCgOTMEKabbroJF5tyFCMCIgJY68StaJG7Nwy7DTOTggDJmdhLmbZEgOQsKX7KjzpDgOQsGl7RYvxaej4FAo/Atm3bEJsmd+7c0ToP871EgOQs8B4XXwNJzrx0T37LJQIkZ9GAw4Xou3btiu+gENraVqxYMWbMGHCaICGAtc4hQ4aULl06WhdivjcIkJwFya08sMUX5Gz37t07bf9buHDh32L4t3z58k38pxoCCJLugTN4/wnXoTTEn5O6det6r3nwvti1a9cIqoil3LJly4AdhkVgvPnz50eL/yx2J6YThADJWfAGjYRa5AtyllALWTkR8C0CcSFn+C1JS0vzrY1KKDZz5kzpJxkXVGBz/cWLF5XQ376S27dvT0lJwcW4kr18TDQCJGf2eyklgQDJGbsBEUgaAvEiZ7jP9NSpU0kzQ/EP79ixQ7wrTPyRLlu27KJFixS3z0B99JahQ4eWKVNGNJbphCLglJzt379/I/+FGAGSM4ORi1lEwBsE4kXO8KPywgsveKNzwL6SmZlZqVIl81/levXqpaenB8xwmIO1TuwQeeCBB8zN59u4IOCUnLVr1y4u32UlqiIQvBGHFhEBVRCIIznDhejr169XxXD/6Il7Ku2M3VmzZu3UqVNQpydBPdu2bcu1Tjs9wbUMyZlr6EJa0D+jJDUhAmFDII7kDONX9erVEeMqbBjGYu+UKVMcjfuFChWaNGlSUG+dP3369LBhw7CS6wgTCttEgOTMJlAU+x8CsQxtLEsEiEAsCMSXnMGlR40aFYs+oSqLkBnuYoBVq1bt73//e1CxAvXEiXgcAcapVf5MxhEBkrM4ghmKqoI6xNAuIuB/BOJOzvLmzXvo0CH/G550Dc+ePXvzzTfHMsQ3bdr0wIEDSTckcQrgnMTzzz8f7ahELNCFsyzJWTjb3b3VifNt1kwEiIA5AnEnZxgImjRpYv5RvgUCzZo1cz9o/lYSF4337dv3/PnzAYYUa50jRozARWG/Gc3/u0SA5MwlcKEtFuBhhaYRAZ8jkAhyhqFsyZIlPjc8ueqNGzcujgN+qVKl5s6dm1yLEv31X375BSFFcGqVa52uew7JmWvoQlow0V7N+okAEYiGQILI2Y033sgL0aNhvnnz5pw5c8Z9uK9Vq9aWLVuifTQw+bjJBUFbuNbpov+QnLkALdRFAjNq0BAioBwCCSJnGNF69eqlHBoeKHzmzJnErdBdddVVIC7Hjx/3wJDkfgIwvvnmm+XKlQv1b6dD40nOHAIWevHkOjm/TgTCjEDiyFn27NkxyRFmbA1tx4a8RA/5BQsWHD169KVLlwwVCFIm1joXL17coEEDrnXa6VQkZ3ZQosz/EQjSYEFbiIBaCCSOnMHDEfldLTQSre3IkSP/P/AlOFW5cuXly5cn2iKf1L9r164OHTrgpHCCQVW7epIztdvPe+194t5UgwiEEIGEkjMMJjNmzAghqoYm4/oEzCZ6PMA2atTom2++MdQneJk//vgjwuyVL1/eY5BV+RzJmSot5Rc9gzdG0CIioAoCiSZnRYoUCeqNQ46aOCMjo3Tp0kkZc3PkyNGzZ89z5845UlhdYax1fvTRRw8++CDXOqX+RnImAcJHCwTUHQWoORFQHYFEkzM4P65PVh2l2PV/7LHHLMbBBL8uXrx4WloaiEvstqhSw+7du1988cV8+fIlGFplqic5U6apfKKoKq5OPYlA8BDwgJxhAmPt2rXBg86+RW+88YZPBtu77rprw4YN9jUPgCTWOt96660YL2PwSfPFqAbJWYwAhq54APyfJhABRRHwgJxhRMNdkKG9EH3NmjVZs2b1z7AOrty6desjR44o2mPdqY0pw48//vihhx4K81onyZl/3FANTdw5G0sRASIQOwLekDOMRIhKFbu2ytVw4sSJEiVK+HAgxmLfsGHDfvrpJ+UgjVHhPXv2dOzYMX/+/D5slESrRHKWaISDVn+MzsbiRIAIuEbAM3KGkO7ff/+9az1VLIjZmoYNG/p5vEYQV+ydVxHbGHXGrfNjxoypUKGCn1sn7rqRnMUd0oBXGKObsTgRIAKuEfCMnGEUa9y4sWs9VSzYv39/JcZuLPaFM1ww2PPSpUsffvjhLFmyKNFSMSpJchYjgKErruKwS52JQDAQ8JKcYWgLzzzNypUrcZmSKqN5tmzZunXrhjuRgtGrnVqxd+/eTp06BX6tk+RMFX/0i55OHYnyRIAIxAsBj8lZmTJlzp8/Hy/lfVvP0aNHr7vuOr+MsLb1QFC61NTUn3/+2bfAJlQxhIIbO3ZsxYoVbQOmmCDJmWINlnR1E+pvrJwIEAETBDwmZxhtevToYaJPAF6B3NStWzfp46prBWrUqIETpgFoCNcmLFu27JFHHgneWifJmWunCGlB1y7EgkSACMSIgPfkDCto6enpMart5+K9e/f2+VCOgLQpKSmTJ0/GjqtVq1bh7nBc+tmsWTNxXe+pp54K2wEOqVPh2qsuXboUKFDA561pXz2SM/tYUfI/CEguEbzHDz/88Br+IwLxRiAukV29J2dw+Vq1agXPzSMW4a7xuMy45M6dOxHR0XBRwYoVK6ItXCKyxuzZs2+//fbIL1OePHkGDBhw4cKFoDaWHbuw1jl+/PhKlSpFMFH6vyRnSjdfEpS34yFKyyxYsCAJsPKTQUcAcx6x+0VSyBla5p133oldeb/VcOjQIWzbctfvQOnq1asHHrBly5aLFy9GTDt+/Dgmtzp37lysWDF31Wql7r777o0bN9pEDBTt+uuvj5TFNkGMYDYLBlgMtPvRRx+NC/PWGsXjBMmZx4Ar/7kA+3PENJIz5fuoLw1QmpwVLlwYd4EHyfdxBULNmjVd9BQc6sT1ozgwaIIG5rrmzZtXtWpVF/VjBm748OHRZsuifRTX1T/99NPa5+rUqbNt27ZowuHJ37dvX9euXQsWLKgho1CC5EyhxvKFqoF3bJIzX/SzwCmhNDlDa2DbU5B8v3v37i662D333LN9+3abOCAu15QpU8SdYZZfBI3AOqbN+vViEyZM0FZXQSJxj/jJkyf1YmHLyczMBDJVqlSxxN9XAiRnvmoO/yvzbwAAAP//2zNFqgAAN/lJREFU7V154FfD2nezpLTRon1DuxY3Uq5ERNkqaZWURCHRJloQRUq2QiVbSlmivUhcKYpyc4VK0SaRkpZL8r4f7/d17pizzdnmnO/5fvqDOXOeeWbOZ555vp/fLM8c8T9p//faa68dwX9EIGwE/vnPfwYfOm+88UbY7VLV97e//W3ZsmXBPyEJGubOnYvPUf3y/5M78sgj77vvvsOHD3tt/9atWxs0aKBSV5EiRVatWmWp/4svvrjnnnsuvfTSc889t2PHjuPGjduxY4el5IIFC/Lnz29UV7Ro0fHjx//222+WwrmW+dZbb11++eXoSgOfJCcuvvhiTx10/fXXJ/lz2LbIEfBkLtkoTHIWuQ3lZAXZTs7QabVr1z506FA2DmqxzZs3bwZl8WSDhQsXXrRokahESv/yyy8OBOjgwYNt27Z1rvHYY49dvny5pBaPK1asaNasmbnsMcccc91113377bfmIosXL86bN69Y5NRTT3377bfNkrmZ8/XXX/fr1+/4448XIUpgmuQsgZ2S6CalfjyTnCXa/rK2cSkgZ8B+9OjRWe0Bfv3114YNG3oyouLFi3/yySeWXz1//vz27duXLFkSCvPkyQMOdNddd33zzTdmYVC3K664wqHeiRMnSqXAg/v27Qu1DqUw2TZlyhSpIB7vuOMOc6k2bdps2rTJLJybOQcOHHjqqafQZWagEpJDcpaQjsiaZqR+JJOcZY0tZlVD00HOChQogJmn7HUCffr08WQ1+fLlW716tfl7MctVv359S1WY0+rVq9eePXukUuCF55xzjmURzKtJwnv37j3vvPMshc2ZN954ozRvt3DhQrMYcjA/N2TIkP3790vV5fIj5hRbtmyZwLVOkjNLG2amLQKpH8YkZ7Z9zxcBEEgHOQMArVq1ylIn4GNo33zzzdLHYtsZ5qWcJ7SAUpkyZT744AOp7Pbt20888UTJiMB3pdVJLIPa0TiprPHYunVrkD+jupEjRxqvzImyZctOmzbNEGYCCGC+s3///ieccIIZrrhySM7iQj5b6039SPbhwbO1L9lujQikhpwBs9mzZ2edH9i4cSMWAb12+IMPPih+6b59+5o3b66oBNNUc+bMEYsjbT7S0bVrV0mmW7duilWIYmgYJjUxhfbqq6+KZwJEGTH9j3/8w+78gdSe3HnEWueECROwt1IEKq40yVlcyGdrvakfqCRn2WqayW53mshZxYoVs2tpDBv27VYhna2mUaNGxgnNn3/+2et+NfCzd999V/KZF1xwgVhpjx49RIFXXnlFfOs17Wl5DvN/OFWwc+dOsQFMA4ElS5ZgMtITmF57ylU+yeQMB5B/4r+EIXBE6ocuyZmr16CADwTSRM7w+QMHDswiV4BNYD66LFMEM1KY7po6dWqtWrV8KMFxAWnVEvvQRT1YSjMEvv/+e6yHim81pHEWdezYseKqaBb1bKRNxUwk7Nzr2d6wuizJ5Ax/tESKPJX7QIDkLKyhRz25hUDKyNnRRx/92Wef+fAg+ovMmDEjXlPDMUnxq4cNGya1p3Tp0sOHD8dJz3LlykmvtD1Wq1YNMdLEdjKdQQBbACdNmlSnTh1tfZGpiOSMFugJAZIzzSOU1aUEgZSRM/TK2Wef/fvvv3tyH/qF161bV7Bgwehs6Kyzzho6dCjWIt9//33EJEOY0yeffBLbyEqUKCFWitmyzLevWbPG065znP2sV68eQqci9iwOFZ522mlSDDOxluBpxLldv369/m7KihqxQo14KEcddVRwnFU0kJxlhVUkp5EkZyrDijJEQEYgfeQMXzh58uTk+CZzSzDnEdGEB3YjYa8YmB9WAy3/YWs5VkKrV69u2AFWRRHfH2TLyHFOXHbZZTNnzsS2Fkk/tr7hQAZYgtdLDpyrM96ihQMGDEAsDzOezAECW7Zsuf3224sVK2YgFlGC5Iz25gkBkrOIRiLVphyBVJIz/ETt2rXLkwfRKXzttddGYVUnnXTSRx99JHEmy0dQtEGDBnltwxlnnLFy5UpLhWImQuPixKVX5Yry2Cr37LPPJn9mVKc5iXWB9z/99NOY1FTE04cYyZkIONOuCJCc+RhlLEIEjkglOUO/ggC5eo1YBJ5//vkozK5Jkya411IkSa5pTKHh2KZiY7p3746TsK46MwKgCDfddJOiZh9ioInmaG2x9GZiK8W4vvLKK6NY6yQ5S2ynJ7NhJGc+XByLEIHUkjMsri1dujRp3mrt2rXHHXdc6GZXt25dhP5XZE6i2IsvvqjSmMGDB4ulFNPOIWdV6nWQQf9effXVCJ+btC5OVHuw1onQxLjpywFJr69IzhLVxclvDMmZ1yFGeSLwBwJpnTnDp+GCwkRdiI6Zpxo1aoRudljD/eqrrxQJk1kMG5Wcm4RLnMylFHOuueYaZ+UB3+Iag/vvvx/h4pL/ExVjC//zn/8888wzOLQREO1McZKzGLsyG6smOQtl3FFJziGQYnKGvhw1alRy3FmXLl2iMC/EU1CkSpZioIwIV2HXMMT19Tcnl6kLVxdUrVrVTnlY+dhsh5BvyenoxLYEc8mg2gHXOknOEtu/yWwYyVlYjo56cguBdJMzrCHidsIk+KyJEydGYVinnHIK5kUsWZd65pQpU+zahlfqeiwlcbTTTnm4+c2aNcOqcRL6OuFt2Lp165133ikFVVHvC5KzhPdv0ppHcqY+uChJBP6LQLrJGb4Tsbhi91aIIpYvX77/gh5e6oEHHshQIoTPePjhh3HlAIKZ4af3zTffxK58S7ZkzsTkmWUIhvLlywdnfqjOYWYuPCT+0IQ5oVtuuWX37t2x93jyG4CexbnXv//97167gOQs+Z2bqBaSnHkdYpQnAn8gkHpyhm+Md80LAcCqVKkSkbV9/vnnmAjp0KGDWT8m1WbNmmWmYpY5nTt3Nmvo168fhLGsOWbMGIT2BYFDpNkKFSp06tRp4cKFlnosM3HHgFl5dDnY/447KBP1+5TkxiBMcfv27XG1hmKPkJwluTcT2DaSM8WRRTEi8BcEcoGcYQYIk0NxuS388v0F8fAewJY+/fRT57uVsF/ekjBJmePHjze3a968eSBhCC1mfoUc/EjjbnJJj+UjQthbaogus1SpUsbF8HH1e3bVu23btiFDhqisdZKcZVfPxt5akrPoHB01pxmBXCBn6D8El4/FSY0bNy4668GaFPbCu+qfNm2aJWcSM3F/pVkPwpniygFzvpGDw6e4E13UY5lGQGCjiLYETrDG0uNZXSnOvSIOX/369R26ieQsq7tYf+NJzhxGE18RAVsEcoScYTcSJpk0O6aPP/440hsncZTStl+FF5j6cj1xuXz5cqHE/ycRqMKcKeXgviZLQiZlqke7lfT7fsRqnebuTlN1sAeslVuudZKcpamjNXwLyZlvJ8aCOY1AjpAz9DHuFNJ57Q/4UOXKlSO1LfWYCLhsVGJL0uOHH34oNTVPnjxSjt0jAjRI2syPxx9/vF3xiPJfe+01DT886a4CMX6HDh164oknin1EcpbuTg/960jOxOHDNBFQRSB3yBkQQUiw0F2PncLWrVur9kH0cpgFMRMmMcdMztQbhSOioirLtH5y9sQTT9h1DfM9IYC1zhdeeAFXZmVMguTME3oUJjlT96WUJAL/RSCnyFnRokV/+OEHDe5y7Nix/4U4AalGjRpZciYjMwg5Q+h5Q49dQj85wxFRDR2dU1XgPlMc1G3VqpWnr77++uu1jQDYuae2UVgDAuknZ6+++qo2E8/Sik444QREBEjrv4girecUOYNhIwxY1P4IRMdys06Mw6pOnTp2tCmTH4SclSlTxlk53uonZz179oy6o3NTPwKkefpwkjNPcKVPmOTMj+fHLl1si8EVwlLh0qVL16xZU8r09whVLVq0aN68ecOGDWvXrl29enWccvenyrUUDq+lz7KNL7r55ptdEfAhkGvkDNYeyicb/SIlfvzxR0QC89ERkRaJlJxhT1ICyRmWlaWu4WMsCJCcxQJ7ciolOfPm2/Pnz48fe4dfERyhR9gbRFf3pvdPaRxSw5hctmwZ9pPiCp0/s////5jiAl0L/Y9pkjMJZ5XHUJgKoryq1JUQGfzhATIRhfPCgYNLLrkkIZ8pNiMHyRnOf0TRxdTpFQGSM6+IpUye5Ex0xS5peOpVq1a5Ri3HHMOcOXMQJAlMzkXjX19fdNFFGzZswCUqjRs3/uubvzxdeeWV3bt3/0tWsAeSMx/45SA5A0oIzRqFB8RlSj66QEORHCRnuCAhii6mTq8IkJx5RSxl8uknZzgZHooTv/TSS3HIHzE5VbQhkNK+fftWrlypvhY5cOBAnO7BzETLli1dq+jfvz8uhA5rgw7JmSvgZoHcJGf4e+Prr78O1wkiooR6bAtzR0Sak4PkrFChQuH2L7X5Q4DkzB9uqSmVfnIWyoGABg0a4KY/XGUDz6X4Y5AJL45zOio3N48aNSqz++T1119X0Y9YSv/6179mz56tHlTJQS3JmQM4dq9yk5wBDfyVEqL7Q6B87Iu3Azn2/BwkZ8Dc69b1EO2BqgwESM4MKHIzkX5yFnzmDPembdmyBeRp8eLF6r8WON2W4Vu41sO5VJcuXTKS+G/Tpk2dhY23N9xwA+RHjx5t5PhOkJz5gC5nyRmwmjlzZijuElvNLrzwQh/gaytSr149Y2xaJoKc1kzmgQBgG/rkaCjWkmtKSM5yrcel700/OQs+c2ZECceVeeq/Cueee67hzZs1a2ZXEKc+9+7dm5HE5Jz6+g4uB8yUclBuV6mUT3ImAaLymMvkDFeGY9VeciU+HocPH64CdYwy559/vjGKLRNByBlgtNQpZoZ++kcFTHyUj95kkXARIDkLF8+s05Z+chZw5uzkk0/GJH/GXT7yyCMqri0jg0twDSeLXTV2BadMmWKIrVu3zk7MnI8NZ5k9aqtXrza/9ZRDcuYJroxwLpMzINCvX7+Azm7JkiXOt4P76JTQi9xyyy3G8LRMvPvuu74rPf300y11ipkFCxb0rd93wVmzZgXsXBYPjgDJWXAMs1pD+slZwJmze++91/CVvskZNNSqVcvsKxHMzGB+kMF9z2YZhxzMtGXahr/vHcRcX5GcuUJkFshxcoYp3jVr1vj2fTt27MC14mZU9eSo36oO7mUMf8vE9OnTfbf51ltvtdRpZB48eNAcTNF3deoFdd7W5duEUl+Q5Cz1Xez8geknZwFnzhA7w/CVQciZ5TFPhEwzlCPhm5xNmDBB3fOaJUnOzJi45uQ4OQM+uPLF34Xohw8fPu+881wRjk5A8RgNNoCKw9MyPXjwYH/tBOvCeLfUaWR+9tln/pQHLHXfffc5/2zwrQYESM40gJzkKkjOnPwYDlpmlg4z7jIIOXv55ZfNNb3yyiuGI0bCNzn78ssvzcrVc0jO1LEyJEnOAAXiufjwbr4JjQF+8ITrTk3s1t+4caM4PC3TuMDD3Bhz+GizTLdu3SwVipmYwTIX1JDTu3dvH93KIuEiQHIWLp5Zpy395CzIsiYuTRJ9ZRByhphnZq+Kv4xF/b7JGRgkbpQy61fMITlTBEoUIzkDGriyArEwPHm9hQsXKk5ciWiHnm7fvr3Djx92mmLRVhyblmm7P4pwxtM5gM6ZZ55pHAOy1JzJRNSS0D9cRWG7du089SmFo0DAwT5VOtGTDC8+j6IHA+pMPzkLsqwp7dgNQs4s/fjOnTtF1+ybnEFJkGBRJGeeHFlGmOQsgwMCwaj7oK1btxYvXtwH2qEXwZ659evXL1q0CAus4rmESpUq4Qwpwk2LA9MubffziSXLGTNmVKtWzdxsvMKcmQozw2Uh6me3zRUFyWnSpIl6n1IyIgTsrCtIz9qVJTmLqBODqCU5szPXP/IRe1b0y48++qiT9F/fScQOvwR/ff/H008//STqD0LOEJLDrF8xh+RMEShRjOTMQAO75lV80KFDh3Bvo1Eq9gSu4siMvh9++AET2+hQ8CFxPDqnMVpFVid9zhVXXIHt/DiL3aZNGxwGKl++/BlnnIETACjlrNZ4e/XVV0s6tT1ixUClQykTKQIkZ5HCm3zl6SdnQZY1JXL24osvqvtH3FBu+FkkVMgZDh+o64ckYk0ZVSDsmaeyojDJmYiGYprkzACqRo0asENXZ2d5JsZQEkvioYceMkaQpwTmvLH06dzmqVOnetIpCiveFOLcAN9vsVrt2psUiBoBkrOoEU64/vSTsyDLmhI58zSzJZ2TVyFnmzdvVvenBQoUEL05yZndSMORWHVU1SVJzkSsRowYYYd/Jn/OnDmxRIUQG2lOY/cbaJA4jlTSWJRUOW2KbWeI5qqiUJLBzWyxhDcz8EFPYZrTuUP5NmoESM6iRjjh+tNPzkKcOYMPVedA+PEWfa4KOUPMs6JFixou0jmBVRJRv3rDzGo5c2bGxDWH5EyECERk06ZNds7um2++wWSMKJ+cNPjZ2LFjxaHknMbVRhgviu3HzW+SH3BWjre4jTfI/lHFhrmKYXegXW8yXw8CJGd6cE5sLeknZyHOnMF14ionV78GAfxhLXlhFXKGItgsrKIfMiNHjhSrIDmzG2OcOVO0qIBiF198sWUXwEoxAx1QedTFsUXMNXAGzkTjnlzwLU+NOeaYY8aPHy8OVYc03IvzMU9PVQcRxhYLy95kpjYESM60QZ3MikjOnDyYtKwJrwoHDT/uVOaII4oVKwYqJrlgRXL2xRdf4F4mZ/14iyNvP/74o1gFyZndACM5czWnsAQwS23uhT59+oSlP1I9GHc4eYpIH+JWzswQA297+OGHa9as6bsBiJ0xe/ZsMWiiOHiRv2DBgsaNG/vWH3rB+fPnm7uSOToRIDnTiXYC6yI5c3JrZnIGl4odJwiSZFesQoUK4qUChgtWJGeQf/zxx+2UZ/KxEDNv3jxDcyZBcmY3ukjOnM0pxLdly5bFlWJiRwSZtw6xYZ5UIWRgnTp1LrzwQoQZA2EKcZER+Fx77bWYSJs7dy6uFsUofuKJJ6677jqc5fTUQg3Czz77rNiPTOtHQCc5q1ix4mj+SxgC6Sdn4e45MygRIv7janPRS+KuwEGDBu3atcuQERPq5Ayl8Ge63Sl9BB+X7hXI1EJyZuc9Sc5EK406fdtttxkd8dVXXxUuXDjqGqk/CgQeeOABox+ZiAUBneQsChOizqAIxGJ2OisN8re75cyZSLlwvjLzF/Ann3xit2aRkfdEzlBk9erVWD8V4/7jdw5/dtvtjCE5szMqkrOgPsJLefxRgbGAvsDpFvWN815qoKwOBESSbTeymB8pAiRnOgw9yXVEal5JUB4pOROJmnPaKznLaEOkcoTHRBxzcLX9+/c7VEFyZmdsJGea/Q/2V+FC9F69emmul9WFiECnTp3sBhTz9SBAchaiPWelKj12FmMtES1rOvAky1fr1q0z24d0Q4BlQcVMkjM7GyM5Mxte1DmtW7eOugo9+hFvDLca4BL0UKrLnz9/9+7d8aOb/NXe888/325AMV8PAiRnoQy6LFaix85irCXIzNmpp56qyI1cxbDWY7YSXBrtWlBRoFy5cmb9ijmMc6YIlCjGOGciGqlM4/7yHTt2wHcdOHCgbdu2Ab8Ru0U/++yzjCfEbrwiRYoEVBhpcbi+GJ02qwYCJGeRWngWKE/9MAhCzo4//nhFbuQqhvNZZmvAdJprQUUBXBhg1q+YQ3KmCJQoRnImopHK9DvvvGO4R2wwUIlx44AD5swMbUjgBhEH4dhfYbJQbC3T+hEgOYt9FMTcAP02p7nGIMua6Bv8jatIj5zFcL+NuacR+si5lOJbRC03K1fPITlTx8qQJDkzoEhrAhcbiM6qcuXKQb5UuuFq2LBhQbRFXRYHOw4fPix+PtOaESA5i9rIk65fs8Hpry7IzBk6D4GIFBmSs1iTJk3MpjBw4EDnUopvp0yZYlaunkNypo6VIUlyZkCR1gRuxhT9VbNmzYJ86fTp00VtCG8WRJuGst99953YYKY1I0BypsHIE12FZoPTX11ActaoUSNFhuQghqOalrc+V6lSxaGU+quA+69JznwMUZIzH6BlV5EXXnhB9FeIYhik/ZiDF7Wdc845QbRpKLtmzRqxwUxrRoDkTIORJ7oKzQanv7qAy5roPEQyU+dJlpI9evSwMwLsRbMsop65YcOGo446yk6/Sj7JmQpKkgzJmQRI+h5vvPFG0V+9+eabvr+xdOnSoqpDhw7hfIBvbXoK4nvFNjOtGQGSMz12ntxaNBuc/uoCzpyh5+rWrWu+bk+dPH3wwQe4cMnOAk477TSE61TXZpa86qqr7JQr5pOcKQIlipGciWikMl29enXRXx08ePCEE07w96XSaYAVK1b406OzFDZLiJ/PtGYESM50WnsS69JscPqrC07O0G3dunUzsyKVnG3btrlenDd48GAVVZYyL730UnCrIjnzgSHJmQ/Qsq7I9u3bRZfl+4jl0qVLRT0PPvhg8qEYM2aM2GamNSNAcpb8MRJtCzUbnP7qQiFn6IN77rnHkh45ZCKMWcOGDVX6b9q0aQ567F7h7+98+fKp6HeWITlzxsfyLcmZJSwpy3zkkUdEl4V9Y8ccc4zXb8TUu6gEaVyi4FWJfvkBAwZIzeajTgRIzvTbfLJq1GltsdQVFjlDt11zzTXq65sIOFm1alXFzkYIpQkTJtiRMMt87AjxvcgitYrkTAJE5ZHkTAWlbJepVauW5LVuv/12rx81f/58UcnatWu9aohFHu5ObDbTmhEgOYvF7BNUqWaD019diOQM3YbA2WBFlmzJyPz555+xbIGLX7x2c8+ePXfv3m3osUtg78uoUaMQiMirfjt5kjM7ZBzySc4cwEnTqzlz5ohe65dffmnQoIH6B3bt2lUsjrTD8SB1tRokmzdvLrWcjzoRIDnTYOSJrkKntcVSV7jkLNOXZ5xxxmOPPfb555+L/AkXk+MHu3///kFu4itVqhSIHS6NETUbadzFOXny5Bo1aoRrUiRnPvAkOfMBWjYWwXDD4UrRd+3atat+/foq39K0aVMc9xHLIqpOwLPVKvWGIoOzSmLLmdaMAMlZKGacxUo0G5z+6oKH0nDoXZyHR6wyTKdVqFAhxKksrHKeffbZt912Gygg2Nj48eMRYwl/yEZ0/J7kzKGL7V6RnNkhk758TFRLjguz1xieDvvP4A169+6NaTap4GWXXZYt+JQpU0ZqPB91IpBr5AyhQPFLhPg1Y8eOnTRpEsK/Dx06tFWrVrhEMcQhU6JECfy5hYPYSISoNhJVOq0tlrqimDmLpCfiU0py5gN7kjMfoGVpEfyxhIA4Zve1efPm++6779xzzzV2f4Ku1a5du2/fvphWN8tPnDgxixDAt5g/gTnaEMgdclaoUKE77rhj48aNxhqRmDhw4ABmWIKcocH4bdeu3csvv4xLL0TNWKFCZtu2bRM6ma3N1OKqiOTM9feA5MwVIrMAyZkZkxTnYK+CJd8y3BqcPqbTjEdzAgYTytlqnSBjAdf8IczRg0COkLPLL79869atImeySz/33HMFChTwav+dOnWyo31GRTiF3bFjR6+aI5fXY2cx1kJy5mpDJGeuEJkFSM7MmKQ7B1H+P/roI3+u7P3338f0QNbhg4Ol/r6XpYIjkAvkDDE+sfRvkCTXxOrVqytVqqQ4jo499tjnn3/eVachAPKXN29eReU6xILbUMI16CFn2J/RuXNnT8e4HHq3WLFid9999+jRo6tVq+YgFtYrkjMfSJKc+QAt24vAd2Mb6G+//abu9H7//ffHH38cCyvZ+O24uU79SykZLgKpJ2f33nuvQYzUE5s2bSpZsqTraEK0hGXLlqmrzUjijygfk3OujfEpEK49JVCbBnIGToZzlJlvB6ny2RN/FoNVGXckI6waDhz8+Saq/5Oc+UCW5MwHaOkoguBn2AQDb+7q7t5++22c7M7er54+fbrrN1IgIgTSTc7at2/vlTkZ8qBQDmdxMNxwtuCNN94w5D0lZs6cieKJGLMRGVZy1GogZ4sWLTK+9/Dhw5j3CtK1HTp0MLQhgZNiQbSplCU5U0FJkiE5kwDJtUec9sLPJ9ZNsPa3Z8+ezJhF0I0tW7YgNNqdd9558sknZzsmjz76qOiLmNaJQIrJGaa+sJ3RE2eShIcNG+YwuDD6JHlPjzid4KBc3yud1hZLXZGG0sj006pVq8RPC7i4iagZoraHHnooamsgOfOBMMmZD9BSXAR/ymfp2qVDp+BHTvRFTOtEIMXkDLeieWJLZmH8OWQ3CVKuXDnEgTcXUc/Zu3dv2bJlHcaFplc6rS2WujTMnEkxxK+66qognff000+LQCGcUhBtKmVJzlRQkmRIziRA+Jg+BLp37y76IqZ1IpBWcoZonaBWGaqEc5RTpkwZM2YMbi/EYqWnwwGI92454qBNnYfZSSIUvKVyrZk6rS2WujSQM+zcFz8NW4aDdOGnn34qamvRokUQbSplSc5UUJJkSM4kQPiYPgQuvfRS0RcxrROBtJKzNm3agBIhcCACBEpDBufqEHFdkaJBg1Qcj3ny5Pn222+hH+E5EMwWcTpOP/30Jk2aXHfddbNnz1bUjOLbtm2DKrN+rTk6rS2WujSQM0QxFj8NV5777kJEs8TxLkMb0kZ8S986XQuSnLlCZBYgOTNjwpyUIYDTDIYvYkIzAmklZw8//PCTTz7pEPcVfxLggJ3dnJaRj4vRzFfm4F41MLARI0aYX2Fs4kayjz/+2NDgnIBwzMNZs8Hpr07DnjPcL4FzAOKn+Q6Bcc0114h6MIumwT5IznyATHLmAzQWyS4EcCud6I6Y1olAWsnZkCFDXEdBZnbNmTzhbd26dSVV1157LaLOSpniI0jb4sWLXTVDAKrEgjGkdVpbLHVpmDlDt61cuVL8Ovxx4K8v33nnHVEPTkv50+OpFMmZJ7gywiRnPkBjkexCAFcaiO6IaZ0IpJWcKa4FzZgxw5VCmRdGmzZt6jrEihYt+s0337gqVyGRrnUFEtBpbbHUpYec4TY98eswK1u8eHGvHYO/A8Q1TShs3LixVyU+5EnOfIBGcuYDNBbJOgRwck30bExrQyCt5ExxCGCvmCt/atasmaRN8cQ0Trq4Kr///vsl5boftZlaXBVpWNZEn+FYL+5nFb/xmWee8dqXCxcuFDV8+eWXeqLhkZx57SnIk5z5AI1Fsg6B9evXi06JaW0I5Dg5w0jBrnxnCmUmZ4rjK3/+/Pirw1n5Aw88oKgtKjFtphZXRXpmztA9iN0ifePVV1+t3m09e/aUimsbnCRn6t1kSJKcGVAwkWIEli5dKvklPupBQJv/T6z1vvfee878yTc5wye77jwjOYvczrWRM6xk//jjj+L3IFy4Ysyzli1bQlgsiz9YHc6zhDucSM584Ely5gM0Fsk6BLDyIPolprUhEBE5wyEPnFeL+h8OyQU3dQS/iI6cud6JTnIWuanrWdbMGGLXrl3N34P1zVKlStlZKuZXcf+rdNgTSpo3b25XJPR8kjMfkJKc+QCNRbIOAcSdMvs05mhAICJyhvBgGhrfu3fv4KYeKTlD2Ftn5kdyFrmdaJs5y9gi+Lj5kxB5BRyxR48eiBtUvnx5cLVTTjnlkksuQZS8HTt2mOXHjRsX0LLz/t8/RSUkZ4pAiWIkZyIaTKcVAVxiaHZQzNGAAMkZyZkGM4uzCs3kDFfszZs3L8gHL1iwAEqC+Hr4U5xO2L9/P+bkVPSQnKmgJMmQnEmA8DGVCNxwww1BvBnL+kaA5IzkzLfxZEdBy2XNAgUK1KpVy/V+hjPPPBP3Zi5atOjiiy9W97zHHnvsiy++6A+duXPnYqFTvS6z5HnnnSdWfdFFF5llpBySMwkQlUeSMxWUKJPtCLRu3Vr0J0xrQ4DkjORMm7HFU5F55qxdu3b79u1DaxCrolKlSnbeE/fSZ8QgiT1hWJG0k7TMx2Veu3btUv9mLH1ixsuVL1rWJWaOGjVKrHTw4MHiW8s0yZklLM6ZJGfO+PBtOhA466yzRH/CtDYESM5IzrQZWzwVSeQM81J79uwxmjJ9+nQ7H3rzzTcbYkj4CEmHEytYWESwFlGPOY07wjDTdtJJJ9m1xFM+pvrEKjp37uxanOTMFSKzAMmZGRPmpA+Bk08+WfQnTGtDgOSM5EybscVTkbSsWbFiRbEd3333nZ0/xWkOUXLgwIF2ks75Rx55JG6ZwF0Qb7zxBm5d3bJly9atW9etW/f+++/jJFSXLl0KFy7srMHT202bNonNVrm9leTME8IZYZIzH6CxSNYhULBgQdGfMK0NAZIzxNhzPlAZJM7Z5MmTnZXztGbkpi7NnCF4GBYQxVoR3N/SY4I8iWI4XGkplqhM3OoqXgCF1VjcjufaQpIzV4jMAiRnZkyYk0oEpLtPRK/IdHQIkJxt377dmT+dc845vkcclpiclY8YMcK38nAKRmdbCdEskTOgtmbNGrFtdh28e/duUaxy5crhIB6llvr164tt3rBhg0ptJGcqKEkyJGcSIHxMKwLSZLzoYZiODoEcJ2c4jedMnvC2Zs2a/gYd7kV0vRtq0KBB/pSHVio620qIZmlZE8BNmzZNbBvuTTKjWbp0aVEGYSkcrrlE5ItChQqZlejPwSKp2OxZs2aptIHkTAUlSYbkTAKEj2lFQE/YUtFxMQ0E0krOsFCuMlIwq+JMzvCjjGieKqrMMk2bNnVWjrcdOnQwF9Sak/phYJ45w/Yv8asff/xxM+IXXHCBKPPRRx+ZZTI5CC2LOTYsJmJTf8D4ZHZVqOdLRzVHjhypUpbkTAUlSYbkTAKEj2lFAJtlRWfItB4E0krOsLzjMNORGUQgRq7kacmSJeYRp8L8UPuyZctc9eMojFm/1hw9dhZjLWZyJkXusezjPn36iG1+7rnnLHsFkTh+++03Q/LWW2+1FPORCZ6HDWReCyJMmtEYJBRv9iQ584oz5EnOfIDGItmIwMSJE0WvwrQeBNJKzo4++ujHHnvMYSLjiiuu2Lt3ryt5Qnhky9HkGvcAUxiuyv/9739bKteaqcfOYqzFvKxZtWpVsT07d+40I654VPPyyy8XVc2fP9+sykcOzO6nn37CVeiTJk2CKatrkHaH1KtXT6UsyZkKSpIMyZkECB/TigDiAYlejmk9CKSVnGGY3HXXXViMuvDCC6UpNEx24CcPJ/ZcyROuPbSbvOjVqxcisVsORsSHR4QEV+UQCHGexbIlSpl67CzGWswzZ4ht4XpgU/GoZp06dcRP27x5sxLojkK4dlM8cWn394FZh8NRTURcQ8wYfPWKFStOPfVUqSzJmQSIyiPJmQpKlEkBAlLER9HjMR0dAikmZ9iijd9KcCAElnr55ZexswhrU/htUqFNGRkwMLuRhet/oPaVV15p0aKFQeBwnu+WW27ZuHGjShVom1HQrhYd+dHZVkI0m2fOAKvrgU3pqKbdRQJg4uKyJj5ZZc3buV9vvPFGETqsKTjLG29PP/10saB4VPOZZ54xXsFApb9XSM4MDNUTJGfqWFEyqxFo27at4T2Y0IZAiskZhkOTJk0QokWFKpllMOEi/YRJ4wsLo0Yp3NODowPGo0qiVatWksJ4HrWZWlwVmWfOALR0YFOi4WXKlBFb63xUE+FkReEGDRoE7EjMu4oK77zzTkWFHTt2FAtiG69RcP369eIrTM4Zr5AgORPRUEyTnCkCRbFsRwDBhkTvwbQeBNJNzjAocJZOhSpJMgjkjrkx1zE1fPhwqaDi49133+2qXJOAHjuLsRbLmTPcOCk2ady4cSLc559/vvh25cqV4lspPXPmTFG4W7dukoD02LBhQ1y7WaFCBSnfeHznnXdEheosXroFr1+/foZOiUlI97iTnBlAqSckSMUuU0+LBFq9akoSAZ0IVKtWTd2kKRkWAqknZ7DhTp06/fzzz4q0CWLz5s0rUqSIovFjtVRdc0by0UcfVVSuQywsS0qsHsuZM+nAJrpcxBq3LYmf89RTT4lvpbS0W3bMmDGSgPgIVp7RfPDgQbtNi99//71Ye5UqVUQNzmnspsyUXb58uRh6zcjPvJV2O5KcOaNq+ZbkzBIWZqYPAexYFT0S03oQyAVyhsGCXx9MhrmyqH379g0bNixPnjyexhd2bKOgq3IIQCwiwD01+C/CeuwsxlosZ86KFy+O/jBahbO1IiiwACPuIvqsdu3a4lspjYgshh4kHA5s4g5N3HFuCL/77ruSKjyWKFHCEEAC8ji+YBZzyAGZw9KqtCQ/YMAAUe0TTzwhaiA5E9FQTJOcKQJFsRQgIHpL0ZMwHR0CEXEF46ctupZDc+/evdXNHj+47du3f/vtty3PaeIq6gcffLBcuXLqCkVJ/CBOnz7dUnOGtOHVSy+95GkSRNQfYTrSHkqCcktyBkC7du2KYCpoIWwCfxpKEOOwBkg3qLprJDr1A5uNGzcWAfnkk0+kSvEoTdpZyphLuea0bNlSrHrx4sViEZIzEQ3FNMmZIlAUSwEC+IEUHQjTGhDIHXJmDJCiRYs2b978pptuwk7r/v37I04nfl6Nt0ES5cuXR+xSLKOtXbsWwbPwDwk8ItM37QvSHqWyGows3irsyBnQAQPD3n8lmOyFcIOEdGBTXE8Uy+GeKBEKsHXxbSYtHdWcOnWqWcZHTo0aNcSqcdJYVEJyJqKhmCY5UwSKYilAAAtPogNhWgMCOUjOUjBSwvwEDUYWbxUO5CwsHL/88kvxG3Fjq6VmhEUWxXCLlFnM91FNsyoxBwzy8OHDRu2Io5Y/f35DAIdfcCArrf9wTarxpSEmSM5CBJOqEo4AduUa3oMJPQiQnCV8UETePD12FmMtlgcCwoUVVYgfaHdgE+unohhisZibgY1oogyWIyUZLLNiQRbnp6R810eENxM1hzVd7FpvWgVIztLas/wuMwJioETRjTAdHQIkZ2Y7zK2c6GwrIZo1kDPpwKbddePfffediIklwZKOakoBydq1a5dZQsXUV5cuXTxZ6oIFC8Tar7zySk/FKSwhEAo5kzb/SVXwkQgkBAGSM9F56kmTnOH+zXz58oU1BHAKEEGyTjzxxLAURq5Hj53FWIsGctasWTPxAxFbz9xtxYoVE2VwQuSoo46SxKSjmgi3IR3VFGe/tm/fLh3JlLRJj4jgIjbgjjvukAT46AmBUMgZzKBmzZqe6qUwEdCMQMmSJaUbU0RPwnRECOQ4ObvtttvwC4jJiEceecTTL53l6MA9ngiohp5C+IWLLrrIUiZxmREZVnLU+iZn1atXx3lgcG3XPsNJYFwNlvnkRYsW4U4ncxHpqCbujzLLSFFkzUc1JReJEyhmJXY50gV5zz77rJ0k81UQCIWcwWZwu8jAgQMxJ8p/RCCBCOBGQlw1mBx/njstyWVyJl1F2LlzZxWf7CAjzmsg7SCZoFept3V/5Kxp06ZGTLKxY8eqdFjVqlVr1aplJ6lyVBNkS9y2b45WLIWowalju+rM+fjTQezr9957zyzDHHUEwiJnYqcwTQSIABHIIJDL5AxxrEQzwDk5dc9slsR0iajt0KFDZpkk5oiNTmXaHzkT71DCBi8sVwfsPFwlIcJreVQTVSA/w88+//xz8zHDp59+WlQiXtDk2rxSpUrhQ4zir7/+umsRCjggQHJm2BITRIAIhI5ALpOztm3binhiw7SDK3Z9hUUwUdu2bdtciyRCQGx0KtP+yBnu0xTRaNKkScDeWrJkiagQ90fZKURMvEaNGkm7zTLCWIYXlWCXrp0Sy3ycVMgU37Nnz2mnnWYpw0xFBEjORFNkmggQgXARyGVyVq9ePRHMb7/9VtEtW4rh+gFRG2ZeLMUSlyk2OpVpf+Rs2rRpIhqWe/w99aV0VBNroJ6KZ4SlpckPP/zQq5L69et37NgRs2heC1JeQoDkTBwgTBMBIhAuArlMznBO88CBAyKemP2SPLD647hx40RVuAxKvWyckmKjU5n2R86MG8ozmIwePTpIJ2EaTMTWx42ZmdoxqSbqwfGT4MdYgnxXLpclORNNkWkiQATCRSCXyRl+Wd58800RT6wa+f65Wbdunajqsssu861Ka0Gx0alM+yNnnTp1EtGYNWtWwF5BRCtD4YwZM3xr++mnnww9SFSoUMG3KhYMggDJmWiHTBMBIhAuAjlOzgYNGiTiiU3Y/tx1w4YNRT3Y0o0bPP2p0l1KbHcq0/7ImXSU94svvgjYMYgVhJZg7RwLpuZ71tWVL1++XOymFi1aqJelZIgIkJyJdsg0ESAC4SKQ4+QMQdolPBFP1IcDl+In4w4eH0riKSJ9f/oe/ZGzwoULi1D8+uuv5pixzh2GBUdctYRrK53FLN+iIJbJcTbTvG1/0qRJYsM8Hdi0rIuZ/hAgORPtkGkiQATCRSDHyRnc8rJly0RIEfjT668w9nZn7tQx9HTv3t2fw4+hlNHotCb8kTP0hLSFH4RJvXsQBQOWBEj37t3boUMH9YKQxG79nTt3Zrpj//79VapUEYsHPLApqmI6CAIkZ2n1GPwuIpAEBEjOENNA6ogRI0aoO20Eh1+6dKmo4YcffsifP7+6hpglxaanMu2bnCFMqwiIpwVEMSAZ7qAoUqSIejf36tVLrBdHE8Sy0oHNFStWiG+Z1oYAyZlopUwTASIQLgIkZ3Dm0uQZEMbvo6KTlw5pomyQUwWKlYYpFq49JVCbb3ImEix8V58+fdRxl0zKU5i0xx57TIRRugSzbNmy4tstW7aot4qSISJAcibaIdNEgAiEiwDJGdw1Yj8hoL8ELFhXwYIFHZw5dhNNmTJFKrVhwwbLmxUd9MT8SvqA9D36Jme48VBEw9MNEkHCpC1cuFCst02bNpKJiAFycaen9JaPehAgOROtlGkiQATCRYDkLOPJBwwYYAYWm46GDx9eu3ZtrF2KDh/hC/r3779jxw6pCHae4dimKJkFaekb0vfom5y1atVKROOtt95S784gYdI2bdok1gv7k+qtVKkSrlfH8jn+OChUqJD0lo96ECA5E62UaSJABMJFgOQs48lxtA5zEHbYItgnomwgiMGqVauMvdpmYU8LX3p+QdxrMX9GynJ8k7OaNWuKUKxZs8YdzT8lfIdJy5s3r3j3OW7DzKYNjH9+fi78n+RMHB1MEwEiEC4CJGfG7wimxyZPnuwb3qFDhxqqsinh+4OzpaBvcgaD+Oqrr4zPHDNmjHq/+g6TJjHCzZs3q1dKSZ0IkJwZQ4MJIkAEQkeA5Ezy5z179pTudHLFHOEOcF2hpCdrHl0/L9sFfJMzdCGipOASCVxiP2HCBE97Cc1h0iwvMjdbSZC1VLM25kSHAMlZtnsGtp8IJBkBkjOz98YFhljiRNhR147DotPUqVMrVqxoVpI1Oa4fme0CQchZkF70FyYtyCmEIK1lWa8IkJxlu2dg+4lAkhEgObPzySVKlMAesrlz5+7atUvqQUytffjhh0OGDDnppJPsimdNvvRt6XuMi5z5C5MWJH5H1thcKhpKcpY+X8EvIgLJQYDkzPWHAmcFcFFmjRo1GjRogJNzlStXVlyhctWcCIHk2GJELYmLnPmjWf4oXSIsKccaQXIW0YClWiJABIAAyVmO/aSYPjf1wyAucuZvgdLfYqipV5kROQIkZ6l3HfxAIhAjAiRnkTvxhFcQo/HpqTouciZt7cfBAldLQNAyERMft627VkGBsBAgORNtlWkiQATCRYDkLCxfna16wrWnBGqLi5zVqlVLRGP9+vWuJlKyZEmxyNq1a12LUCAuBEjORFtlmggQgXARIDmLy7cnpd5w7SmB2uIiZ0cfffTu3bsNQObNm6fS5RAzivTt21elCGViQYDkzDBUJogAEQgdAZKzWBx7gioN3aSSpjAucoY+btmy5d69ewEIYslWr15dpdePO+64kSNHos1du3ZVkadMXAiQnCVtpLM9RCBNCJCcxeXbk1JvmqzZ8ltiJGfoY2wjq1u3Li5lSkp/sx0hIUByZjncmEkEiEAoCJCcheSqs1ZNKGaUZCXxkrOstQs23AUBkrMkj3q2jQhkOwIkZy4uOPWvs92CXdtPcpZ6G47lA0nOXIceBYgAEfCNAMlZLI49QZX6Np1sKUhyliBrS1FTSM6yxQOwnUQgGxEgOUvRz4WvT8lGq/XU5ldffdUXMBaFevTo0Z//kooAeseizyLLIjnzNAwpTASIgCcESM4ic95ZotiTuWSjcIjkbOPGjdmIQI60Gb2jc8yRnOWIXfEziUAsCJCc6fTnSawrFrPTWSnJmU60Y6yL5CxG8Fk1ESAC4SJAcpZEwqSzTeHaUwK1kZwlsFOiaBLJWRSoUicRIAKxIEByppMIJbGuWMxOZ6UkZzrRjrEukrMYwWfVRIAIhIsAyVkSCZPONoVrTwnURnKWwE6JokkkZ1GgSp1EgAjEggDJmU4ilMS6YjE7nZWSnOlEO8a6SM5iBJ9VEwEiEC4CJGdJJEw62xSuPSVQW4hxzvDzn8APZJMyCGzatEnnwOFpTRoeESAC0SFAcqbTnyexruhsKyGaOXOWkI6IuhmcOYsaYeonAkRAGwIkZ0kkTDrbpM3U4qqI5Cwu5DXXS3KmGXBWRwSIQHQIkJzpJEJJrCs620qIZi5rJqQjom4GlzWjRpj6iQAR0IYAyVkSCZPONmkztbgqIjmLC3nN9ZKcaQac1REBIhAdAiRnOolQEuuKzrYSopnkLCEdEXUzSM6iRpj6iQAR0IYAyVkSCZPONmkztbgq4p6zuJDXXC/3nGkGnNURASIQHQIkZzqJUBLris62EqKZM2cJ6Yiom8GZs6gRpn4iQAS0IUBylkTCpLNN2kwtropIzuJCXnO9JGeaAWd1RIAIRIcAyZlOIpTEuqKzrYRo5rJmQjoi6mZwWTNqhKmfCBABbQiQnCWRMOlskzZTi6sizpzFhbzmejlzphlwVkcEiEB0CJCc6SRCSawrOttKiGaSs4R0RNTNIDmLGmHqJwJEQBsCJGdJJEw626TN1OKqiOQsLuQ110typhlwVkcEiEB0CJCc6SRCSawrOttKiGbuOUtIR0TdDO45ixph6icCREAbAhGRs44dO/aN/l+9evWSSHeyq03aTC2uijhzFhfymuvlzJlmwFkdESAC0SEQETnLLn6S062NzrYSopnkLCEdEXUzSM6iRpj6iQAR0IYAyVlOM7Mjjvhf+z9jqCPUq6AAAAAASUVORK5CYII="}}},{"cell_type":"code","source":"import shapely.wkt\nfrom shapely.geometry import Polygon\n\nwkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10))\"\npolygon = shapely.wkt.loads(wkt_string)","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.1522Z","iopub.execute_input":"2022-06-28T06:48:41.152679Z","iopub.status.idle":"2022-06-28T06:48:41.251535Z","shell.execute_reply.started":"2022-06-28T06:48:41.152586Z","shell.execute_reply":"2022-06-28T06:48:41.250667Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"In this article, we will first look at how to visualize a polygon with the Shapely package or the Matplotlib library. \nThen we will go over five feature engineering ideas from polygons in WKT format.\n\n# How to Visualize a Polygon\n\nThe first thing, you might want to do with the polygon is to visualize it to get a better intuition about it. \nYou can either plot the polygon directly via the Shapely package or you can plot the polygon via its coordinates using the Matplotlib library.\n\n## Visualization via Shapely Package\n\nTo visualize the mere shape of the polygon, you can display the Shapely polygon after loading it.","metadata":{}},{"cell_type":"code","source":"wkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10))\"\npolygon = shapely.wkt.loads(wkt_string)\npolygon","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.25293Z","iopub.execute_input":"2022-06-28T06:48:41.253416Z","iopub.status.idle":"2022-06-28T06:48:41.268229Z","shell.execute_reply.started":"2022-06-28T06:48:41.253383Z","shell.execute_reply":"2022-06-28T06:48:41.266991Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"wkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"\npolygon = shapely.wkt.loads(wkt_string)\npolygon","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.270005Z","iopub.execute_input":"2022-06-28T06:48:41.270469Z","iopub.status.idle":"2022-06-28T06:48:41.278656Z","shell.execute_reply.started":"2022-06-28T06:48:41.270426Z","shell.execute_reply":"2022-06-28T06:48:41.277747Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"While this is a quick option, its disadvantage is that you don’t get an intuition about the coordinates.\n\n## Visualization via Matplotlib Library\n\nTo visualize the polygon by its coordinates, you can use the Matplotlib library in addition to the Shapely package.\n\n","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nplt.rcParams.update({'font.size': 14})\nplt.rcParams[\"figure.figsize\"] = (4,4)","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.280263Z","iopub.execute_input":"2022-06-28T06:48:41.280785Z","iopub.status.idle":"2022-06-28T06:48:41.289084Z","shell.execute_reply.started":"2022-06-28T06:48:41.280752Z","shell.execute_reply":"2022-06-28T06:48:41.288157Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"From the Shapely polygon, you can retrieve the polygon's x and y coordinates from the `xy` attribute of the exterior (`polygon.exterior.xy`) and interiors (`polygon.interiors[i].xy`). The 'exterior' is the outer shaper of the polygon. Additionally a polygon can have none, one or more 'interiors', which are smaller polygons within the exterior.\nwith their `xy` attributes to plot the polygon as follows:\n","metadata":{}},{"cell_type":"code","source":"def plot_polygon(wkt_string, ax=None):\n    polygon = shapely.wkt.loads(wkt_string)\n    \n    # Retrieve and plot x and y coordinates of exterior\n    x, y = polygon.exterior.xy\n    ax.plot(x, y, color = 'black')\n    \n    # Retrieve and plot x and y coordinates of interior\n    for interior in polygon.interiors:\n        x, y = interior.xy\n        ax.plot(x, y, color = 'black')\n        \n    ax.set_title(wkt_string.replace(\"),\", \"),\\n\"), fontsize=14)\n    ax.set_xlim([0,100])\n    ax.set_ylim([0,100])","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.290298Z","iopub.execute_input":"2022-06-28T06:48:41.29075Z","iopub.status.idle":"2022-06-28T06:48:41.304731Z","shell.execute_reply.started":"2022-06-28T06:48:41.290721Z","shell.execute_reply":"2022-06-28T06:48:41.303836Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(12, 6))\nwkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10))\"\nplot_polygon(wkt_string, ax[0])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"\nplot_polygon(wkt_string, ax[1])\n\nplt.tight_layout()\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:41.306786Z","iopub.execute_input":"2022-06-28T06:48:41.307194Z","iopub.status.idle":"2022-06-28T06:48:41.681034Z","shell.execute_reply.started":"2022-06-28T06:48:41.30716Z","shell.execute_reply":"2022-06-28T06:48:41.679724Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 1. Find the Area of a Polygon \nAfter you have visualized the polygon, you might want to know how to calculate the area of the polygon from its given coordinates.\nInstead of writing your own function to calculate it, you can simply retrieve the polygon’s area from the Shapely polygon’s attribute `area`. \nLet’s plot a few polygons and verify their areas. \n\nBelow, you can see a quadratic polygon with an edge length of 80 units. The Shapely polygon’s `area` attribute returns a value of 6400, which corresponds to 80 times 80. And is, therefore, correct.","metadata":{}},{"cell_type":"code","source":"polygon.area","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:41.683344Z","iopub.execute_input":"2022-06-28T06:48:41.683968Z","iopub.status.idle":"2022-06-28T06:48:41.69049Z","shell.execute_reply.started":"2022-06-28T06:48:41.683877Z","shell.execute_reply":"2022-06-28T06:48:41.689431Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_polygon_with_area(wkt_string, ax=None):\n    polygon = shapely.wkt.loads(wkt_string)\n    \n    # Retrieve and plot x and y coordinates of exterior\n    x, y = polygon.exterior.xy\n    ax.plot(x, y, color = 'black')\n    \n    # Retrieve and plot x and y coordinates of interior\n    for interior in polygon.interiors:\n        x, y = interior.xy\n        ax.plot(x, y, color = 'black')\n    \n    ax.set_title(f\"Area: {polygon.area:0.0f}\", fontsize=14)\n    ax.set_xlim([0,100])\n    ax.set_ylim([0,100])\n    \nfig, ax = plt.subplots(nrows=1, ncols=2, figsize=(12, 6))\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10))\"\nplot_polygon_with_area(wkt_string, ax[0])\n\nwkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10))\"\nplot_polygon_with_area(wkt_string, ax[1])\n\nplt.tight_layout()\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:41.691757Z","iopub.execute_input":"2022-06-28T06:48:41.692605Z","iopub.status.idle":"2022-06-28T06:48:42.012531Z","shell.execute_reply.started":"2022-06-28T06:48:41.692571Z","shell.execute_reply":"2022-06-28T06:48:42.011841Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"However, not all polygons are closed shapes. Sometimes, polygons can have holes, which are called interiors in the Shapely package.\nIf we plot and verify their area, we can see that the area of the polygons with interiors is smaller than the same polygon without any interiors because the area of the interior is subtracted from the area of the exterior.","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(nrows=1, ncols=3, figsize=(18, 6))\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10))\"\nplot_polygon_with_area(wkt_string, ax[0])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"\nplot_polygon_with_area(wkt_string, ax[1])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20), (60 60, 80 60, 80 80, 60 80, 60 60))\"\nplot_polygon_with_area(wkt_string, ax[2])\n\nplt.tight_layout()\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:42.014106Z","iopub.execute_input":"2022-06-28T06:48:42.014902Z","iopub.status.idle":"2022-06-28T06:48:42.587739Z","shell.execute_reply.started":"2022-06-28T06:48:42.014866Z","shell.execute_reply":"2022-06-28T06:48:42.586781Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 2. Find the Perimeter of a Polygon\nNext, you might want to know how to calculate the perimeter of the polygon from its given coordinates.\n\nLet’s plot a few polygons again and verify their perimeters. \nBelow, you can again see the quadratic polygon from our previous example with an edge length of 80 units. The Shapely polygon’s `perimeter` attribute returns a value of 320, which corresponds to four times 80. And is, therefore, correct. \n\nAgain, some polygons have interiors. If we retrieve the perimeter for a polygon with interiors, the perimeter increases, because the perimeter of the interior is added.\nYou can create new features for the outer and inner perimeters as follows:","metadata":{}},{"cell_type":"code","source":"perimeter = polygon.length\nouter_perimeter = polygon.exterior.length\ninner_perimeter = perimeter - outer_perimeter","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:42.588923Z","iopub.execute_input":"2022-06-28T06:48:42.589301Z","iopub.status.idle":"2022-06-28T06:48:42.595181Z","shell.execute_reply.started":"2022-06-28T06:48:42.589269Z","shell.execute_reply":"2022-06-28T06:48:42.594019Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_polygon_with_perimeter(wkt_string, ax=None):\n    polygon = shapely.wkt.loads(wkt_string)\n    \n    # Retrieve and plot x and y coordinates of exterior\n    x, y = polygon.exterior.xy\n    ax.plot(x, y, color = 'black')\n    \n    # Retrieve and plot x and y coordinates of interior\n    for interior in polygon.interiors:\n        x, y = interior.xy\n        ax.plot(x, y, color = 'black')\n    \n    perimeter = polygon.length\n    outer_perimeter = polygon.exterior.length\n    inner_perimeter = perimeter - outer_perimeter\n    \n    ax.set_title(f\"Perimeter: {perimeter:0.0f} \\nOuter Perimeter: {outer_perimeter:0.0f} \\nInner Perimeter: {inner_perimeter:0.0f} \", fontsize=14)\n    ax.set_xlim([0,100])\n    ax.set_ylim([0,100])\n    \nfig, ax = plt.subplots(nrows=1, ncols=3, figsize=(18, 6))\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10))\"\nplot_polygon_with_perimeter(wkt_string, ax[0])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"\nplot_polygon_with_perimeter(wkt_string, ax[1])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20), (60 60, 80 60, 80 80, 60 80, 60 60))\"\nplot_polygon_with_perimeter(wkt_string, ax[2])\n\nplt.tight_layout()\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:42.596752Z","iopub.execute_input":"2022-06-28T06:48:42.597895Z","iopub.status.idle":"2022-06-28T06:48:43.064767Z","shell.execute_reply.started":"2022-06-28T06:48:42.597848Z","shell.execute_reply":"2022-06-28T06:48:43.063661Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3. Get the Number of Interiors of a Polygon\nAs you have already seen, polygons can have so-called interiors. \nThese are the holes in the exterior polygon. \nThe Shapely package provides an attribute for the number of interiors a polygon has.","metadata":{}},{"cell_type":"code","source":"num_interiors = len(list(polygon.interiors))","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:48:43.066007Z","iopub.execute_input":"2022-06-28T06:48:43.067045Z","iopub.status.idle":"2022-06-28T06:48:43.072327Z","shell.execute_reply.started":"2022-06-28T06:48:43.067007Z","shell.execute_reply":"2022-06-28T06:48:43.071112Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_polygon_with_num_interiors(wkt_string, ax=None):\n    polygon = shapely.wkt.loads(wkt_string)\n    \n    # Retrieve and plot x and y coordinates of exterior\n    x, y = polygon.exterior.xy\n    ax.plot(x, y, color = 'black')\n    \n    # Retrieve and plot x and y coordinates of interior\n    for interior in polygon.interiors:\n        x, y = interior.xy\n        ax.plot(x, y, color = 'black')\n    \n    num_interiors = len(list(polygon.interiors))\n    \n    ax.set_title(f\"# Interiors: {num_interiors}\", fontsize=14)\n    ax.set_xlim([0,100])\n    ax.set_ylim([0,100])\n    \nfig, ax = plt.subplots(nrows=1, ncols=3, figsize=(18, 6))\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10))\"\nplot_polygon_with_num_interiors(wkt_string, ax[0])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20))\"\nplot_polygon_with_num_interiors(wkt_string, ax[1])\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20), (60 60, 80 60, 80 80, 60 80, 60 60))\"\nplot_polygon_with_num_interiors(wkt_string, ax[2])\n\nplt.tight_layout()\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:43.074225Z","iopub.execute_input":"2022-06-28T06:48:43.074908Z","iopub.status.idle":"2022-06-28T06:48:43.528894Z","shell.execute_reply.started":"2022-06-28T06:48:43.074861Z","shell.execute_reply":"2022-06-28T06:48:43.52787Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4. Check if a Polygon is Invalid\nPolygons can be invalid.\nFor example, if a polygon’s interior intersects with the exterior or if the interior lies outside of the exterior. \nWhen you plot a Shapely polygon, the package indicates whether the polygon is valid or invalid with the polygon’s coloring. \nA valid polygon is filled with green color, while an invalid polygon is visualized in red.\nA new feature can be created from the validity of a polygon.\nFor this, you can use the boolean attribute `is_valid`.","metadata":{}},{"cell_type":"code","source":"polygon.is_valid","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:43.530521Z","iopub.execute_input":"2022-06-28T06:48:43.531173Z","iopub.status.idle":"2022-06-28T06:48:43.538625Z","shell.execute_reply.started":"2022-06-28T06:48:43.531127Z","shell.execute_reply":"2022-06-28T06:48:43.537397Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"wkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10))\"\npolygon = shapely.wkt.loads(wkt_string)\nprint(f\"Polygon is valid: {polygon.is_valid}\")\ndisplay(polygon)\n\nwkt_string = \"POLYGON ((10 10, 20 10, 20 80, 90 80, 90 90, 10 90, 10 10),(60 60, 70 60, 70 70, 60 70, 60 60))\"\npolygon = shapely.wkt.loads(wkt_string)\nprint(f\"Polygon is valid: {polygon.is_valid}\")\ndisplay(polygon)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-28T06:48:43.541368Z","iopub.execute_input":"2022-06-28T06:48:43.542376Z","iopub.status.idle":"2022-06-28T06:48:43.559342Z","shell.execute_reply.started":"2022-06-28T06:48:43.542328Z","shell.execute_reply":"2022-06-28T06:48:43.557668Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 5. Create a Mask of the Polygon\nAside from creating new features from the polygon’s attributes, you could also create a mask from the polygon’s coordinates.\nIf you want to apply some computer vision models to it.","metadata":{}},{"cell_type":"code","source":"from PIL import Image, ImageDraw\nimport numpy as np \n\ndef create_mask_from_polygon(wkt_string):\n\n    width = 100\n    height = 100\n    img = Image.new('L', (width, height), 0)\n    \n    # Load polygon\n    polygon = shapely.wkt.loads(wkt_string)\n\n    # Draw exterior\n    ImageDraw.Draw(img).polygon(polygon.exterior.coords, \n                                fill = 1)\n    \n    # Draw interior\n    for interior in polygon.interiors:\n        ImageDraw.Draw(img).polygon(interior.coords, \n                                    fill = 0)\n\n    return np.array(img)\n\nwkt_string = \"POLYGON ((10 10, 90 10, 90 90, 10 90, 10 10), (20 20, 50 20, 50 50, 20 50, 20 20), (60 60, 80 60, 80 80, 60 80, 60 60))\"\nmask = create_mask_from_polygon(wkt_string)\nplt.imshow(mask, cmap='binary')\nplt.show()\n\ndisplay(mask)","metadata":{"execution":{"iopub.status.busy":"2022-06-28T06:49:12.392252Z","iopub.execute_input":"2022-06-28T06:49:12.392688Z","iopub.status.idle":"2022-06-28T06:49:12.52911Z","shell.execute_reply.started":"2022-06-28T06:49:12.392655Z","shell.execute_reply":"2022-06-28T06:49:12.528115Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Conclusion\n\nWKT format is a simply way to describe a polygon. With the help of the Shapely package, you can convert the WKT string to a Shapely polygon object and take advantage of its attributes. In this article, you have learned how to visualize a polygon with Matplotlib and/or Shapely. Additionally, we have discussed five ideas to create new features from the polygon:\n1. Area of a polygon\n2. Perimeter of a polygon\n3. Number of interiors of a polygon\n4. Validity of a polygon\n5. Mask of a polygon","metadata":{}},{"cell_type":"markdown","source":"# References\n\n[1] S. Gillies, “The Shapely User Manual.” shapely.readthedocs.io. [https://shapely.readthedocs.io/en/stable/manual.html](https://shapely.readthedocs.io/en/stable/manual.html) (accessed June 20, 2022)","metadata":{}}]}