Monday, March 4, 2013

Sobel and Prewitt edge detector in C++: Image Processing

Sobel and Prewitt are used extensively for detecting edges in image processing.

Sobel Operator

The operator calculates the gradient of the image intensity at each point, giving the direction of the largest possible increase from light to dark and the rate of change in that direction. The result therefore shows how "abruptly" or "smoothly" the image changes at that point, and therefore how likely it is that part of the image represents an edge, as well as how that edge is likely to be oriented. The Sobel kernels are given by

Here the kernel hx is sensitive to changes in the x direction, i.e., edges that
run vertically, or have a vertical component. Similarly, the kernel hy is sensitive to changes in y direction, i.e., edges that run horizontally, or have a horizontal component. But how do we combine the results of convolution with these two kernels to give a single measure of the presence of an edge?
. The two gradients computed at each pixel (Gx and Gy) by convolving with above two kernels can be regarded as the x and y components of gradient vector. This vector is oriented along the direction of change, normal to the direction in which the edge runs. Gradient magnitude and direction are given by

Note: To reduce the complexity of above square root calculation, gradient magnitude is some times approximated by absolute summation i.e G = abs(Gx) + abs(Gy).
Prewitt Operator
Prewitt operator also works same as Sobel operator. But the kernel data is little different
Other operations are same as Sobel.

Source Code: C++
The following source code is written for Sobel but is similar for Prewitt.


using namespace std;
using namespace cv;

// Computes the x component of the gradient vector
// at a given point in a image.
// returns gradient in the x direction
int xGradient(Mat image, int x, int y)
    return<uchar>(y-1, x-1) +
                2*<uchar>(y, x-1) +
       <uchar>(y+1, x-1) -
        <uchar>(y-1, x+1) -
                   2*<uchar>(y, x+1) -
          <uchar>(y+1, x+1);

// Computes the y component of the gradient vector
// at a given point in a image
// returns gradient in the y direction

int yGradient(Mat image, int x, int y)
    return<uchar>(y-1, x-1) +
                2*<uchar>(y-1, x) +
       <uchar>(y-1, x+1) -
        <uchar>(y+1, x-1) -
                   2*<uchar>(y+1, x) -
          <uchar>(y+1, x+1);

int main()

      Mat src, dst;
      int gx, gy, sum;

      // Load an image
      src = imread("lena.jpg", CV_LOAD_IMAGE_GRAYSCALE);
      dst = src.clone();
      if( ! )
      { return -1; }

        for(int y = 0; y < src.rows; y++)
            for(int x = 0; x < src.cols; x++)
      <uchar>(y,x) = 0.0;

        for(int y = 1; y < src.rows - 1; y++){
            for(int x = 1; x < src.cols - 1; x++){
                gx = xGradient(src, x, y);
                gy = yGradient(src, x, y);
                sum = abs(gx) + abs(gy);
                sum = sum > 255 ? 255:sum;
                sum = sum < 0 ? 0 : sum;
      <uchar>(y,x) = sum;

        imshow("final", dst);

        imshow("initial", src);


    return 0;

Original Image
After applying Sobel Operator


  1. I am afraid that this line is wrong: sum = abs(gx) + abs(gy);

    You should have written:
    sum = sqrt(gx*gx + gy*gy);

    Subadditivity for absolute value (Real numbers) looks like:
    |a+b| ≤ |a| + |b|

    1. |a| + |b| is more then twice faster, and give an approximate result

  2. Zbigniew you are 100% right. The resultant of two orthogonal vectors is sqrt(gx*gx + gy*gy). But this is computationally expensive operation, above program calculates gradient of each pixel and calculate the resultant value. So to reduce the complexity in very large amount many researchers and authors have suggested that summing absolute value instead of calculating square root is a good approximation. Many books including book by Gonzalez and woods have written about this approach. See the following link for more detail

  3. How can I find the local maximum and local minimum for edges. So that I can calculate the width of the edge from the difference of local maximum and minimum.