3_homogeneous_relief_zone.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "comfortable-ground",
   "metadata": {},
   "outputs": [],
   "source": [
    "# import libraries\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import geopandas as gpd\n",
    "import matplotlib.pyplot as plt\n",
    "from rasterstats import zonal_stats\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "funny-impression",
   "metadata": {},
   "outputs": [],
   "source": [
    "# load sota basin shapefile\n",
    "sota = gpd.read_file('./data/input/sota_basin.shp')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "found-breath",
   "metadata": {},
   "outputs": [],
   "source": [
    "# load sota sub basins shapefile\n",
    "subs = gpd.read_file('./data/input/sub_basins_sota.shp')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "indie-preview",
   "metadata": {},
   "outputs": [],
   "source": [
    "# find minimum and maximum elevation of each sub basin\n",
    "subs = subs.join(pd.DataFrame(zonal_stats(\n",
    "    vectors=subs['geometry'],\n",
    "    raster='./data/input/dem_sota.tif',\n",
    "    stats=['max','min'])), how='left')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "hungarian-fraction",
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate height difference\n",
    "subs['D'] = subs['max'] - subs['min']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "talented-today",
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate perimeter and area\n",
    "subs['P'] = subs.length\n",
    "subs['A'] = subs.area"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "parental-disco",
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate shape coefficient\n",
    "subs['C'] = 0.282 * subs['P'] * (subs['A'] ** -0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "placed-imperial",
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate equivalent rectangle length\n",
    "subs ['L'] = (subs['A'] ** 0.5) * (subs['C'] / 1.128) * (1 + np.sqrt(1 - ((1.128 / subs['C']) ** 2)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "short-applicant",
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate specific height difference\n",
    "subs['Ds'] = (subs['D'] / subs['L']) * (subs['A'] ** 0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "diagnostic-demand",
   "metadata": {},
   "outputs": [],
   "source": [
    "# classify relief by specific height difference\n",
    "subs.loc[subs['Ds'] < 10, 'relief'] = 'R1'\n",
    "subs.loc[(10 < subs['Ds']) & (subs['Ds'] < 25), 'relief'] = 'R2'\n",
    "subs.loc[(25 < subs['Ds']) & (subs['Ds'] < 50), 'relief'] = 'R3'\n",
    "subs.loc[(50 < subs['Ds']) & (subs['Ds'] < 100), 'relief'] = 'R4'\n",
    "subs.loc[(100 < subs['Ds']) & (subs['Ds'] < 250), 'relief'] = 'R5'\n",
    "subs.loc[(250 < subs['Ds']) & (subs['Ds'] < 500), 'relief'] = 'R6'\n",
    "subs.loc[subs['Ds'] > 500, 'relief'] = 'R7'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "addressed-relief",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['R3', 'R4', 'R2'], dtype=object)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# check relief classes in sota\n",
    "subs['relief'].unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "immediate-tourism",
   "metadata": {},
   "outputs": [],
   "source": [
    "# create new dataframe of relif classes of sota\n",
    "relief_sota_class = subs[['PolygonId', 'relief', 'geometry']].to_crs('EPSG:4326')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "ideal-interstate",
   "metadata": {},
   "outputs": [],
   "source": [
    "# categorize by sensibility of relief classes to water\n",
    "r_conditions_ = [(relief_sota_class['relief'] == 'R2'),\n",
    "                  (relief_sota_class['relief'] == 'R3'),\n",
    "                  (relief_sota_class['relief'] == 'R4')\n",
    "                 ]\n",
    "r_choices_ = [3, 2, 1]\n",
    "relief_sota_class['W_Type'] = np.select(r_conditions_, r_choices_, default=None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "acquired-chaos",
   "metadata": {},
   "outputs": [],
   "source": [
    "# dissolve relief data with water type column\n",
    "relief_sota_class = relief_sota_class.dissolve(by='W_Type').reset_index()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "tender-poetry",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<AxesSubplot:>"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAALAAAAD4CAYAAABMvGxfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAtNklEQVR4nO2dd3hUVd7HP2eSSU9IDyEJECD0EiA0owgCFuyra9m1bkF91bWtruura33Xvrqvu6+Ka921rYsNBQVBQDoBAqETanrvyUwyM+f9IxMIqVPunTuT3M/zzJOZW879XfjOmXNP+f6ElBIdHV/FoHUAOjruoAtYx6fRBazj0+gC1vFpdAHr+DT+WgfgDLGxsXLo0KFah6GjAdu3by+XUsZ13O5TAh46dChZWVlah6GjAUKIE11t15sQOj6NLmAdn0YXsI5PowtYx6fRBazj0+gC1vFpdAHr+DS6gHV8Gl3A/YzKpkZu//Yr6pubtQ5FEXQB9yOqTU3cvfxbVhzJ5S+bN3CiulrrkNzGp4aSdbqmsqmR4vp64kPDWHb4IA0tzewuKWFoZCSjY+O4bORohBC8unkjm/JPAvBe9g4+ytnF/154MecPT9P4DlxHF7CP02y18srmjXyYs6vbY9afPEGz1crSQwc6nfvfq1cyL3U4fgbf/DHuNWohxDtCiFIhxJ52234uhNgrhLAJITK6OW+UECK73atWCHGvfd8TQoiCdvsWKnZH/YznNqzrUbwAS/bv7STeNiqamrj0k39xqKJcjfBUx5Gv3XvAhR227QF+Bqzr7iQp5UEpZbqUMh2YCjQCX7Q75JW2/VLKZU5FrXOKnJJit8s4UF7GHd9+rUA0nqfXJoSUcp0QYmiHbfsBhBCOXmcecERK2eWUOB3HkVKytSCfD/fs4nBFhWI154maavaVlTI2Ll6R8jyFp9rA1wEfd9h2lxDiJiALeEBKWdXViUKIRcAigMGDB6sapLdyy5dLyK+rwV8YsErJkapKxa9hk5L/7N/LY7FxzlRMmqN6y10IEQBcBnzWbvPrwHAgHSgCXu7ufCnlYillhpQyIy6u04T8Po/FZmNXSTFHq6o4VFmhinjbeC97B18e2K9a+WrgiUfPi4AdUsqStg1SyhIppVVKaQPeAqZ7IA6fwmqz0dDczKpjR6gxmzx23ZzSYnzJ7MYTTYjr6dB8EEIkSimL7B+vpPWhUKcdz/y0hg9zdmGx2Tx63fd27WRkTCzXjJuAwQeaEo50o30MbAJGCSHyhRC/FkJcKYTIB2YB3wohvrcfO0gIsazduSHAAuDzDsW+IITIEULsBuYC9yl0P32GhNAwj4u3jUdWr+RIpXpNFSVxpBfi+m52fdFxg5SyEFjY7nMjENPFcTc6EWO/oqqpiU/35vDK5g2axtFstWh6fUfRR+K8AKvNRmVTE7EhITy/YR3/3qdtiyrI35/o4BBNY3AUXcAaU1RXx0M/fMeGvJNcODyN744c1jok/A0GIoOCtA7DIXxzALyPYLXZ+HRvDhvyWifYeIN4AcIDAvn64H7MFu9vRugC1pB/5WTzv1s3aR1GJ4rq6/jz+nVUmzzXfecquoA1pM7cTKjRqHUYXVLXbOajPbsobajXOpQe0QWsIXdNn8kTc+ZpHUa3vLZ1M3cv/0brMHpEF7CGFNfX8c/d2VqH0SO7Soo5UF6mdRjdogtYI2xS8tLG9exWYDqkmjRbrewqLur9QI3QBawRL2z8ic8P7NM6DIf4z/69WDUaFewNXcAacf6wEVqH4BAZg5KYEJ+gdRjdog9kaERvy4C0ZmRMLLdPncZFI0YS6O+9MvHeyPo4QV4sCoMQ3D/zLJ9Yraw3ITSgxWplw0nvXV01fVCyT4gXdAFrwo6iQk7W1mgdRrfMTU3VOgSH0QXsYcwWCy9u/EnrMLrljozp3DRxstZhOIwuYA/z921b2OHF/aof7NrJnrKS3g/0EnQBe5h7ZswiwM9P6zC6paGlhZ9/9gm3f/sV+V7czGlDNWce+3HH7UuHsoUQWe22RwshVgohDtv/Rrl/K75Bs9XKgEDvn2u74kguN3+5hJxS766NVXPmacdcu/tOe6E/DKySUqYBq+yf+wUf5uyirLFB6zAc4lh1Fb9d+gVHKiu0DqVbehWwlHIdUNlh234p5UE3rns58L79/fvAFW6U5VP4ws9ye0obGnhx03rNFpj2htptYAmsEEJstzvstJHQtqze/rdbPyMhxCIhRJYQIquszHtnRTnCtsJ8Ptmbo3UYTrPiSC6Prl7plfMh1BZwppRyCq3mJncKIWY7W0BfceaRUvLq5o00W61ah+ISXx7Y75Uz51QVsH2ZPVLKUlqX4bc58JQIIRKh1eQEKFUzDm9ACMFlo8ZoHYbL/HneAiYnDtI6jE6oJmAhRKgQIrztPXA+px14vgZutr+/GfhKrTi8iWvGjmfRlAzOGzqMqKAgkiMiiAgM1DqsXjEaDFzgpUPLvc4osTvzzAFi7W48j9P6UPcaEEerM0+2lPICIcQg4B9SyoVAAvCF3enQH/hISvmdvdjngH8LIX4NnAR+ruxteSdCCB4++9wztp2oruaBFcs4XFlJXbNZo8h6JnPwEEIDArQOo0uELxm5ZWRkyKysrN4P9DH2lJbwm6VfUNrgnd1rn1x1LdOTkjWNQQixvUNXLKCPxHkFVaYmxcSbFBLGoJAwRcpqY0R0tKLlKYkuYI05Xl3F8+sdGQ/qnYzAaEJ+t4aEV3LIJIpJIdEY3Uzekp6Q6NU2U7qANeblTevZ58aqXwGkhkYQ4m8kdGcFNpOF2v0lFN+7nPpFyzirwfHaONjfn8TgUEaGRp7aNj0pyeXYPIH3LgvoJxTV1bl03vSAaIw/nESarBQv28yo8QMpye88ypf3xxXMfu581gXV9lqeYclhSlcfxubvx+xrJ2GICWZERDO2s2wYvDQNl/4QpyFHqypZ+OEHNNscH9zIJArbt7mU/HAY4cR/XcLCMTRekcoxSwP+BgMjrMEENNqwGQUNb2yn7kD3XfELbjqXaRdOZu51mY5fUGG6e4jTa2ANOVZd5ZR4AUz/yqE6Kx9nvdNLlu2H7w4QZ69IK4QB2eLYtVd+sJY1n2wg88rpBAR6lxWWd/4u9BPm2gc0HGVcaBRVWfmuX9AmkRb7y0HxttHSbKE83/tmpekC1pCDFeVUOeEAGbujyumaV0kMft4nF++LqB/xHyec2GODgil8V7v2f2BwANGJ3rfuQBewhmy0G1v3RnJIOAlPbsdm0s5wutnUQkWh9yV+0QWsEXVmM8equ0xO2olhVQZMRT13g6mNlJKcdd6XBFEXsEYcra7qdW5wuDGAzDV1lDzzo4ei6pnPXv6amnJtv0gd0QWsEWaLhbBeZnhNPWKj+Mu9WOqaPRRVzxzfk8f95/6JZnOL1qGcQhewRkxPSmZmUjIGIQjy69wdf7YtkvLPvC+B6cn9BTx6ybNsWpqF1QtWl+gC1oBmUzPlhZUEPZ/F5LePM31r0xlpXTNFFCWP/EBTXrV2QfbAzlU5/Ony5ynL075fWB+J8zBl+RU8eumzHN112tyvencRU89ayDZz61N+xRNrNO1xcITkkYnEJXdKwupx9BrYw3z6/JdniPcU7+ZwbkM452w101LV5PnAnCRqYCR+/to7DKnmzCOESBFC/CiE2G8/9p52+54QQhTYHXuyhRALuyqjr/HNmytZ/vaqLvdVbs0j/79XUvjRTg9H5RqT507QOgRAXWceC/CAlHIMMJPWZfVj2+1/xe7Yky6lXNZ1EX2Lwtwimk3e8wTvKimjBnHBrXO0DgNQ0ZlHSlkkpdxhf18H7Ae8e3a0ysy7YTb+Ru1/dt1h1mUZPL30YeIHe4dHh0ce4oQQQ4HJwJZ2m+8SQtwEZNFaU3c5LGV39FkEMHjwYJUjVQ+bzcZf71iMxclZYN5CyugkZl81k5ufuhYhtJxSdCaqP8QJIcKAJcC9Usq2YZzXgeFAOlAEvNzd+X3FmSf/UBH7N3tHMm9XuO3FG7nl6eu8SrygsoCFEEZaxfuhlPLztu1SyhIppVVKaQPe4rRjT5+kpryWjV9u1ToMtwiPVnals1Ko1oQQrV/Vt4H9Usq/dNiX2GbuB1zJaceePofVauW+cx4j72Ch1qG4jBCCrct3MnbWKK1D6YSazjyZwI1AjhAi217cI/YehxeEEOm0ulceB25T8qa8ieN78sg/5L0pBXoiPCqUxz57gJTRSVQUOjZzztP0KmAp5fXd7Pqii2MLgYX29+uh6wUEUsobnYjRp1n21g/40sLZ9tRVNTAwNZ7YQdHEDvJOcxN9JE5lZlw81esefBzF4GfA1OCdfm1t6AJWmekXTWZkxjCtw3CJUdOGk5SWqHUYPaIL2APED1Gv+y8yI5mYs4ZgCFL+efzGx6/xumX0HdFno6nMvk0HCQpVzwM4+JcTWC+qiL1pGDYpifIPJOVgE5ZBITT9M4fqnQXEzE4lZO4w8p7seh5GRyJiwrnwV+cxPtP7eh06ogtYRdZ+tolnrv1L7we6iCHQn71BTWCGclPrDLZKs4kjyQDN+N0ymMH/NY5DlmbMVjOp8aGYS3t3wRw2aQi/ff4G1eJWEl3AKmBqNHN09wn+dvfbql5n0BPncbS5jpn+0QTtrcRS1ogxNRLbgADyE/w40lDDsYbWwc85x/zIr2h0qFxvmCbpKLqAVeCT577gw2eWqHqNxL9cyNaAekKEkcYXN1F2vHM/bfroeKIuG4MlPpi8v37XRSldM3XBJCVDVRVdwL3w1y0b+d30WQ51hR3afoS//PYNju3uYsK6Qhgjgki4ZgKbjHU0W1onBkXNH0HDP7Z1OrbuQGmPpn0diYwfwFNf/YHh6UOVCld1dAF3QEpJjdlEtcnE2zu382HOLg5VVHDfzLMYEd3zEppXb1/MkezjqsWWcu/ZFI0LY21DNbRbULkvI5zEH6Jo6KIWdoah41MYPX2ET/Vb6wLuwDeHD/Lalk2EBway055VfnnuIY5XV/GbyRkMCAqixmTicGUFF6WNZEJ8Ak31Tbx46985vP2oanFFZiSzLtWKraG6074yUyNj5wyn4T33rKeyV+9h/RdbOednM9wqx5PoAu7Ad7mHya3qbKG0v7yMB1YuP2Pbx3t288ezZ1O19hg/LdnS6RwlCbphAraup0wDsGVqEDMCZpC/2L04asq8y7ikN/SBjA7EBAc7fGyN2cTDq1awoVndFLiRExPZaug5x7LJaqFgWiQpj53n1rV2rvatVLi6gDtw1EG/svbUv7FdhUhOE3bteFocyFOc21BDdYrjfsNdMWWedyzWdBS9CdGOZ9at4Uilcw6McUEh1B1SvgYOHhRB3PWTqE0LZ63J8ZjcXbA0LnO0myV4Fr0GtvPTyeO8k72dkoZ6p84bV6BOPBEZyfyY0sJ2J8QLkNNYSdPi+aQ8NR/h5ALStCmpJKUNdOocrdEFbCewC3+y3ogODKL8HXWaD8VL9xHvYn62wsZ61kTUE3tOqlPnzb56FsYA75680xFdwHZO1FQ7fU6l2UTBo5NJeWq+4vHEzByCyQ3zvGB/f6q2nQSD4326S179lqKjJS5fUwtUc+axH3ehEOKgECJXCPFwu+3RQoiVQojD9r+ae9fnVrpmVFdlNrEhspHBD88h5fF5BA+KUCSe0HnDqHUj+XeTxULYi/NpemMeQYmOxVRdWsMDcx7njfvfw2rxjeX/qjnzCCH8gL8DFwFjgevbOfM8DKySUqYBq+yfNcXZtm97Wmw2fhxoYk1UA7GLpikST/m7O9xOE7vdVElhYz0Rj5zt8Dll+RUsefVbaitcS8DoaVRz5qF1qXyulPKolLIZ+AS43L7vcuB9+/v3gSucCVoNpgwcpEg5DQmBRIxNQPi7NxwbmTnUoa4zR9hkrSR+QZrDxxsD/LHZfGMdn5pt4CQgr93nfE5bSyW0Lau3/41XMQ6H2FOmTNtvm7mSnYtSka+fT2BCuMvllHy9l9kmZZojAEXXpJL0y8n4hfbsCg9w1X2XEOOFGYm6Qk0Bd1UFOf21FkIsEkJkCSGyysrUG/FK62WijrMcaaih5ukZhI927btpM1nIf3gFU4OUWQ18vKGWn6YF0vTKuYQOjUIY/TAEdt3zMn3hFEWu6QnUHMjIB1LafU4G2tw9StrMTYQQiUC3c/6klIuBxdCaK1mNQGtMJpbs26t4uQWN9dT9bgzjv42h5FvnM/wIoObu7xn1fws42MUkHmeRQF5jHXn3jiLUaCQ5IBR5x/cMvGIC/kMjIcQfJLxWdYDFjHH7ep5ATQFvA9KEEKlAAXAd8Av7vq+Bm4Hn7H+/UjGOXqloauSQi70QvVHbbGbbBVGMuuIiGrAyqMGPyqfX0uzg6gjZYiXyX4c457KR5IZZKGpy/WGzPQ0tLRxsqWbk6xewvqEaaLa/IF6h5pQncKQb7WNgEzBKCJEvhPi1EOJKu0vPLFqdeb63HztICLEMQEppAe4CvqfVWvXfUsq2au45YIEQ4jCwwP5ZM3JK1f0Pa7HZ2NNQxbGGWjZQRdn/zHCoLdpG+dqjFD7wHcOalB9kONRFzX7FKN+ofUFFZx7752VAJ/NqKWUFMM/xMNVlWJRnXWfKTU0Mf2E+ZXcucyr3sV+LDVQeKBsZE8ulI31nPoQ+EgdMiE/grBTPeg9vaakk5bfOmXIWPbrK5eFlR1k4YiTj4hNUvYaS6AK2838LLyUlYoBHr1ny6W6njjdGBpMi3Zsu2R2hRiNn26I48dJaVcpXC13AdiICg1h7y29YNKXbkXFFEUBAjHO1qbmkjsCv1Fm2NEGGU3T/cnYv2c7iBz+g3AsTe3eFLuAOLJo6jaRw5QYQuiMqMIiGY86LxBBiJNhfmc6jtgyhMUHBiK9zT23/7OWl/HLIHSx55RtFrqMmuoA7EB0cwusXX8agcNdH0RyhRdoIjA116pzEK8ezfe4AmizKJEGcme/HiAe3EXP/ekpXHDpjn81qo6nepMh11EQXcBdEBwdT1aRessHowCCmbG3EVNG7zRO0DkAMui6d7PnRDos3IiCQOdWhDAjo7MvmJwSZRFHyj23IFivWxq5TfwUGO97VpxW6gLsgPjSMGyamq1b++K315L+XhXBwXHHwbTP4aWYQ9S2OZ62fsreZvCdWMWZJIXPyA0gKCTuVj/msLSaK713e62DK8MnOTYjXAl3AXXCypprsYuXTAoQZA5i9s4WCdx33bzAE+bNtonM14ZziQIo/zAagbNVh8l5aR/B/reKc+jDig0MoW37AoXJO7sunxmSi2uS9qW/1RZ0dKG9s5NKP/6lYO7M9U7c3kf++c0uQjGGBJASEcLSl52X1AEaDgczyIPKeW91pn7RIqt7eQVKQP7UO5mKOiA3nvhXLKKit5fWLL/P4gI8j6DVwB2JDQrhvZqby5QYFU73umFPnBMaH0vRsJkcbehcvwMjgAeQ91Vm8bTQcraB2n2PD5lLA3yzHWHP8GIcrK7jxi884UV19ar/VZqO+2fEmjVroNXAXFNYp604jgBGf5VOWW+7wOWFpcVgeyuBYg2M+Fee0RFL+sHKDEAPGJLCz6vQkwaL6eub98x1uTZ/C2Nh4Hlm9kiB/f+6flcmNKj4v9IYu4A4U1NayLPdQ7wc6gRCC2hzn8sQF3zONzQ2O9RMnBodS/PAP2EzKNXvEPVOh6czr26Tk7Z2nm0Bmq4XH16xif1kpT82dj7+bS6BcQRewnTqzmec3/kRWYQGlDY51bznC7JZIDJUm8sodmz4JkHj1RLbYqh0+PtDgh9VkcWpiUE/EL0gjy+z49T/Zm8PGvJOMi0/ghfkXEBrgue43vQ0MmC2tNclHObs4VOH4z3xvRAYGUfr0GvKeXQNOrDFrPFiG0eC4KcnxhloSLlAun4XNYnN6Pd7J2hqW5x5i7YnjisXhCP1ewGaLheuWfMqXB51fMdEbSf7BtNQ6P5pVk1PE9Fzp8Krk+OAQ6vcrN6fZ1uD6w9ldy5dy//fL2JR30iMJHvu1gK02Gy9s/IldJcWqlB9T4Ya3Qoi/w7XgmC11NJ6sdv1aHQid7d4AxpcH9/PLLz5jQ95JhSLqnn4t4K8PHuDd7B2qlB1uDKD46TUun2+LdCw116CQMPLf6ZxewC1alDE1OVKlzjKt9vTrh7isIpWc+YCMo5KTbvwU97Z+e2JINMEfH6Qh9wBKj5M1Dx8ANvenU8aFODdZyRVctZbq1RpKCDFKCJHd7lUrhLjXvu8JIURBu30LO56vJjYpOVBeRp1ZnTzAEQGBlPwr2+Xz484bwcbgnp1x8qyN1O4uoimv2uXrtGfAhERS/nAuiVdPZHMPTvDOEGxU3yjQVWupXq2hpJQHpZTpUsp0YCrQyJnr6F5p229fO+cRrDYbme8sZuFHH/DN4d7MhZxnWmA0yX/KwlziujVTdVY+qcE9T+esMpuInJLU4zGOYowMpuTucaxJNLP+7BBsCj18Oeu17AouWUvhvDXUPOCIlFK9/FMOsvbEcbd80Hoj7FgD5nL3+pFbak0kH+n91yFwkjKJuBMXTaegUdl/E6PBwLlD1J/N5upDnLPWUNcBH3fYdpcQYre9idKtj5HSzjzHXUgh4AzWMmUGQaQD3mp7R/gTNS2l1+N6IvHqiaxLUn5OQ4vNxs/+/SF7VLYsUL0XQggRAFwGfNZu8+vAcCAdKAJe7u58KeViKWWGlDIjLs79rO82KRUbseqK4k93uV2GX4iR/Sm9R1lmauTwrcNIfmYByb+eRtJz5xP41kUk/9oxh8ywtDis5yZjVam/tqGlhXUnjmO2WCiqq1Nl8o+rvRAOW0PRaq+6Q0p56qvY/r0Q4i3AY4uvbk2fwmtbN1PnhvduT0TPHELZ6tzeD+yBgddO4lCTY0PPNc1m1oaZYYIfUAsNED41mKnJ86DWDM1W8v6+CQmEp8UhBET8ajLCKtkV1Ux1s3rt1ACDH7tLinlgxXKW5R5i8sBEPrnqWox+yuVidlXAzlhDXU+H5kOb+O0fr6TVb9gjFNbVYbJ0vYTGXaYHRFO7x700VSlPzmfNAPfao3XNzayJaoao1iQ0yW8uxN/W6pzpbzBgsdln26k8G7LZZmXF0dNf5p3FRdzy1eeMjInhwbPOIUSBXgqXrKXoxhqqvbWU/XOIff/nHYp9QQiRI4TYDcwF7nP7ThzAbLHw8ub1ivnutic2KJjmV7dhLnWvDSz9nG3gSPxE9/dTZmpkZ1Ml28ytNa1FhXt3hk35J3l/105OKPQs4o61VCdrqC6spRqBTr6lUsobnYhREZ5dv5alBw9QrFIPRKIhWJF0W0UJftDDd8BokFw9wsLoqNYmUEZcGSarkauWJ9G1o6138k72Dp6ZO59ANy0C+s1InEEI1cQLEFlrw/EJk10TmBDOMXPPMV4+zMoTk7/GT56eNSfxJ9jvVzT5RloLAJbs30ut2cT5w9NYe+IYD2fOZpALfhz9Zi7Ez0aPU63scGMAJY93v5THUUKSInpdi/efXH8+Pn56XEli5GjxHIL8HU+R6y2sPHqEB1d+xzeHDrLgn++6NC+l39TAe8t66ihxjyHGMJrcmfdgJzDJMW+257NCmFg7nxM5EBoOz99WxqTbmlgzxneaEB1psljIcWFWYL8QsNliYekhx5aSu0LcwTqUmDgogh17Km+yCp59LIjKradTkOS9uYU5T86jIT7w1AObrzEhwfksof2iCdFstTIr2b0Rq56oSVPGhqpkgWNDwwJoONpZpHmPr6LyjmXMLfB+R52OjI9P4Ooxzjfz+oWAwwMDuWHiJBLDwlQpv0Whf8Uwm2MFSSD+hsnd7m/YknfKuM9XuGliOuGBjs2Bbk+/EDBAkL+RD664mtTIKEU60NsTc9T9+Q/GqGCKcXx0cF2qtdusoBXrjjHrpO/8116cNpKrXKh9oR8JGGB4dAwrb7yVbb+5g/FxyqWmqx8S6rTXb0fCR8bRZHO8H8wmJXG/TO92/4mX1jEqNNKtmDzFpIREhHDtAbRfCRha+4ODjUaWXPMLzh8+QpEys8yV1D+bScR45x5CghIjiMxIJvzNhWy7fhDVZucWgFYt7f7BVEgIX7yXOSf8ObdRXatYdxDAhSMczyLaEd9qKCmI0c+PACeWrvdGXmMdJbcNJy1oMv7Pb6O+BxeegJgQ6p/N5ISpAQG0NLnWaxB5QRp1B7rvHqzalkfVtjwMQf6cfdk4LLOT2Gzxrh6K1Kgokt1I7dDvauD2xJRbGBRy5oOdOwm2m61W9jZUwUM9J28xVzRSb23BYnPef6E9DaMd+4+3mSwU/XsXZXctI3NVLekh3mHSNyE+gTcvvrz3A3ugXwvY+N1xQu78kYwP8piTH8CcyhCmfllKiL9rD3kGITjHGknAuz1n/TQY/RiowMjZNnMlcfOc+/ktXrqPgPeUz0rqLJMSBvL4uecx3M0Uv/22CbF73T52rtqDbLFStaOAqh2nVyhPviyNDTg/W2p8cBSFi3pf3pf8h9mscdC0ryeC/f1pLnV+fkd1diGTbx7HThebLkrw2ynTmJI4yO1y+mUN3Gxq5uXfvE51ade2pcX3f8ecg5ApohxuUszJD6DeAfECTtlM9cQoYwTVOc4bcVsbW7A+v4WILtIPeIKpiYOYlzpMkbL6pYARgsLcHsbdbZK81zdTfM9yMg/ITv/RA4NDGR8STUZgNDO+KmPaRwXkvbTOocmMwYMiODhEmR++7KZKBj92Xm8WEl3ScLSCyRXqL3vviu1FhYoZyvTLJkRjreMTH/P+bxMJgf5Mve9sNiS1ML0lnMLbvqMJXDIUMQT6U2pSzv1yTUwjM88b4dIypoInVzH7+ik0npVALg1Uq+ST0RE/IZilUGbUfingwzucc0q3mS2cfG4NEyYmUrS7yL1p43+cibVRubZnTFAwZWuOuHSutEgK/rkd/gmj/7aQzU6MBLpDsNGo2EBSrwIWQrwDXAKUSinH27dFA58CQ4HjwDVSdrZzEUIcB+oAK2CRUmY4c75a7F67z6Xzqne7l/gl+ZYMcqzKGkFJCcIA0s2VQo5mTHKXB2ZlMn+YMgNIoKIzTzvm2t132udwdeZ8xYkeGElAkOfbf5Wrc2mRyi6bmHjchrS4r77AKvXzXSSGhXHVmHGMionFTyE3d08583TE3fPd4srfLWTGxVM8eUkAGk9WE+qn3Bcn0M+P8g92KlOYVf3FnguGjWBgmLLD2mo780hghRBiuxBikQvnK+7M04apwTPtvY7ECeW6rib4D6CpSJmENCefWs1IlSf/nDtUeasptbvRMqWUU2g1N7lTCDHb2QKUduZpI2ed8o7sjhCioCVFWKFyuYyFhFAH5yO7QkJoGCOjYxUv19WIS+yOPPTkzGNfZo+UspRWZ8rpzpyvFuUFFZgaPV8DB8aGssNarVyBCltC+an4IHfhiDSSIpxfddwbrgq4zZkHunHmEUKECiHC294D53PagafX89XCZrPx0Z+/6P1AFYi/dqJixiKDQ8KxHatWpCxo9WPb06Jsfrw2Avz8eGz2XFXKVtOZJwFYL4TYBWwFvpVSfmff1+X5amC1nn7qP3mggBuH3cnS179X63LdIgXsHh/Y46hZekg0s5sHMOeIgejAoG6PC/LzJ+jeNRR8qNADHK3DyxM/K2BoqPK15HlDh51KNK40qjnzSCmPApO6KbOiq/PV4E+XP0/+wdYkg9WltTTWaZO4OmbmEIqsXXs+pISEM7TWj8rH11NgfygbMn4g3DmGSnPneP0NBmzNyudyrlh3jJTkSI5PV+ZBMyY4mIcyZ3PO4CGKlNcVfXokbtt3O9m2PNsj6Z56w3j9OJraTSYPDwhgcl0whrImTr60ko55PGv3FDP0TQOTFoygZGQoFbZmSppah6D9DQa3By66Q0yKBxzLzdwbLy64iDkq9Dy0p08L+K93vOUV4gUIKDMxISWaMmlmdDHUfLqHAnvi7e5+XGuyC6nJbpX24BmD4VcjiDAYST7QSHFUMM0V7ppZnUlQYgTbwxpBgZ6S4VHRTE9Kdr+gXujTAq6vVm7SjLvkPd1qPRUe6E+e2fmf/8otJ4k6VEZzVZMiJipdEXPfLPa2uF/7TkwYyOJLLld89XdX9NnplFJKBqYqt/JYKWwuiLeNlqomVf0n/WqUGU6+ZdIU4kPV8eDoSJ8UcF1VPQ/Nf5KjuzTPKeNTNG5Upm4/UeOxeVl9U8CNtU1k/7jXa9q/vkLwXPcfuIwGg0s2qa7S5wR8dPcJfj/3ca3D8ElaIlxrs46OjWPIgEgEcFHaSH4+dryygfVAn3uIMzc1U1lcrXUYPomx3AROTlcwCMHTc+cxeEAkLVbPO2z3OQHHpcTQ4saDUn8leFAEe5KBLuYHBfj58av0qQwMC2N4dDRWm+SR1SsorKvj97POZmqiMhlDXaHPCTgyLoLJ88az4wf3sgX1N5qK6xjgF0B5Fyv9bpiQzkOZ55yx7cebfk1BXS1DI7vNUekR+lwb2N/oT/rcCVqH4XvYJMaH1hEbdNpwJTYkhKvGjOPemWd1Otzo56e5eKEP1sBNDSY++vMSrcPwSVqqmvATrXXa3KHD+J/z5iu+gkJp+lwNHBQSSOYVPXuT6XTPyJZgnpk7n7cvu9LrxQt9UMBCCBqc8H3QOZO0nyr4xYQuJxF6JX1OwKDdejdfYWBqPA+9fxdPfvEQySNP5+UwBho5/yZ1Jp6rRZ9rA9tstlPzf/szd7xyC4nDEljz6QZWf7T+1PbBY5L447/uYcTk1lE3P38Dm7/ZzrQLJ2NqNDP9ou5zb3gjfU7ADTWNlBd4l4mzpzn7ZzO44u6LMBgMzLxkKuffPAdjoJGAICNhUWEkp52udWdcPJUZF0/VMFr3cGRJ0TtCiFIhxJ5226KFECuFEIftfzv1pwghUoQQPwoh9gsh9goh7mm37wkhRIEQItv+WtjxfFcJHRDCgpvPJTLedddvX6emrPZUzgkhBFMXTGLi7LGMnp52hnj7Amo681iAB6SUY4CZtC6rH9tu/yt2x550KaWDvqS9YzAYeOAfdzDvF2crVaTPUXC4qN9MZFLNmUdKWSSl3GF/XwfsBzwy5mhqMPPtWz944lJeyZiZaRgUsm7ydtR25gFACDEUmAxsabf5LiHEbnsTpdshHVeceQ5vP9qveyKmzPedbjB3Uf1rKoQIA5YA90op24wHXgeGA+lAEfByd+e74syz4oM17oTs08y74RwuuX2B1mF4DFd7IUqEEIlSyqKenHWEEEZaxfuhlPLztu1SypJ2x7wFfONiHJ2oLK5izccblCrOZwiPDuPxJb9n0rmuZbz0VdR05hHA28B+KeVfOuxr/yh8Jacde9xm23fZtKjgmeDt3P7yzf1OvKCuM08mcCNwXhfdZS8IIXKEELuBucB97t5Is7mFbxevZOUHa90tyueIHxxL5pX9c/6Hms486+nG8kBKeaMTMTrE+s+38Orti5Uu1id45KN7CY1wL1ezr9Jn+lp2ruqfE9h/+ehVjDtrlNZhaEafGEq2Wq3s+GG31mF4hNikaP70n9+z+sOfaKxv4panrtM6JE3pEwIuy6ug9GT3ybX7Eg01jYycOowxM1zP8N6X6BNNiPBoz7jAeAMvrn4CP38/rcPwGvqEgINCtEmZqgW1FXVah+BV9AkBIyAqoW/PPoseGMkjH93L2FkjtQ7Fq+gTAvbz8+Osy6ZpHYaqpE4cwtzrMvttd1l39AkBAxza7lq6VV8gKDSQq+69WOswvJI+0QsBEBHrOUM5TyGE4Kmv/sCE2WP0mrcb+kwNXFnkOUtPTyGlZPSMEbp4e6BPCDg3+xjH9+RpHYbbCCEIiwwFIDwqlHvfWERkXN9+OHUXn29CNJtbeO3Of/SJJTQD4iL417G/s+HLbUTGRTBl/kStQ/J6fF7ABYcKqSo5ndfB4GfA5oHE1WowaHgCgcGBnHd9/13P5yw+L+DUCUP4IPdv1FbWYW5sZkBcBNUl1az/fCu52cfI+Wk/xcc8msnWZSbP000JncXnBdxGRHQ4RLe+jx8cx8/s3U652cd48da/e22+jLY27/D0oUxfOEXjaHyPPiPg7hiRnsqbO19ixw+7efTS52gxK5guXgGe+/5RhoxL6VfD4UrSJ3ohHGHK/Il8WriYe15fxLBJ6qU+dYYFN53LsElDdPG6gWrOPPbjLhRCHBRC5AohHnb2fKUJjwrjktsW8ObOl3hx1eMMHZ+CMcAfg59nv8dnXT6NPy//bx58906MAeonA+zLqObMI4TwA/4OXASMBa5v58zjiLOPqqTPHc9rm5/lk4LF/M+3jzD1fM95KcQmRTPtgvRT9k86rqOaMw8wHciVUh6VUjYDn9jPc/R81QkKCSQiJpyM8yfx1Fd/YMhY9XP7AuQfLvLIdfoDajrzJAHth8fyOW0t5bCzjyvOPK4QEGjkrxueUXWyePTASJ766g/816u3qnaN/oaavRBd/T46PVwmpVwMLAbIyMhQdbgtIDgANX/Vr3nwcmZdmqHeBfohajrz5AMp7T4nA23O0w45+3ia/IOFWFqUTdY3cGgc1zx4Ofs2H+KKuy9StGwdFZ15gG1AmhAiVQgRAFxnP8/R8z1ORVG14mXe9OS1XHrHBfzh/bv1tWwqoJozj5TSAtwFfE+rteq/pZR77cV2eb7WjJw6jIwLJrnUOxCXEnNGf27alFRSRieRmNqjcaeOmwhfmsWVkZEhs7KyVL9ObvYxdqzczdI3VvQ4j8LgZ+CGx67G1GDm2ocuJ/9QIb+f+wQtzRZe3/7CqQGT/uLVqyZCiO1Syk4PEH1+KNkVRqSnMiI9lZITZaz8YC1N9Z0TCBsDjTz5xYNMu/B0UpSxs0bxWenbrP10I7HJ0bpwPYBeAzvAvs2H+P6d1WxamkXyqEFcctv5DJs4hKHjUno/WUcR9BrYDcbOHMnYmSPdt9DUURz9N07Hp9EFrOPT6ALW8Wl0Aev4NLqAdXwaXcA6Po0uYB2fRhewjk+jC1jHp/GpoWQhRBngnQYP3RML9JUEHlreyxApZadcwz4lYF9ECJHV1Ri+L+KN96I3IXR8Gl3AOj6NLmD16Uv5b73uXvQ2sI5Po9fAOj6NLmAdn0YXsAIIIYKEEFuFELuEEHuFEE92cYwQQvyv3ehwtxDCK82AHbyXOUKIGiFEtv31Jy1iBX1JkVKYgfOklPVCCCOwXgixXEq5ud0xFwFp9tcM4HX7X2/DkXsB+ElKeYkG8Z2BXgMrgGyl3v7RaH91fDq+HPjAfuxmINLuSuRVOHgvXoMuYIUQQvgJIbJptclaKaXc0uGQnswOvQoH7gVglr2ZsVwIMc6zEZ5GF7BCSCmtUsp0Wj3gpgshxnc4RBGzQ0/gwL3soHVuwiTgNeBLz0Z4Gl3ACiOlrAbW0NkUvCezQ6+ku3uRUta2NTOklMsAoxAi1uMBogtYEYQQcUKISPv7YGA+cKDDYV8DN9l7I2YCNW0eyd6EI/cihBgo7AZyQojptOqowsOhAnovhFIkAu/b0yoYaDUy/EYIcTuAlPINYBmwEMgFGgFvdbl25F6uBu4QQliAJuA6qdGQrj6UrOPT6E0IHZ9GF7COT6MLWMen0QWs49PoAtbxaXQB6/g0uoB1fJr/B4GCxNy7QlH/AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot categorized relief\n",
    "relief_sota_class.plot(column='W_Type')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "respected-resource",
   "metadata": {},
   "outputs": [],
   "source": [
    "# export to homogeneous relief zones of sota\n",
    "#relief_sota_class.to_file('./data/output/homogeneous_relief_zone.shp')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.9.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}