add interactive rasterio, matplotlib tutorials

Author: Sara Safavi

jupyter-notebooks/inspecting-sat-imagery/Inspecting Satellite Imagery.ipynb renamed to jupyter-notebooks/getting-to-know-sat-imagery/Inspecting Satellite Imagery.ipynb

@@ -5,21 +5,26 @@
    "metadata": {},
    "source": [
     "# Inspecting Satellite Imagery using Rasterio\n",
-    "## A first look at analyzing satellite data with Python\n",
+    "## A first look at satellite data with Python\n",
     "\n",
-    "At this point, you've explored different ways of searching for, filtering, and downloading satellite imagery. Now let's use one of these acquired datasets and dig into it a bit with Python.\n",
+    "At this point, you will have learned different ways of searching for, filtering, and downloading satellite imagery. Now let's use one of these acquired datasets and dig into it a bit with Python.\n",
     "\n",
-    "Here we're going to use a Python library called `rasterio`: you may be familiar with it already, or perhaps with the related C library, `GDAL`."
+    "Here we're going to use a Python library called `rasterio`: you may be familiar with it already, or perhaps with the related C library, `GDAL`. If you've used `numpy` before, working with `rasterio` will feel very familiar."
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 50,
+   "execution_count": null,
    "metadata": {},
    "outputs": [],
    "source": [
     "import rasterio\n",
-    "satdat = rasterio.open(\"example.tif\")"
+    "\n",
+    "# feel free to replace 'example' with your own GeoTIFF:\n",
+    "# note that this Notebook will assume we're working with PlanetScope 4-band imagery\n",
+    "image_file = \"example.tif\"\n",
+    "\n",
+    "satdat = rasterio.open(image_file)"
    ]
   },
   {
@@ -32,40 +37,22 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 51,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "BoundingBox(left=540759.0, bottom=3754401.0, right=568047.0, top=3767985.0)"
-      ]
-     },
-     "execution_count": 51,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
+   "outputs": [],
    "source": [
     "# Get the bounding box of this GeoTIFF\n",
+    "\n",
     "satdat.bounds"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 63,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "Width: 27288.0, Height: 13584.0\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
-    "# Get dimensions, in map units (here, that's meters)\n",
+    "# Get dimensions, in map units (using the example GeoTIFF, that's meters)\n",
     "\n",
     "width = satdat.bounds.right - satdat.bounds.left\n",
     "height = satdat.bounds.top - satdat.bounds.bottom\n",
@@ -75,40 +62,22 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 64,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "Width: 9096, Height: 4528\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
     "# Get dimensions, in pixels\n",
+    "\n",
     "px_width = satdat.width\n",
     "px_height = satdat.height\n",
     "print(\"Width: {}, Height: {}\".format(px_width, px_height))"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 65,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "(3.0, 3.0)"
-      ]
-     },
-     "execution_count": 65,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
+   "outputs": [],
    "source": [
     "# How many meters to a pixel?\n",
     "\n",
@@ -120,70 +89,51 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 66,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "CRS({'init': 'epsg:32611'})"
-      ]
-     },
-     "execution_count": 66,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
+   "outputs": [],
    "source": [
     "# Get coordinate reference system\n",
+    "\n",
     "satdat.crs"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 59,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "Top left corner coordinates: (540759.0, 3767985.0)\n",
-      "Bottom right corner coordinates: (568047.0, 3754401.0)\n"
-     ]
-    },
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "/opt/conda/lib/python3.6/site-packages/IPython/core/interactiveshell.py:2910: FutureWarning: The value of this property will change in version 1.0. Please see https://github.com/mapbox/rasterio/issues/86 for details.\n",
-      "  exec(code_obj, self.user_global_ns, self.user_ns)\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
     "# Get coordinates of top-left & bottom right corners\n",
-    "# NOTE: how to do this depends on your Rasterio version:\n",
-    "# the below example may generate a FutureWarning, which is safe to ignore here\n",
     "\n",
