grass.jupyter: timeseries visualizations

doc/notebooks/temporal.ipynb

@@ -0,0 +1,318 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Spatio-Temporal Analysis with grass.jupyter\n",
+    "\n",
+    "As part of a GRASS mini grant, we've been adding visualization functions for time space datasets (strds and stvds). You can find out more about the project and follow the progress on the [GRASS wiki page](https://trac.osgeo.org/grass/wiki/GSoC/2021/JupyterAndGRASS/MiniGrant2022).\n",
+    "\n",
+    "This interactive notebook is available online thanks to the [https://mybinder.org](Binder) service. To run the select part (called a *cell*), hit `Shift + Enter`."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Download Dataset"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import subprocess\n",
+    "import sys\n",
+    "import shutil\n",
+    "\n",
+    "# Download zip\n",
+    "!curl http://fatra.cnr.ncsu.edu/temporal-grass-workshop/NC_spm_temporal_workshop.zip -o ../../data/NC_spm_temporal_workshop.zip\n",
+    "\n",
+    "# Unpack zip to grassdata\n",
+    "shutil.unpack_archive(\"../../data/NC_spm_temporal_workshop.zip\", \"../../data/grassdata\", \"zip\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Start GRASS GIS"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Ask GRASS GIS where its Python packages are.\n",
+    "sys.path.append(\n",
+    "    subprocess.check_output([\"grass\", \"--config\", \"python_path\"], text=True).strip()\n",
+    ")\n",
+    "\n",
+    "# Import GRASS packages\n",
+    "import grass.script as gs\n",
+    "import grass.jupyter as gj\n",
+    "\n",
+    "# Start GRASS Session\n",
+    "gj.init(\"../../data/grassdata\", \"NC_spm_temporal_workshop\", \"climate_2000_2012\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Set computational region to the elevation raster."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\"g.region\", raster=\"elev_state_500m\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Create empty space time datasets"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\n",
+    "    \"t.create\",\n",
+    "    output=\"tempmean\",\n",
+    "    type=\"strds\",\n",
+    "    temporaltype=\"absolute\",\n",
+    "    title=\"Average temperature\",\n",
+    "    description=\"Monthly temperature average in NC [deg C]\",\n",
+    ")\n",
+    "\n",
+    "gs.run_command(\n",
+    "    \"t.create\",\n",
+    "    output=\"precip_sum\",\n",
+    "    type=\"strds\",\n",
+    "    temporaltype=\"absolute\",\n",
+    "    title=\"Preciptation\",\n",
+    "    description=\"Monthly precipitation sums in NC [mm]\",\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Create list of rasters to be registered to empty space time datasets"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "tempmean_list = gs.read_command(\n",
+    "    \"g.list\", type=\"raster\", pattern=\"*tempmean\", separator=\"comma\"\n",
+    ").strip()\n",
+    "\n",
+    "precip_list = gs.read_command(\n",
+    "    \"g.list\", type=\"raster\", pattern=\"*precip\", separator=\"comma\"\n",
+    ").strip()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Register the rasters to the space time dataset created above"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\n",
+    "    \"t.register\",\n",
+    "    input=\"tempmean\",\n",
+    "    type=\"raster\",\n",
+    "    start=\"2000-01-01\",\n",
+    "    increment=\"1 months\",\n",
+    "    maps=tempmean_list,\n",
+    "    flags=\"i\",\n",
+    ")\n",
+    "\n",
+    "gs.run_command(\n",
+    "    \"t.register\",\n",
+    "    input=\"precip_sum\",\n",
+    "    type=\"raster\",\n",
+    "    start=\"2000-01-01\",\n",
+    "    increment=\"1 months\",\n",
+    "    maps=precip_list,\n",
+    "    flags=\"i\",\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Extract a small subset for visualization"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\n",
+    "    \"t.rast.extract\",\n",
+    "    input=\"precip_sum\",\n",
+    "    output=\"precip_sum_2010\",\n",
+    "    where=\"start_time >= '2010-01-01' and start_time < '2011-01-01'\",\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Set the color table for all rasters in timeseries"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\"t.rast.colors\", input=\"precip_sum_2010\", color=\"precipitation_monthly\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Temporal Visualizations\n",
+    "\n",
+    "The `TimeSeries` class contains visualization functions for GRASS space time dataset (strds or stvds). The `time_slider` function allows users to interactively view the evolution of the dataset through time using IPython and Jupyter Widgets."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "img = gj.TimeSeries(\"precip_sum_2010\", fill_gaps=False, use_region=True)\n",
+    "img.d_legend(color=\"black\", at=(10,40,2,6)) #Add legend\n",
+    "img.overlay.d_vect(map=\"boundary_county\", fill_color=\"none\")\n",
+    "img.overlay.d_barscale()\n",
+    "img.time_slider()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "We can also display the TimeSeries as an animation with `animate`."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "img.animate(duration=500, label=True, text_size=16, text_color=\"gray\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Here is the same example but with two of the rasters unregistered, creating a dataset with variable timesteps."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gs.run_command(\n",
+    "    \"t.unregister\",\n",
+    "    type=\"raster\",\n",
+    "    input=\"precip_sum_2010\",\n",
+    "    maps=\"2010_02_precip,2010_08_precip\",\n",
+    ")\n",
+    "print(\n",
+    "    gs.read_command(\"t.rast.list\", input=\"precip_sum_2010\", columns=\"name,start_time\")\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "img2 = gj.TimeSeries(\"precip_sum_2010\")\n",
+    "img2.d_legend(color=\"gray\", at=(10, 0, 30, 0))  # Add legend\n",
+    "img2.time_slider()  # Create TimeSlider"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "By default, `fill_gaps=False`, so there is are blank images where we removed rasters. By setting `fill_gaps=True`, we will see the gap filled by the previous time step."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "img3 = gj.TimeSeries(\"precip_sum_2010\", fill_gaps=True)\n",
+    "img3.d_legend(color=\"gray\", at=(10, 0, 30, 0))  # Add legend\n",
+    "img3.time_slider()  # Create TimeSlider"
+   ]
+  }
+ ],
+ "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.10"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}

