Raspberry Pi
Edge Detection
Contents |
0 Edge Detection
The algorithm we will look at in this tutorial is an edge detection algorithm, specifically an edge detection algorithm based on the Sobel operator. This algorithm works by calculating the gradient of the intensity of the image at each point, finding the direction of the change from light to dark and the magnitude of the change. This magnitude corresponds to how sharp the edge is.
1 The Sobel Operator
To calculate the gradient of each point in the image, the image is convolved with the Sobel Kernel. Convolution is done by moving the kernel across the image, one pixel at a time. At each pixel, the pixel and its neighbours are weighted by the corresponding value in the kernel, and summed to produce a new value. This operation is shown in the following diagram.
~~~~~~~~~~~~~~~~~~~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
| k00 | k10 | k20 | | i00 | i10 | i20 | i30 | ..
|~~~~~+~~~~~+~~~~~+ |~~~~~+~~~~~+~~~~~+~~~~~+~~~
| k01 | k11 | k21 | | i01 | i11 | i21 | i31 | ..
|~~~~~+~~~~~+~~~~~+ |~~~~~+~~~~~+~~~~~+~~~~~+~~~
| k02 | k12 | k22 | | i02 | i12 | i22 | i32 | ..
|~~~~~+~~~~~+~~~~~+ |~~~~~+~~~~~+~~~~~+~~~~~+~~~
| i03 | i13 | i23 | i33 | ..
Applying the kernel to the |~~~~~+~~~~~+~~~~~+~~~~~+~~~
to the first pixel gives us | .. | .. | .. | .. | ..
an output value for that
pixel
o00 = k11 * i00 + k12 * i01 + k21 * i10 + k22 * i11
we treat it as though the kernel has been overlaid onto the image, with the centre pixel of the kernel aligning we the first pixel in the image. Then we multiple each entry in the kernel by the value beneath it, and sum them to produce a single single output value from that pixel.
For the pixels on the boundary, we just ignore any entries in the kernel that fall outside.
Edge detection using the Sobel Operator applies two separate kernels to calculate the x and y gradients in the image. The length of this gradient is then calculated and normalised to produce a single intensity approximately equal to the sharpness of the edge at that position.
The kernels used for Sobel Edge Detection are shown below.
~~~~~~~~~~~~~
|-1 | 0 | 1 |
|~~~+~~~+~~~|
|-2 | 0 | 2 |
|~~~+~~~+~~~|
|-1 | 0 | 1 |
~~~~~~~~~~~~~
x gradient kernel
~~~~~~~~~~~~~
|-1 |-2 |-1 |
|~~~+~~~+~~~|
| 0 | 0 | 0 |
|~~~+~~~+~~~|
| 1 | 2 | 1 |
~~~~~~~~~~~~~
y gradient kernel
2 The Algorithm
Now the algorithm can be broken down into its constituent steps
1 Iterate over every pixel in the image
2 Apply the x gradient kernel
3 Apply the y gradient kernel
4 Find the length of the gradient using pythagoras' theorem
5 Normalise the gradient length to the range 0-255
6 Set the pixels to the new values
This illustrates the main steps, though we miss out some specifics such as what we do when we meet the boundary of the image. Now we will try implementing this in Python, and we will fill in any gaps in the algorithm as we go along
Create a new file for the python script, you can call it what you like. We will use the code from the first tutorial as a template for this one.
import imgproc
from imgproc import *
# import the maths module
import math
# open a webcam to take pictures, or load a sample image
camera = Camera(160, 120)
# Open a viewer window to display images
viewer = Viewer(160, 120, "Edge Detection")
# take a picture from the camera
img = camera.grabImage()
# display the image in the viewer
viewer.displayImage(img)
# iterate over each pixel in the image
for x in range(0, img.width):
for y in range(0, img.height):
# our code will go here
# display the image again
viewer.displayImage(img)
# delay for 5000 milliseconds to give us time to see the changes, then exit
waitTime(5000)
We need to do one tweak so that this will work for our edge detection. When we are iterating over the image, we will be changing value however when we move on to the next pixel, the just changed value will still be in the area of influence of our kernel. For this reason we need to make sure the image we read pixels from and the image we write pixels to are different, and we can do this by just created an empty image, which we will populate with the new values.
img_edge = Image(img.width, img.height)
Now we will read pixels from the original image, but write the edge values to
img_edge
Before we start fiddling with the pixel values, we need to decide how we will
handle the image boundaries. There are a few different ways to do this, one
way is to check each pixel to see if it is on the boundary of the image, and
ignore the non-existant pixels on that side. An alternative way would be to
wrap the pixels, so that a non-existant pixel on the left edge actually maps
to the pixel on the same row on the right edge. For simplicities sake I am
simply going to ignore the boundary pixels, instead choosing to start one pixel
in from the starting point, and 1 pixel away from the finishing edge. To do
this, change the iterators to start at 1, and finish at
img.width - 1
and img.height - 1
respectively. Now
we can safely index adjacent pixels without falling off the edge of the image.
Right, lets move on to applying the gradients. For the x gradient, we need the column of pixels to the left and right, for each of the 6 pixels we index we will add up the red, green and blue intensities, multiply it by the kernel, and accumulated it into Gx. We will also do the same for the y gradient. As Gx and Gy both use the same pixels on occasions, it is best to apply the two kernels simultaneously.
# initialise Gx to 0 and Gy to 0 for every pixel
Gx = 0
Gy = 0
# top left pixel
r, g, b = img[x - 1, y - 1]
# intensity ranges from 0 to 765 (255 * 3)
intensity = r + g + b
# accumulate the value into Gx, and Gy
Gx += -intensity
Gy += -intensity
# now we do the same for the remaining pixels, left to right, top to bottom
# remaining left column
r, g, b = img[x - 1, y]
Gx += -2 * (r + g + b)
r, g, b = img[x - 1, y + 1]
Gx += -(r + g + b)
Gy += (r + g + b)
# middle pixels
r, g, b = img[x, y - 1]
Gy += -2 * (r + g + b)
r, g, b = img[x, y + 1]
Gy += 2 * (r + g + b)
# right column
r, g, b = img[x + 1, y - 1]
Gx += (r + g + b)
Gy += -(r + g + b)
r, g, b = img[x + 1, y]
Gx += 2 * (r + g + b)
r, g, b= img[x + 1, y + 1]
Gx += (r + g + b)
Gy += (r + g + b)
Now that the kernels are applied, Gx and Gy contain un-normalised values, which equal the relative length of the gradient in their respective axis. We wish to calculate the length of the gradient, and normalise it to a range suitable for displaying.
Lets calculate the relative length of each gradient. We can do this with pythagoras' theorem:
Pythagoras' theorem
sqr(a) + sqr(b) = sqr(c)
c = sqrt( sqr(a) + sqr(b) )
And in Python:
length = math.sqrt((Gx * Gx) + (Gy * Gy))
This gives us an un-normalised length, to get the length normalised, we need to
know the range of length
We know that each pixel's intensity value has a range of 0 to 765, due to it being the addition of 3 variables with range 0 to 255. Also, from looking at the kernels, we can see the maximum amount you can accumulate is +/- 4 intensities. Therefore the total range of Gx and Gy is between 0 and 3060. These values are then squared, added and square rooted, giving us a final range of 0 to 4328 ( sqrt(2 * sqr(3060)) ). We can renormalise by dividing by the upper bounds and multiplying by 255.
# normalise the length
length = length / 4328 * 255
# convert the length to an integer
length = int(length)
# set edge image to a grayscale of the normalised length
edge_img[x, y] = length, length, length
Right at the end of the code we need to change displayImage
to
display our new image rather than the original. Other than that, the script is
finished, and we can test it to see the results. The script should produce a
new grayscale image, highlighting edges from the original image.
This concludes this tutorial. Below are the code listings for this tutorial. The C code will actually continually capture images from the webcam, and you will find that on the Raspberry Pi it is capable of displaying the edge detected image in realtime.
3 Code Listing
3.1 Python Code
import imgproc
from imgproc import *
# also import the maths module for the square root function
import math
# open a webcam to take pictures
camera = Camera(160, 120)
# Open a viewer window to display images
viewer = Viewer(160, 120, "The first step")
# grab a couple of frames, just to flush the camera
camera.grabImage()
camera.grabImage()
# take a picture from the camera
img = camera.grabImage()
# create an empty image to show the edges
img_edge = Image(160, 120)
# display the image in the viewer
viewer.displayImage(img)
# iterate over each pixel in the image
for x in range(1, img.width - 1):
for y in range(1, img.height - 1):
# initialise Gx to 0 and Gy to 0 for every pixel
Gx = 0
Gy = 0
# top left pixel
r, g, b = img[x - 1, y - 1]
# intensity ranges from 0 to 765 (255 * 3)
intensity = (r + g + b)
# accumulate the value into Gx, and Gy
Gx += -intensity
Gy += -intensity
# now we do the same for the remaining pixels, left to right, top to bottom
# remaining left column
r, g, b = img[x - 1, y]
Gx += -2 * (r + g + b)
r, g, b = img[x - 1, y + 1]
Gx += -(r + g + b)
Gy += (r + g + b)
# middle pixels
r, g, b = img[x, y - 1]
Gy += -2 * (r + g + b)
r, g, b = img[x, y + 1]
Gy += 2 * (r + g + b)
# right column
r, g, b = img[x + 1, y - 1]
Gx += (r + g + b)
Gy += -(r + g + b)
r, g, b = img[x + 1, y]
Gx += 2 * (r + g + b)
r, g, b= img[x + 1, y + 1]
Gx += (r + g + b)
Gy += (r + g + b)
# calculate the length of the gradient
length = math.sqrt((Gx * Gx) + (Gy * Gy))
# normalise the length of gradient to the range 0 to 255
length = length / 4328 * 255
# convert the length to an integer
length = int(length)
# draw the length in the edge image
img_edge[x, y] = length, length, length
# display the edge image
viewer.displayImage(img_edge)
# delay for 5000 milliseconds to give us time to see the changes, then exit
waitTime(5000)
# end of the script
3.2 C Code
To compile the C code type this command in the terminal.
sudo gcc -O2 -o edge edge.c -limgproc
Now here is the C source
#include <stdlib.h>
#include <stdio.h>
#include "imgproc.h"
#include <math.h>
int main(int argc, char * argv[])
{
// initialise the library
init_imgproc();
// open the webcam, with a capture resolution of width 320 and height 240
Camera * cam = camOpen(320, 240);
// create a new viewer of the same resolution with a caption
Viewer * view = viewOpen(320, 240, "First");
while(1){
// capture an image from the webcam
Image * img = camGrabImage(cam);
// Create an empty image
Image * img_edge = imgNew(320, 240);
// iterators
unsigned int x, y;
for(x = 1; x < img->width - 1; x++){
for(y = 1; y < img->height - 1; y++){
// initialise Gx and Gy to 0
int Gx = 0;
int Gy = 0;
unsigned int intensity = 0;
char * pixel;
// Left column
pixel = imgGetPixel(img, x - 1, y - 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += -intensity;
Gy += -intensity;
pixel = imgGetPixel(img, x - 1, y);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += -2 * intensity;
pixel = imgGetPixel(img, x - 1, y + 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += -intensity;
Gy += +intensity;
// middle column
pixel = imgGetPixel(img, x, y - 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gy += -2 * intensity;
pixel = imgGetPixel(img, x, y + 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gy += +2 * intensity;
// right column
pixel = imgGetPixel(img, x + 1, y - 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += +intensity;
Gy += -intensity;
pixel = imgGetPixel(img, x + 1, y);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += +2 * intensity;
pixel = imgGetPixel(img, x + 1, y + 1);
intensity = pixel[0] + pixel[1] + pixel[2];
Gx += +intensity;
Gy += +intensity;
// calculate the gradient length
unsigned int length = (unsigned int)sqrt( (float)(Gx * Gx) + (float)(Gy * Gy) );
// normalise the length to 0 to 255
length = length / 17;
// draw the pixel on the edge image
imgSetPixel(img_edge, x, y, length, length, length);
}
}
// display the edge image
viewDisplayImage(view, img_edge);
// now we will free the memory for the various objects
imgDestroy(img_edge);
imgDestroy(img);
}
// finally quit
viewClose(view);
camClose(cam);
// finally we uninitialise the library
quit_imgproc();
return 0;
}