DelaunayTriangulation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Delaunay Triangulation and Voronoi Diagrams\n",
    "\n",
    "- Overview\n",
    "- What is Delaunay Triangulation?\n",
    "- What is a Voronoi Diagram ?\n",
    "- Delaunay Triangulation Code\n",
    "- Delaunay Triangulation Animation Code"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# What is Delaunay Triangulation?\n",
    "\n",
    "1. Cover entire space\n",
    "2. No point should be inside the circumcircle of any triangle\n",
    "\n",
    "![Delaunay Triangulation](data/images/opcv4face-w4-m3-dealunyTriangle-1-small.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "![Delaunay triangulation favors small angles](data/images/opcv4face-w4-m3-angles.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# What is a Voronoi Diagram?\n",
    "\n",
    "![Voronoi Diagram](data/images/opcv4face-w4-m3-voronoiDiagram-small.png)\n",
    "\n",
    "> If you connect the points in neighboring Voronoi regions, you get a Delaunay triangulation!\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "# Subdiv2D : OpenCV implementation of Delaunay Triangulation\n",
    "\n",
    "\n",
    "## OpenCV Output \n",
    "\n",
    "> x10 y10 x22 y22 x30 y30\n",
    "\n",
    "## We want \n",
    "\n",
    "> 10 22 30\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 1: Import required module"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import cv2\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 2: Find index of closest landmark point\n",
    "\n",
    "The function `findIndex` finds the index of the closest landmark point to the given vertex from the vector of points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def findIndex(points, point):\n",
    "  diff = np.array(points) - np.array(point)\n",
    "\n",
    "  # Find the distance of point from all points\n",
    "  diffNorm = np.linalg.norm(diff, 2, 1)\n",
    "\n",
    "  # Find the index with minimum distance and return it\n",
    "  return np.argmin(diffNorm)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 3: Write Delauanay Triangles to a file\n",
    "\n",
    "The function `writeDelaunay` gets the Delaunay triangles using the `getTriangleList` method of the `subdiv` object.\n",
    "\n",
    "`getTriangleList` returns a array of 6 floating point numbers. The six numbers are the `x` and `y` coordinates of the three vertices of a triangle.\n",
    "\n",
    "We are interested in storing the output file as a list of indices of the vertices and not the actual vertices. So, we loop through all the triangles returned by `getTriangleList` and find the index of the vertices.\n",
    "\n",
    "Each row of the output file contains three number that represent a triangle. The numbers are the indices of points in the input file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "# write delaunay triangles to file\n",
    "def writeDelaunay( subdiv, points, outputFileName ) :\n",
    "\n",
    "    # Obtain the list of triangles.\n",
    "    # Each triangle is stored as vector of 6 coordinates\n",
    "    # (x0, y0, x1, y1, x2, y2)\n",
    "    triangleList = subdiv.getTriangleList();\n",
    "    filePointer = open(outputFileName,'w')\n",
    "    # Will convert triangle representation to three vertices pt1, pt2, pt3\n",
    "    for t in triangleList :\n",
    "        pt1 = (t[0], t[1])\n",
    "        pt2 = (t[2], t[3])\n",
    "        pt3 = (t[4], t[5])\n",
    "\n",
    "        # Find the landmark corresponding to each vertex\n",
    "        landmark1 = findIndex(points,pt1)\n",
    "        landmark2 = findIndex(points,pt2)\n",
    "        landmark3 = findIndex(points,pt3)\n",
    "\n",
    "        filePointer.writelines(\"{} {} {}\\n\".format(landmark1, landmark2, landmark3 ))\n",
    "\n",
    "    filePointer.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 4: Read image\n",
    "\n",
    "We specify the colors for drawing lines and points on the image. Then we read the image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# Define colors for drawing.\n",
    "delaunayColor = (255,255,255)\n",
    "pointsColor = (0, 0, 255)\n",
    "\n",
    "# Read in the image.\n",
    "img = cv2.imread(\"data/images/smiling-man.jpg\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "# Step 5: Calculate Delaunay Triangles\n",
    "\n",
    "Delaunay triangulation and voronoi diagrams are calculated using the class `Subdiv2D`.\n",
    "\n",
    "`Subdiv2D` is initialized with a rectangle containing all the points. The landmark points are read from the file.\n",
    "\n",
    "The `insert` method is used to insert points.\n",
    "\n",
    "`writeDelaunay` function uses `getTriangleList` method to get the delaunay triangles."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "# Rectangle to be used with Subdiv2D\n",
    "size = img.shape\n",
    "rect = (0, 0, size[1], size[0])\n",
    "\n",
    "# Create an instance of Subdiv2D\n",
    "subdiv = cv2.Subdiv2D(rect);\n",
    "\n",
    "# Create an array of points.\n",
    "points = []\n",
    "\n",
    "# Read in the points from a text file\n",
    "with open(\"data/images/smiling-man-delaunay.txt\") as file :\n",
    "    for line in file:\n",
    "        x, y = line.split()\n",
    "        points.append((int(x), int(y)))\n",
    "\n",
    "outputFileName = \"results/smiling-man-delaunay.tri\"\n",
    "\n",
    "# Insert points into subdiv\n",
    "for p in points:\n",
    "    subdiv.insert(p)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "# Step 6: Write Delaunay Triangles to file\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Writing Delaunay triangles to results/smiling-man-delaunay.tri\n"
     ]
    }
   ],
   "source": [
    "writeDelaunay(subdiv, points, outputFileName)\n",
    "print(\"Writing Delaunay triangles to {}\".format(outputFileName))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# Delaunay Triangulation Animation Code\n",
    "It is instructive to see an animation of Delaunay triangulation and Voronoi diagram as the points get added. Let’s dive into the code."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 1: Import required modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "#!/usr/bin/python\n",
    "import cv2\n",
    "import numpy as np\n",
    "import random"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 2: Check if a given point falls inside a rectangle\n",
    "\n",
    "This function checks if a given point falls inside a rectangle. This is the python counter-part of the OpenCV `rect.contains()` method provided in C++."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# Check if a point is inside a rectangle\n",
    "# Rect is an array of (x, y, w, h)\n",
    "def rectContains(rect, point) :\n",
    "    if point[0] < rect[0] :\n",
    "        return False\n",
    "    elif point[1] < rect[1] :\n",
    "        return False\n",
    "    elif point[0] > rect[2] :\n",
    "        return False\n",
    "    elif point[1] > rect[3] :\n",
    "        return False\n",
    "    return True"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 3: Draw a point on image\n",
    "\n",
    "The function `drawPoint` draws a point on an image using a specified color."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# Draw a point on the image\n",
    "def drawPoint(img, p, color ) :\n",
    "    cv2.circle( img, p, 2, color, -1, cv2.LINE_AA, 0 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 4: Draw Delaunay Triangles\n",
    "\n",
    "The function `drawDelaunay` draws Delaunay triangles on an image using a specified color.\n",
    "\n",
    "Points are inserted using an instance of `Subdiv2D` before calling this function. Delaunay triangles are obtained by calling the method `getTriangleList`.\n",
    "\n",
    "Once we get the triangles, we draw the lines by going over each vertex of the triangle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "# Draw delaunay triangles\n",
    "def drawDelaunay(img, subdiv, delaunayColor ) :\n",
    "\n",
    "    # Obtain the list of triangles.\n",
    "    # Each triangle is stored as vector of 6 coordinates\n",
    "    # (x0, y0, x1, y1, x2, y2)\n",
    "    triangleList = subdiv.getTriangleList();\n",
    "    size = img.shape\n",
    "    r = (0, 0, size[1], size[0])\n",
    "\n",
    "    # Will convert triangle representation to three vertices pt1, pt2, pt3\n",
    "    for t in triangleList :\n",
    "        pt1 = (t[0], t[1])\n",
    "        pt2 = (t[2], t[3])\n",
    "        pt3 = (t[4], t[5])\n",
    "\n",
    "        # Draw triangles that are completely inside the image\n",
    "        if rectContains(r, pt1) and rectContains(r, pt2) and rectContains(r, pt3) :\n",
    "            cv2.line(img, pt1, pt2, delaunayColor, 1, cv2.LINE_AA, 0)\n",
    "            cv2.line(img, pt2, pt3, delaunayColor, 1, cv2.LINE_AA, 0)\n",
    "            cv2.line(img, pt3, pt1, delaunayColor, 1, cv2.LINE_AA, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 5: Draw Voronoi Diagrams\n",
    "\n",
    "The function `drawVoronoi` draws Voronoi diagrams on an image.\n",
    "\n",
    "Points are inserted using an instance of `Subdiv2D` before calling this function. The Voronoi diagram is represented as a collection of facets (polygons) which can be obtained using the the method `getVoronoiFacetList`.\n",
    "\n",
    "The facets so obtained are drawn on the image using `fillConvexPoly` and the boundaries of the facets are drawn using polylines."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "# Draw voronoi diagram\n",
    "def drawVoronoi(img, subdiv) :\n",
    "\n",
    "    # Get facets and centers\n",
    "    ( facets, centers) = subdiv.getVoronoiFacetList([])\n",
    "\n",
    "    for i in range(0,len(facets)) :\n",
    "        ifacetArr = []\n",
    "        for f in facets[i] :\n",
    "            ifacetArr.append(f)\n",
    "\n",
    "        # Extract ith facet\n",
    "        ifacet = np.array(ifacetArr, np.int)\n",
    "\n",
    "        # Generate random color\n",
    "        color = (random.randint(0, 255), random.randint(0, 255), random.randint(0, 255))\n",
    "\n",
    "        # Fill facet with a random color\n",
    "        cv2.fillConvexPoly(img, ifacet, color, cv2.LINE_AA, 0);\n",
    "\n",
    "        # Draw facet boundary\n",
    "        ifacets = np.array([ifacet])\n",
    "        cv2.polylines(img, ifacets, True, (0, 0, 0), 1, cv2.LINE_AA, 0)\n",
    "\n",
    "        # Draw centers.\n",
    "        cv2.circle(img, (centers[i][0], centers[i][1]), 3, (0, 0, 0), -1, cv2.LINE_AA, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 6: Find index of closest landmark point\n",
    "\n",
    "The function `findIndex` finds the index of the closest landmark point to the given vertex from the vector of points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def findIndex(points, point):\n",
    "  diff = np.array(points) - np.array(point)\n",
    "\n",
    "  # Find the distance of point from all points\n",
    "  diffNorm = np.linalg.norm(diff, 2, 1)\n",
    "\n",
    "  # Find the index with minimum distance and return it\n",
    "  return np.argmin(diffNorm)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 7: Write Delaunay Triangles to file\n",
    "\n",
    "Save output Delaunay triangles to file. The function `writeDelaunay` takes as input an instance of `Subdiv2D`, the list of input points, and saves it to disk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# write delaunay triangles to file\n",
    "def writeDelaunay( subdiv, points, outputFileName ) :\n",
    "\n",
    "    # Obtain the list of triangles.\n",
    "    # Each triangle is stored as vector of 6 coordinates\n",
    "    # (x0, y0, x1, y1, x2, y2)\n",
    "    triangleList = subdiv.getTriangleList();\n",
    "    filePointer = open(outputFileName,'w')\n",
    "    # Will convert triangle representation to three vertices pt1, pt2, pt3\n",
    "    for t in triangleList :\n",
    "        pt1 = (t[0], t[1])\n",
    "        pt2 = (t[2], t[3])\n",
    "        pt3 = (t[4], t[5])\n",
    "\n",
    "        # Find the landmark corresponding to each vertex\n",
    "        landmark1 = findIndex(points,pt1)\n",
    "        landmark2 = findIndex(points,pt2)\n",
    "        landmark3 = findIndex(points,pt3)\n",
    "\n",
    "        filePointer.writelines(\"{} {} {}\\n\".format(landmark1, landmark2, landmark3 ))\n",
    "\n",
    "    filePointer.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 8: Create a bounding box\n",
    "\n",
    "We specify the colors for drawing lines and points on the image. Then we read the image and create a bounding box named `rect`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# Define colors for drawing.\n",
    "delaunayColor = (255,255,255)\n",
    "pointsColor = (0, 0, 255)\n",
    "\n",
    "# Read in the image.\n",
    "img = cv2.imread(\"data/images/smiling-man.jpg\");\n",
    "\n",
    "# Rectangle to be used with Subdiv2D\n",
    "size = img.shape\n",
    "rect = (0, 0, size[1], size[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 9: Draw Landmark Points\n",
    "\n",
    "The image bounding box is used to initialize an instance of `Subdiv2D`. The points are subsequently added to this object for calculating and displaying Delaunay triangulation and Voronoi diagrams."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create an instance of Subdiv2D\n",
    "subdiv = cv2.Subdiv2D(rect);\n",
    "\n",
    "# Create an array of points.\n",
    "points = [];\n",
    "\n",
    "# Allocate space for voronoi Diagram\n",
    "imgVoronoi = np.zeros(img.shape, dtype = img.dtype)\n",
    "\n",
    "# Read in the points from a text file\n",
    "with open(\"data/images/smiling-man-delaunay.txt\") as file :\n",
    "    for line in file :\n",
    "        x, y = line.split()\n",
    "        points.append((int(x), int(y)))\n",
    "\n",
    "outputFileName = \"results/smiling-man-delaunay.tri\"\n",
    "\n",
    "# Draw landmark points on the image\n",
    "for p in points :\n",
    "    drawPoint(img, p, pointsColor )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 10: Draw Delaunay Triangulation and Voronoi Diagram\n",
    "\n",
    "We calculate and draw the Delaunay triangulation ( on `imgDelaunay` ) and Voronoi Diagram on ( `imgVoronoi` ) every time a new point is added thus creating an animation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Writing Delaunay triangles to results/smiling-man-delaunay.tri\n"
     ]
    }
   ],
   "source": [
    "# Insert points into subdiv\n",
    "plotPoints = []\n",
    "count = 0\n",
    "for p in points :\n",
    "    count += 1\n",
    "    writeFileName = \"delaunayAnimation/delaunayAnimation\" + str(count).zfill(4) + \".png\"\n",
    "    \n",
    "    subdiv.insert(p)\n",
    "    plotPoints.append(p)\n",
    "\n",
    "    imgDelaunay = img.copy()\n",
    "\n",
    "    # Draw delaunay triangles and voronoi diagrams\n",
    "    drawDelaunay(imgDelaunay, subdiv, delaunayColor);\n",
    "    drawVoronoi(imgVoronoi,subdiv)\n",
    "\n",
    "    for pp in plotPoints :\n",
    "        drawPoint(imgDelaunay, pp, pointsColor)\n",
    "\n",
    "    # Create a horizontal stack\n",
    "    imgDisplay = np.hstack([imgDelaunay, imgVoronoi])\n",
    "    # Write image to file\n",
    "    cv2.imwrite(writeFileName,imgDisplay)\n",
    "\n",
    "writeDelaunay(subdiv, points, outputFileName)\n",
    "print(\"Writing Delaunay triangles to {}\".format(outputFileName))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "## Step 11: Create animation from images\n",
    "\n",
    "To create animation from saved images, we will use ImageMagick's `convert` command."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "convert: unable to open image 'delaunayAnimation/delaunayAnimation*': No such file or directory @ error/blob.c/OpenBlob/3457.\n",
      "convert: no decode delegate for this image format `' @ error/constitute.c/ReadImage/512.\n",
      "convert: no images defined `results/delaunayAnimation.gif' @ error/convert.c/ConvertImageCommand/3275.\n",
      "rm: delaunayAnimation/delaunayAnimation*: No such file or directory\n"
     ]
    }
   ],
   "source": [
    "!convert delaunayAnimation/delaunayAnimation* -loop 0 -delay 100 results/delaunayAnimation.gif\n",
    "!rm delaunayAnimation/delaunayAnimation*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# Animation\n",
    "\n",
    "![Delaunay Animation](results/delaunayAnimation.gif)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "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.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}