4a - Sentinel-1 GRD Batch - Subset.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<font face=\"Calibri\" size=\"2\"> <i>Open SAR Toolkit - Tutorial 4a, version 1.2, June 2020. Andreas Vollrath, ESA/ESRIN phi-lab</i>\n",
    "</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![title](auxiliary/header_image.PNG)\n",
    "\n",
    "--------\n",
    "\n",
    "# OST Tutorial IV-A\n",
    "## How to create near-daily timeseries over Vienna. Introduction to GRD Batch Processing part I.\n",
    "\n",
    "[![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/ESA-PhiLab/OST_Notebooks/blob/master/4a%20-%20Sentinel-1%20GRD%20Batch%20-%20Subset.ipynb)\n",
    "\n",
    "--------\n",
    "\n",
    "**Short description**\n",
    "\n",
    "This notebook provides an introduction to the batch processing of Sentinel-1 GRD data using OST's *Sentinel1Batch* class. This is a subclass of the *Sentinel1* class, and thus inherits all the functionalities of the *Generic* and the *Sentinel1* classes for the generation of a project as well as data search and refinement as presented in the OST Tutorial II notebook. The *Sentinel1Batch* class holds functions for the batch processing of single calibrated backscatter products. Furthermore, time-series processing and the generation of multi-temporal statistics, referred to as timescans, are introduced. \n",
    "\n",
    "Within the given example, time-series for 4 different overlapping tracks are going to be produced over the city of Vienna, Austria. The notebook demonstrates:\n",
    "\n",
    "1. the reduction of data processing by automatically subsetting the data,\n",
    "2. time-series processing and the corresponding ARD types,\n",
    "3. merging the track specific time-series into a unique time-series with almost daily coverage,  \n",
    "4. creation of a timeseries animation for outreach purposes.\n",
    "\n",
    "--------\n",
    "\n",
    "**Requirements**\n",
    "\n",
    "- a PC/Mac with at least 16GB of RAM\n",
    "- about 25 GB of free disk space\n",
    "- a Copernicus Open Data Hub user account, ideally valid for at least 7 days (https://scihub.copernicus.eu)\n",
    "--------\n",
    "\n",
    "**NOTE:** all cells that have an * after its number can be executed without changing any code. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 0\\* - Install OST and dependencies \n",
    "\n",
    "**NOTE:** Applies only if you haven't fully installed OST yet, e.g. on Google Colab, "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# !apt-get -y install wget\n",
    "# !wget https://raw.githubusercontent.com/ESA-PhiLab/OST_Notebooks/master/install_ost.sh\n",
    "# !bash install_ost.sh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1* - Import of Libraries\n",
    "\n",
    "In this step we import some standard python libraries for OS independent path handling as well as the *Sentinel1_GRDBatch* class thta handles the full workflow from search, download and processing of multiple GRD scenes. In addition, the OST helper module *vector* is loaded to create an AOI based on Point coordinates, and the *raster* module for creating a time-series animation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this is for the path handling and especially useful if you are on Windows\n",
    "from pathlib import Path\n",
    "from pprint import pprint\n",
    "\n",
    "\n",
    "# this is the s1Project class, that basically handles all the workflow from beginning to the end\n",
    "from ost import Sentinel1Batch\n",
    "from ost.helpers import vector, raster"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2* - Set up the project \n",
    "\n",
    "Here you going to initialize the *Sentinel1Batch* class by determining the project folder, the AOI and the start and end date. In addition we determine the image product type (i.e. GRD) and the ARD type that we use for processing. In this cas we choose the OST-RTC that produces Radiometrically Terrain Corrected products, i.e. the images are corrected for radiometric distortions along mountainous slopes. This type of ARD is advised when doing land cover and land use studiesover rugged terrain."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# define the project directory\n",
    "project_dir = Path.home().joinpath('OST_Tutorials/Tutorial_4a')\n",
    "\n",
    "# define aoi with a 2 point coordinates and create a buffer with 20km radius\n",
    "lat, lon = '48.25', '16.4' # Vienna\n",
    "aoi = vector.latlon_to_wkt(lat, lon, buffer_meter=17500, envelope=True)\n",
    "\n",
    "# define the start and end date\n",
    "start = '2020-05-01'\n",
    "end = '2020-05-31'\n",
    "\n",
    "# initialize the class to s1_grd instance\n",
    "s1_grd = Sentinel1Batch(\n",
    "    project_dir=project_dir,\n",
    "    aoi=aoi,\n",
    "    start=start,\n",
    "    end=end, \n",
    "    product_type='GRD', \n",
    "    ard_type='OST-RTC'\n",
    ")\n",
    "\n",
    "# do the search\n",
    "if not s1_grd.inventory_file:\n",
    "    s1_grd.search()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3* - Plot refined data inventory\n",
    "\n",
    "The resultant dataframe from the search inventory is visualised. We do not do a refinement step here, since all images are fully overlapping the AOI. This allows us to create a combined, almost daily time-series of all images."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#---------------------------------------------------\n",
    "# for plotting purposes we use this iPython magic\n",
    "%matplotlib inline\n",
    "%pylab inline\n",
    "pylab.rcParams['figure.figsize'] = (19, 19)\n",
    "#---------------------------------------------------\n",
    "\n",
    "# plot the inventory\n",
    "s1_grd.plot_inventory(s1_grd.inventory, transparency=.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4* - Download of GRD scenes\n",
    "\n",
    "As already shown in Tutorial II, you will download the scenes based on the refined inventory dataframe for the respective produckt key."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s1_grd.download(s1_grd.inventory, concurrent=10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5* - Set ARD parameters\n",
    "\n",
    "Similar to the *Sentinel1Scene* class (Tutorial I and III), the *Sentinel1Batch* class handles the defintion of ARD types in a hierarchical dictionary structure. You can use the same types and steps to customize as for the *Sentinel1Scene* class. For our example, we already intialised the class instance with the OST-RTC ARD type, in order to remove the backscatter distortion on mountainous slopes. This applies to all single image processing in the first step. The subsequent time-series processing will bring all the imagery to a common grid and apply a multitemporal speckle filter, that is much more efficient than speckle filtering applied on a single image. Since speckle filters a conceptionalised to work on the raw power data of SAR imagery, the conversion to the dB scale is only applied after the multi-temporal speckle filtering. In addition, all images are clipped to the very same extent, that is defined by the common data coverage of all images per track as well as the AOI.\n",
    "\n",
    "Note that it is possible to change the datatype of the output to either unsigned 16 or 8-bit integer. The backscatter data is therefore linearly stretched between -30 to 5 dB. This has the advantage to reduce the necessary storage by a factor of 2 for 16-bit uint and a factor of 4 for 8-bit uint data.\n",
    "\n",
    "The exact processing steps are as follows and depend on the ARD type:\n",
    "\n",
    "|               ARD Types          | OST-GTC      | OST-RTC         | CEOS           | Earth-Engine |\n",
    "| :-------------                   | :----------: | :-----------:   | :----------:   | -----------: |\n",
    "| **Time-series processing steps**                                                          |        \n",
    "|  Create stack                    |       x      |       x         |       x        |      x       |\n",
    "|  Multi-temporal speckle filter   |       x      |       x         |       -        |      -       |\n",
    "|  dB conversion                   |       x      |       x         |       x        |      -       |\n",
    "|  Layover/Shadow mask creation    |       x      |       x         |       x        |      x       |\n",
    "|  Datatype conversion             |       -      |       -         |       -        |      -       |\n",
    "|  Clip to common extent           |       x      |       x         |       x        |      x       |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('-----------------------------------------------------------------------------------------------------------')\n",
    "print('Time-series processing parameters hold in the configuration file:')\n",
    "print('-----------------------------------------------------------------------------------------------------------')\n",
    "print('')\n",
    "pprint(s1_grd.ard_parameters['time-series_ARD'])\n",
    "\n",
    "# custimze some single scene ARD parameters\n",
    "s1_grd.ard_parameters['single_ARD']['resolution'] = 50  # reduce for processing time and disk space\n",
    "s1_grd.ard_parameters['time-series_ARD']['dtype_output'] = 'uint8'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 6* - Run the batch routine\n",
    "\n",
    "To process all the data, including time-series and timescans is as easy as one command. All the complexity is handled by OST in the back, and you just have to wait, since processing can take quite a while. Note the keywords to aly higher level tim-series and timescan generation. Mosaicking refers toacross track mosaicking, which for this example is not the case. The *overwrite* argument tells OST if it should start processing from scratch (i.e. set to **True**), or process from where it stopped. The latter comes in handy when wokng on cloud instances that crash or automatically shutdown every once in a while. \n",
    "\n",
    "**Note** that subsetting is set automatically to **True** when all tracks hold in the inventory fully overlap the AOI."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s1_grd.grds_to_ards(\n",
    "    inventory_df=s1_grd.inventory, \n",
    "    timeseries=True, \n",
    "    timescan=True, \n",
    "    mosaic=False, \n",
    "    overwrite=False\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7* - Merge single per-track time-series to one single time-series\n",
    "\n",
    "By using a little helper function from OST's raster module, combining the 4 time-series to a unique one is as easy the following command. Within your processing directory, a new foler called combined is created. If multi-temporal statistics for this new time-series should be created, set the timescan argument to **True**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "raster.combine_timeseries(\n",
    "    processing_dir=s1_grd.processing_dir, \n",
    "    config_dict=s1_grd.config_dict, \n",
    "    timescan=True\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 8* - Create a time-series animation\n",
    "\n",
    "Finally, a time-series animation is created. therefore we need to pass the time-series folder to the command. product_list expects a list of 1 to 3 elements. For GRD data this is either a single polarisation, or both bands. OST will calculate the power ratio of band 1 and 2 for a 3-band RGB composite. A shrink factor can be set to reduce image resolution and memory needs.\n",
    "\n",
    "**Note:** This needs imagemagick installed, which is not a default requirement by OST.\n",
    "You can install it on e.g. Ubuntu by typing:\n",
    "sudo apt-get install magick"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# define Time-series folder\n",
    "ts_folder = s1_grd.processing_dir.joinpath('combined/Timeseries')\n",
    "\n",
    "# create the animation\n",
    "raster.create_timeseries_animation(\n",
    "    timeseries_folder=ts_folder, \n",
    "    product_list=['bs.VV', 'bs.VH'], \n",
    "    out_folder=s1_grd.processing_dir,\n",
    "    shrink_factor=3, \n",
    "    add_dates=True\n",
    ")"
   ]
  }
 ],
 "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.8.5"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": true,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}