python/grass/jupyter/Makefile

@@ -11,10 +11,11 @@ MODULES = \
 	interact_display \
 	region \
 	render3d \
-    reprojection_renderer \
+  reprojection_renderer \
+	utils \
+	timeseries \
 	utils
 
-
 PYFILES := $(patsubst %,$(DSTDIR)/%.py,$(MODULES) __init__)
 PYCFILES := $(patsubst %,$(DSTDIR)/%.pyc,$(MODULES) __init__)
 

python/grass/jupyter/__init__.py

@@ -46,4 +46,5 @@
 from .render3d import *
 from .setup import *
 from .utils import *
+from .timeseries import *
 from .reprojection_renderer import *

python/grass/jupyter/display.py

@@ -30,6 +30,7 @@ class GrassRenderer:
     Elements are added to the display by calling GRASS display modules.
 
     Basic usage::
+
     >>> m = GrassRenderer()
     >>> m.run("d.rast", map="elevation")
     >>> m.run("d.legend", raster="elevation")
@@ -39,10 +40,12 @@ class GrassRenderer:
     as a class method and replacing "." with "_" in the name.
 
     Shortcut usage::
+
     >>> m = GrassRenderer()
     >>> m.d_rast(map="elevation")
     >>> m.d_legend(raster="elevation")
     >>> m.show()
+
     """
 
     def __init__(
@@ -56,6 +59,7 @@ def __init__(
         renderer="cairo",
         use_region=False,
         saved_region=None,
+        read_file=False,
     ):
 
         """Creates an instance of the GrassRenderer class.
@@ -71,9 +75,12 @@ def __init__(
         :param int text_size: default text size, overwritten by most display modules
         :param renderer: GRASS renderer driver (options: cairo, png, ps, html)
         :param use_region: if True, use either current or provided saved region,
-                          else derive region from rendered layers
+                        else derive region from rendered layers
         :param saved_region: if name of saved_region is provided,
-                            this region is then used for rendering
+                        this region is then used for rendering
+        :param bool read_file: if False (default), erase filename before re-writing to
+                         clear contents. If True, read file without clearing contents
+                         first.
         """
 
         # Copy Environment
@@ -107,6 +114,8 @@ def cleanup(tmpdir):
 
         if filename:
             self._filename = filename
+            if not read_file and os.path.exists(self._filename):
+                os.remove(self._filename)
         else:
             self._filename = os.path.join(self._tmpdir.name, "map.png")
         # Set environment var for file