ndwi.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "otherwise-capability",
   "metadata": {},
   "source": [
    "# NDWI analysis with Google Earth Engine"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "unauthorized-contents",
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [],
   "source": [
    "# the regulars\n",
    "import geopandas as gpd\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# work with dates\n",
    "from datetime import datetime as dt\n",
    "\n",
    "# allow images to display\n",
    "from IPython.display import Image\n",
    "\n",
    "# earth engine\n",
    "import ee"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "existing-plastic",
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [],
   "source": [
    "# Trigger the authentication flow.\n",
    "ee.Authenticate()\n",
    "\n",
    "# Initialize the library.\n",
    "ee.Initialize()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "documentary-anniversary",
   "metadata": {},
   "source": [
    "## Set some parameters\n",
    "\n",
    "Change these parameters for your area and timeframe of interest"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "right-detroit",
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [],
   "source": [
    "# coordinates of the 2020 Kyushu Region Floods\n",
    "lat = 38.46\n",
    "lon = 140.35\n",
    "\n",
    "# point of interest as an ee.Geometry\n",
    "poi = ee.Geometry.Point(lon,lat)\n",
    "\n",
    "# start date of range to filter for\n",
    "start_date = '2020-06-01'\n",
    "\n",
    "# end date\n",
    "end_date = '2020-09-30'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "computational-repository",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# Define a region of interest with a buffer zone\n",
    "roi = poi.buffer(5000) # meters"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "graduate-occasions",
   "metadata": {},
   "source": [
    "## Sentinel-2 data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "soviet-volume",
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [],
   "source": [
    "# get the satellite data (Sentinel)\n",
    "sentinel = ee.ImageCollection(\"COPERNICUS/S2_SR\")\\\n",
    "            .filterBounds(roi)\\\n",
    "            .filterDate(start_date,end_date)\\\n",
    "            .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "sublime-monroe",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# how many images did we get?\n",
    "print('Total number:', sentinel.size().getInfo())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "reverse-hungary",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# put the images in a list\n",
    "sentinel_list = sentinel.toList(sentinel.size());"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "consolidated-million",
   "metadata": {},
   "outputs": [],
   "source": [
    "# set some parameters for image\n",
    "parameters = {\n",
    "                'min': 0,\n",
    "                'max': 2500,\n",
    "                'dimensions': 800,\n",
    "                'bands': ['B4', 'B3', 'B2'],\n",
    "                'region':roi\n",
    "             }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "furnished-traffic",
   "metadata": {
    "scrolled": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# loop through each image and display it\n",
    "for i in range(sentinel.size().getInfo()):\n",
    "\n",
    "    # when was this image taken?\n",
    "    date = ee.Date(ee.Image(sentinel_list.get(i)).get('system:time_start'))\n",
    "    time = date.getInfo()['value']/1000\n",
    "    date = dt.utcfromtimestamp(time).strftime('%Y-%m-%d %H:%M:%S')\n",
    "    \n",
    "    # cloud cover\n",
    "    cloud = ee.Image(sentinel_list.get(i)).get('CLOUDY_PIXEL_PERCENTAGE').getInfo()\n",
    "    \n",
    "    print('Image #',i,date,'Cloud cover:',cloud)\n",
    "    \n",
    "    display(Image(url = ee.Image(sentinel_list.get(i)).getThumbUrl(parameters)))\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "incredible-morning",
   "metadata": {},
   "source": [
    "## Normalized Difference Water Index (NDWI)\n",
    "\n",
    "NDWI = (green – NIR) / (green + NIR)\n",
    "\n",
    "For Sentinel data, green = B3 and the NIR is B8"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "worth-dietary",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# ndwi palette: red is low water, blue is high water\n",
    "palette = ['red', 'yellow', 'green', 'blue'];\n",
    "\n",
    "# notice, no bands here because we will calculate ndwi per image\n",
    "ndwi_parameters = {'min': -1,\n",
    "                   'max': 0.5,\n",
    "                   'dimensions': 800,\n",
    "                   'palette': palette,\n",
    "                   'region': roi}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "capital-lewis",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# loop through each image and display it\n",
    "for i in range(sentinel.size().getInfo()):\n",
    "\n",
    "    # first, calculate ndwi for each image\n",
    "    ndwi = ee.Image(sentinel_list.get(i)).normalizedDifference(['B3', 'B8'])\n",
    "    \n",
    "    # when was this image taken?\n",
    "    date = ee.Date(ee.Image(sentinel_list.get(i)).get('system:time_start'))\n",
    "    time = date.getInfo()['value']/1000\n",
    "    date = dt.utcfromtimestamp(time).strftime('%Y-%m-%d %H:%M:%S')\n",
    "    \n",
    "    # print some information\n",
    "    print('Image #',i,date)\n",
    "    \n",
    "    # display the image\n",
    "    display(Image(url=ndwi.getThumbUrl(ndwi_parameters)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "latest-marina",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import folium\n",
    "from folium import plugins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "incorporated-going",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# Google function that allows ee layers on folium\n",
    "def add_ee_layer(self, ee_image_object, vis_params, name):\n",
    "    \"\"\"Adds a method for displaying Earth Engine image tiles to folium map.\"\"\"\n",
    "    map_id_dict = ee.Image(ee_image_object).getMapId(vis_params)\n",
    "    folium.raster_layers.TileLayer(\n",
    "        tiles=map_id_dict['tile_fetcher'].url_format,\n",
    "        attr='Map Data &copy; <a href=\"https://earthengine.google.com/\">Google Earth Engine</a>',\n",
    "        name=name,\n",
    "        overlay=True,\n",
    "        control=True\n",
    "    ).add_to(self)\n",
    "\n",
    "# Add Earth Engine drawing method to folium\n",
    "folium.Map.add_ee_layer = add_ee_layer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "retained-clone",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Add custom basemaps to folium\n",
    "basemaps = {\n",
    "    'Google Maps': folium.TileLayer(\n",
    "        tiles = 'https://mt1.google.com/vt/lyrs=m&x={x}&y={y}&z={z}',\n",
    "        attr = 'Google',\n",
    "        name = 'Google Maps',\n",
    "        overlay = True,\n",
    "        control = True\n",
    "    ),\n",
    "    'Google Satellite': folium.TileLayer(\n",
    "        tiles = 'https://mt1.google.com/vt/lyrs=s&x={x}&y={y}&z={z}',\n",
    "        attr = 'Google',\n",
    "        name = 'Google Satellite',\n",
    "        overlay = True,\n",
    "        control = True\n",
    "    ),\n",
    "    'Google Terrain': folium.TileLayer(\n",
    "        tiles = 'https://mt1.google.com/vt/lyrs=p&x={x}&y={y}&z={z}',\n",
    "        attr = 'Google',\n",
    "        name = 'Google Terrain',\n",
    "        overlay = True,\n",
    "        control = True\n",
    "    ),\n",
    "    'Google Satellite Hybrid': folium.TileLayer(\n",
    "        tiles = 'https://mt1.google.com/vt/lyrs=y&x={x}&y={y}&z={z}',\n",
    "        attr = 'Google',\n",
    "        name = 'Google Satellite',\n",
    "        overlay = True,\n",
    "        control = True\n",
    "    ),\n",
    "    'Esri Satellite': folium.TileLayer(\n",
    "        tiles = 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}',\n",
    "        attr = 'Esri',\n",
    "        name = 'Esri Satellite',\n",
    "        overlay = True,\n",
    "        control = True\n",
    "    )\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "previous-summit",
   "metadata": {
    "scrolled": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "# Create a map\n",
    "my_map = folium.Map(location=[lat, lon], zoom_start=10)\n",
    "\n",
    "# Add custom basemaps\n",
    "basemaps['Google Maps'].add_to(my_map)\n",
    "basemaps['Google Satellite Hybrid'].add_to(my_map)\n",
    "\n",
    "\n",
    "# loop through each image and display it\n",
    "for i in range(sentinel.size().getInfo()):\n",
    "\n",
    "    # first, calculate ndwi for each image\n",
    "    ndwi = ee.Image(sentinel_list.get(i)).normalizedDifference(['B3', 'B8'])\n",
    "    \n",
    "    # when was this image taken?\n",
    "    date = ee.Date(ee.Image(sentinel_list.get(i)).get('system:time_start'))\n",
    "    time = date.getInfo()['value']/1000\n",
    "    date = dt.utcfromtimestamp(time).strftime('%Y-%m-%d %H:%M:%S')\n",
    "        \n",
    "    my_map.add_ee_layer(ndwi, \n",
    "                        ndwi_parameters, \n",
    "                        name=date)\n",
    "\n",
    "\n",
    "# Add a layer control panel to the map\n",
    "folium.LayerControl(collapsed = False).add_to(my_map)\n",
    "\n",
    "# Display the map.\n",
    "display(my_map)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "athletic-storm",
   "metadata": {},
   "outputs": [],
   "source": [
    "my_map.save('kumamoto.html')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eleven-henry",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "modular-writer",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}