04_Dune_Migration.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "622ee68a-0283-4df8-9381-624576329238",
   "metadata": {},
   "source": [
    "# Dune Migration\n",
    "\n",
    "_**Caitlin Haedrich and Pratikshya Regmi**, North Carolina State University_\n",
    "\n",
    "In this notebook, we're going to introduce temporal datasets in GRASS to see how Jockeys Ridge has moved over the last two decades. We'll create an animation showing the migration."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "569247bd-4525-4a85-8721-59c26bcbd0ab",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "## 1. Import Python Packages and Start GRASS Session"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a4f240ce-ec1c-4cce-9788-9b4690dd5740",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Import Python standard library and IPython packages we need.\n",
    "import subprocess\n",
    "import sys\n",
    "\n",
    "# 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 the GRASS GIS packages we need.\n",
    "import grass.script as gs\n",
    "import grass.jupyter as gj\n",
    "\n",
    "# Start GRASS Session\n",
    "gj.init(\"nags_head/PERMANENT\");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6645dd05-0032-4e10-9713-36c08b2520b3",
   "metadata": {},
   "source": [
    "Before we get started, we set and check the computational region."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "60bf1724-3727-44f1-83a3-7b6957bacf2a",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "!g.region region=jockeys_ridge -p"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f70a0555-8ff9-4423-8349-9aeec056d191",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "## 2. Getting Started with GRASS Temporal Tools\n",
    "\n",
    "\n",
    "To better handle time series of maps, GRASS has temporal datasets which serve as containers for the individual maps representing a certain point or interval of time. We can further manipulate them all together in the temporal dataset instead of the individual maps. For example, we can easily aggregate data at different time intervals, compute univariate statics, or set a universal color scheme.\n",
    "\n",
    "These temporal datasets can contain rasters (space-time raster dataset or _strds_), 3D rasters (_str3ds_) or vectors (_stvds_). Visit [the temporal data processing manual](https://grass.osgeo.org/grass83/manuals/temporalintro.html) for more info.\n",
    "\n",
    "First, we create an empty dataset of type strds (space-time raster dataset)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "100f09ca-47e0-40ec-9384-5481f61f7405",
   "metadata": {},
   "outputs": [],
   "source": [
    "gs.run_command('t.create', output='JockeysRidge', type='strds',\n",
    "                temporaltype='relative', title=\"Jockeys Ridge Elevation Series\",\n",
    "                description=\"from 1996 to 2020 with gaps\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b4cfc45f-cb81-4678-879b-c5372ecf0a26",
   "metadata": {},
   "source": [
    "After defining our temporal dataset, we can add and remove layers from it with `t.register` and `t.unregister`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "73296491-177b-4b6d-aea5-d2e953bb9acb",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get list of rasters we'd like to add to our dataset\n",
    "DEMs = gs.read_command(\"g.list\", type=\"raster\", pattern=\"JR_????\", separator=\"comma\").strip().split(\",\")\n",
    "DEMs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8a38223f-9f7c-4961-bf6b-dc83b23cb0cf",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get the year of each DEM\n",
    "years=[name[-4:] for name in DEMs]\n",
    "years"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "809b22cd-2d2d-45ff-a10f-b291f808ba9c",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Put DEM name and date together in table to register them\n",
    "table = \"\"\n",
    "\n",
    "for DEM, year in zip(DEMs, years):\n",
    "    row = f\"{DEM}|{year}\\n\"\n",
    "    table += row\n",
    "\n",
    "# Print formatted table\n",
    "print(table)\n",
    "\n",
    "# Write the formatted table to a text file\n",
    "with open('./table.txt', 'w') as file:\n",
    "    file.write(table)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "41896d84-065c-4f14-a2a1-04fa29c68964",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Register rasters to dataset\n",
    "gs.run_command(\"t.register\", input=\"JockeysRidge\", file=\"./table.txt\", type=\"raster\", unit=\"years\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d096ca40-b419-4036-b83d-0633bb19c1af",
   "metadata": {},
   "source": [
    "Check what you have now in JockeysRidge dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6d0afb5a-d441-4368-80ee-388c9c4c1272",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(gs.read_command(\"t.rast.list\", input=\"JockeysRidge\", columns=\"name,start_time\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b49f4fdc-d1cd-4b70-900b-f5ec39a5e39b",
   "metadata": {},
   "source": [
    "Set the same color table for all maps."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8c979376-d9f7-45b0-920e-10ba42d02f41",
   "metadata": {},
   "outputs": [],
   "source": [
    "gs.run_command(\"t.rast.colors\", input=\"JockeysRidge\", color=\"elevation\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "77bfc852-83a3-4c37-9753-e9619091ab80",
   "metadata": {},
   "source": [
    "Animate the time series."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "785fc75c-7086-4fa1-8784-f925aded1349",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "!r.mask raster=JR_2014"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f61bb4e3-05fe-4125-8b7c-12b02e6c50cf",
   "metadata": {},
   "outputs": [],
   "source": [
    "animation = gj.TimeSeriesMap()\n",
    "animation.add_raster_series(\"JockeysRidge\", fill_gaps=True)\n",
    "animation.d_legend(color=\"black\", at=(12,72,0,3), fontsize=12)\n",
    "animation.d_barscale(length=250)\n",
    "# animation.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f844da39-7a35-48c6-b68b-b1449a426764",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "animation.save(filename=\"migration.gif\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0f3347f4-d02e-4756-851f-f00e0f25b9ef",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from IPython.display import display, Image\n",
    "\n",
    "display(Image(filename=\"migration.gif\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "86458af3-1773-4348-8fd8-beafcb814cea",
   "metadata": {},
   "source": [
    "**Try it Yourself!**\n",
    "\n",
    "**Question 1**\n",
    "\n",
    "_The animation reveals that the datasets for the years 1996, 1997, 1998, 2000, 2004, 2005, 2006, and 2007 are incomplete. Modify the code above to adjust the DEM list to include only the years with complete datasets._"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac507053-e999-456a-a517-5522574c6631",
   "metadata": {},
   "source": [
    "<details>\n",
    "    <summary>👉 <b>Hint!!</b></summary>\n",
    "\n",
    "```python\n",
    "DEMs=['JR_2001','JR_2008','JR_2014','JR_2016','JR_2018','JR_2020']\n",
    "``` \n",
    "</details>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "50dbf1c5-2e51-4a66-9e1e-cb2b52eab494",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "30aa9970-9e4d-475d-8c6e-85dc82a28192",
   "metadata": {
    "tags": []
   },
   "source": [
    "**Question 2**\n",
    "\n",
    "_Hurricane Floyd was a very powerful Cape Verde hurricane which struck the Bahamas and the East Coast of the United States. It was the sixth named storm, fourth hurricane, and third major hurricane in the 1999 Atlantic hurricane season._\n",
    "\n",
    "![](img/Floyd_banner.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "030a57ef-bbb0-4f9e-bb34-e826e6032291",
   "metadata": {},
   "source": [
    "_We have lidar datasets that were taken before and after hurricane (JR_1999_0909, JR_1999_0918). Using the dataset provided, calculate the difference in landmass volume before and after the impact of Hurricane Floyd. Also, compute the net volume change._"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b351c1d7-c04a-414f-8105-c62a2925a56e",
   "metadata": {},
   "source": [
    "<details>\n",
    "    <summary>👉 <b>Hint 1</b></summary>\n",
    "\n",
    "```python\n",
    "!r.mapcalc \"elevation_difference = JR_1999_1104 - JR_1998\"\n",
    "```\n",
    "\n",
    "</details>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7bc1ad23-a386-498f-9155-e2f2a5d5461a",
   "metadata": {
    "tags": []
   },
   "source": [
    "<details>\n",
    "    <summary>👉 <b>Hint 2</b></summary>\n",
    "\n",
    "```python\n",
    "!r.colors map=elevation_difference color=difference\n",
    "example = gj.Map()\n",
    "example.d_background(color=\"white\")\n",
    "example.d_rast(map=\"elevation_difference\")\n",
    "example.d_legend(raster=\"elevation_difference\")\n",
    "example.show()\n",
    "``` \n",
    "\n",
    "<details>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a0876c6-7596-45f1-bf7b-e35b2a8fe0cb",
   "metadata": {},
   "source": [
    "<details>\n",
    "    <summary>👉 <b>Hint 3</b></summary>\n",
    "    \n",
    "```python\n",
    "stats= gs.parse_command(\"r.univar\", map='elevation_difference', flags='g')\n",
    "stats['sum']\n",
    "print(f\"Net sand change is {stats ['sum']} cubic meter\")\n",
    "```\n",
    "<details>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "02754376-2c2b-4c26-9e87-f7afc185ad4b",
   "metadata": {},
   "source": [
    "## Further Readings\n",
    "\n",
    "[Mitasova, H., Overton, M., Recalde, J.J., Bernstein, D., and Freeman C., 2009, Raster-based Analysis of Coastal Terrain Dynamics from Multitemporal Lidar Data, Journal of Coastal Research 25(2), p. 507-514, DOI: 10.2112/07-0976.1.](https://meridian.allenpress.com/jcr/article-abstract/25/2%20(252)/507/28283/Raster-Based-Analysis-of-Coastal-Terrain-Dynamics)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "GRASS GIS",
   "language": "python",
   "name": "grass"
  },
  "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.12.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}