-    "try:\n",
-    "    topleft = satdat.transform * (0, 0)\n",
-    "    botright = satdat.transform * (width, height)\n",
-    "    \n",
-    "except TypeError:\n",
-    "    topleft = satdat.affine * (0, 0)\n",
-    "    botright = satdat.affine * (width, height)\n",
+    "# NOTE: how to do this depends on your Rasterio version --\n",
+    "# if you're running against a pre-1.0 release use satdat.affine instead\n",
+    "\n",
+    "topleft = satdat.transform * (0, 0)\n",
+    "botright = satdat.transform * (width, height)\n",
     "    \n",
     "print(\"Top left corner coordinates: {}\".format(topleft))\n",
     "print(\"Bottom right corner coordinates: {}\".format(botright))"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Another way of viewing most of the previous values:\n",
+    "# Get the basic metadata of this GeoTIFF\n",
+    "\n",
+    "satdat.meta"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
     "## Bands\n",
-    "So far, we haven't done too much raster-specific work yet. Since we know we're inspecting a multispectral satellite image, let's see what we can learn about its bands."
+    "So far, we haven't done too much geospatial-raster-specific work yet. Since we know we're inspecting a multispectral satellite image, let's see what we can learn about its bands."
    ]
   },
   {
@@ -193,6 +143,7 @@
    "outputs": [],
    "source": [
     "# Get the number of bands by listing their indices\n",
+    "\n",
     "satdat.indexes"
    ]
   },
@@ -214,7 +165,21 @@
     "blue = satdat.read(1)\n",
     "green = satdat.read(2)\n",
     "red = satdat.read(3)\n",
-    "nir = satdat.read(4)"
+    "nir = satdat.read(4)\n",
+    "\n",
+    "# or:\n",
+    "# blue, green, red, nir = satdat.read()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Pixels\n",
+    "\n",
+    "In a raster dataset, each pixel has a value. Pixels are arranged in a grid, and pixels representing equivalent data have the same value:\n",
+    "\n",
+    "![pixels2.png](pixels2.png)"
    ]
   },
   {
@@ -225,7 +190,7 @@
    "source": [
     "# bands are stored as numpy arrays\n",
     "\n",
-    "type(blue)"
+    "print(type(blue))"
    ]
   },
   {
@@ -234,7 +199,9 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "blue"
+    "# how many dimensions would a raster have?\n",
+    "\n",
+    "print(blue.ndim)"
    ]
   },
   {
@@ -243,56 +210,48 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "# Output min & max pixel values in each band\n",
+    "# take a look at the summarized array\n",
     "\n",
-    "print(blue.min(), blue.max())\n",
-    "print(green.min(), green.max())\n",
-    "print(red.min(), red.max())\n",
-    "print(nir.min(), nir.max())"
+    "print(blue)"
    ]
   },
   {
-   "cell_type": "markdown",
+   "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
-    "## Pixels"
+    "# Output a min & max pixel value in each band\n",
+    "\n",
+    "# we'll need to use numpy directly for this\n",
+    "import numpy\n",
+    "\n",
+    "print(numpy.amin(blue), numpy.amax(blue))\n",
+    "print(numpy.amin(green), numpy.amax(green))\n",
+    "print(numpy.amin(red), numpy.amax(red))\n",
+    "print(numpy.amin(nir), numpy.amax(nir))"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 61,
+   "execution_count": null,
    "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "6382"
-      ]
-     },
-     "execution_count": 61,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
+   "outputs": [],
    "source": [
-    "# Let's grab the pixel 5km east and 5km south of the upper left corner\n",
+    "# Let's grab the pixel 2km east and 2km south of the upper left corner\n",
     "\n",
     "# get the pixel \n",
-    "px_x = satdat.bounds.left + 5000\n",
-    "px_y = satdat.bounds.top - 5000\n",
+    "px_x = satdat.bounds.left + 2000\n",
+    "px_y = satdat.bounds.top - 2000\n",
     "\n",
     "row, col = satdat.index(px_x, px_y)\n",
     "\n",
-    "# Now let's look at the value of Band1 (\"blue\") at this pixel\n",
-    "blue[row, col]"
+    "# Now let's look at the value of each band at this pixel\n",
+    "print(\"Red: {}\".format(red[row, col]))\n",
+    "print(\"Green: {}\".format(green[row, col]))\n",
+    "print(\"Blue: {}\".format(blue[row, col]))\n",
+    "print(\"NIR: {}\".format(nir[row, col]))"
    ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": []
   }
  ],
  "metadata": {
@@ -311,7 +270,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.6.2"
+   "version": "3.6.3"
   }
  },
  "nbformat": 4,