0.OGR.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hello \"OGR\" World!\n",
    "\n",
    "\n",
    "This tutorial is designed to give you an overview of the OGR module.\n",
    "\n",
    "OGR is part of the [GDAL module](https://gdal.org/). \n",
    "\n",
    "Let us start by importing the OGR module"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ogr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Jupyter Notebook allows you to run commands (just the way you would do in your command prompt / terminal) by prefixing the exclaimation mark.\n",
    "\n",
    "The \"ogrinfo --formats\" command lists [all supported GIS-Vector file formats](https://gdal.org/drivers/vector/index.html)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!ogrinfo --formats"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "+ The first thing you need to do is select a **\"driver\"**.\n",
    "\n",
    "+ A driver is an object that knows how to interact with a certain data type (such as a shapefile).\n",
    "\n",
    "+ Need an appropriate driver in order to read or write data (need it explicitly for write).\n",
    "\n",
    "+ Each file type (shapefile, KML etc.) has a unique \"Driver short name\". \n",
    "\n",
    "+ You can get that information from the [corresponding GDAL/OGR page](https://gdal.org/drivers/vector/shapefile.html).  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "driver = ogr.GetDriverByName('ESRI Shapefile')\n",
    "print(driver)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "+ Next, we specify the file path and specify the file open mode (0 = read; 1 = write).\n",
    "\n",
    "+ If ds (datasource) is None, then there is an issue with file opening."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = driver.Open('data/sites.shp', 0)\n",
    "print(ds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "layer = ds.GetLayer(0) #or ds.GetLayer() #the number is optional for shapefiles.\n",
    "layer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#To get the number of features in the file\n",
    "n = layer.GetFeatureCount()\n",
    "print ('feature count:', n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#To get the bounding box of the file\n",
    "extent = layer.GetExtent()\n",
    "print ('extent:', extent)\n",
    "print ('ul:', extent[0], extent[3])\n",
    "print ('lr:', extent[1], extent[2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![Vector Data Types](https://www.earthdatascience.org/images/courses/earth-analytics/spatial-data/points-lines-polygons-vector-data-types.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![Spatial Data Attributes](https://www.earthdatascience.org/images/courses/earth-analytics/spatial-data/spatial-attribute-tables.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![Bounding box](https://www.earthdatascience.org/images/courses/earth-analytics/spatial-data/spatial-extent.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#You can randomly select features, in case you know the FID. \n",
    "#This file has 42 features (0..41).\n",
    "#Let us change the value from 41 to 42, and see what happens.\n",
    "\n",
    "#Let us read the first feature \n",
    "feat = layer.GetFeature(0)\n",
    "print(feat.GetField('id')) #should be a different id\n",
    "\n",
    "feat = layer.GetFeature(41) #as count starts from 0\n",
    "print(feat.GetField('id'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Let us read the next feature (from the beginning)\n",
    "feat = layer.GetNextFeature()\n",
    "print(feat.GetField('id')) #what will it print?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "while feat:\n",
    "    print (feat.GetField('id'))\n",
    "    feat.Destroy()\n",
    "    feat = layer.GetNextFeature() #what happens if I forget this?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Things to try**\n",
    "1. Do a layer.GetNextFeature now. You get None. Why?\n",
    "2. ResetReading, and comment feat = layer.GetNextFeature in the previous cell and run. Why did it hang?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat = layer.GetNextFeature()\n",
    "print (feat.GetField('id'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#If you have reached the end of the file, you need to reset the file pointer to the start.\n",
    "layer.ResetReading()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#You can change the data type of the value being retrieved from the attribute table\n",
    "#Change the GetField function to GetFieldAsString and see the result\n",
    "feat = layer.GetNextFeature()\n",
    "value = feat.GetField('id')\n",
    "print(value)\n",
    "print (type(value))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "+ Feature objects have a method called GetGeometryRef() which returns a Geometry object (could be Point, Polygon, etc)\n",
    "+ Point objects have GetX()and GetY() methods"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geometry = feat.GetGeometryRef()\n",
    "geom_type = geometry.GetGeometryType()\n",
    "x = geometry.GetX()\n",
    "y = geometry.GetY()\n",
    "print('x = ', x)\n",
    "print('y = ', y)\n",
    "print('Geometry type = ', geom_type)\n",
    "#print(ogr.wkbPoint, ogr.wkbLineString, ogr.wkbPolygon)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Classwork**: Read coordinates and attributes from a shapefile\n",
    "+ Loop through the points in sites.shp\n",
    "+ Print out id, x & y coordinates, and cover type for each point to a text file, one point per line\n",
    "\n",
    "**Hint:** The two attribute fields in the shapefile are called \"id\" and \"cover\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat.Destroy() #Only for memory management. All variables get killed on session/kernal expiry."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds.Destroy() #Only for memory management. All variables get killed on session/kernal expiry."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Homework**:\n",
    "1. Complete the [slides Nos.](https://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_slides1.pdf) 27-38\n",
    "2. Try out the [assignments](https://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_hw1.pdf). (there is some overlap with point 1)\n",
    "3. Explore the numerous things you could do with [Geometry](http://pcjericks.github.io/py-gdalogr-cookbook/geometry.html) and [Vector Layers](http://pcjericks.github.io/py-gdalogr-cookbook/vector_layers.html)"
   ]
  }
 ],
 "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}