{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "71b729ba",
   "metadata": {},
   "source": [
    "# Calculate coverage example\n",
    "\n",
    "This examples shows how `lusos` can calculate the coverage percentages of BGT-soilmap unit combinations in each raster cell of a 2D grid. The way the coverage is calculated is shown on the [previous page](./concept.ipynb). This tutorial shows for a small area of 5x5 km how the coverage for the BGT-soilmap combinations can be calculated. We will use sample data that is available from the lusos-repository.\n",
    "\n",
    "We begin with the necessary imports, a definition of the xmin, ymin, xmax, ymax bounding box of the example area so we can indicate the location and load the BGT and soilmap data. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "db2015fd",
   "metadata": {},
   "outputs": [],
   "source": [
    "import geopandas as gpd\n",
    "from shapely import geometry as gmt\n",
    "\n",
    "import lusos\n",
    "\n",
    "# Bounding box for the example area\n",
    "xmin, ymin, xmax, ymax = 111_000, 455_000, 116_000, 460_000\n",
    "study_area = gpd.GeoDataFrame(geometry=[gmt.box(xmin, ymin, xmax, ymax)], crs=28992)\n",
    "\n",
    "soilmap = lusos.data.sample_soilmap()\n",
    "bgt = lusos.data.sample_bgt()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "91a23a01",
   "metadata": {},
   "source": [
    "Let's first checkout the BGT data: ![BGT data](../_static/example_bgt_data.png)\n",
    "\n",
    "Now let's checkout the soilmap data: ![!Soilmap data](../_static/example_soilmap_data.png)\n",
    "\n",
    "The maps show that the BGT data has 9 different units and the soilmap has 13 different units (i.e. soilunit_codes). Lusos groups the units of the soilmap into four main groups based on the \"soilunit_code\": \"peat\", \"moerig\", \"buried\" and \"other\". The BGT and grouped soilmap therefore have 36 unique combinations and calculating the coverage of the BGT-soilmap combinations for a grid of nrows and ncolumns will thus result in a NxNx36 sized grid.\n",
    "\n",
    "As you can see, almost the complete study area is covered by both BGT and soilmap data. In builtup areas (e.g. cities) and in the locations of water bodies, soilmap data is not present. In the locations where the soilmap is missing, the combined result will produce a missing value.\n",
    "\n",
    "Let's now define the grid where we are going to calculate the coverage for. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7501112d",
   "metadata": {},
   "outputs": [],
   "source": [
    "xresolution = yresolution = 25\n",
    "\n",
    "grid = lusos.LassoGrid(xmin, ymin, xmax, ymax, xresolution, yresolution)\n",
    "grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3d3a77c7",
   "metadata": {},
   "outputs": [],
   "source": [
    "coverage = lusos.bgt_soilmap_coverage(bgt, soilmap, grid)\n",
    "coverage"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "default",
   "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.13.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}