# How to Detect Gradients and Edges in Images?

I want to be able to find points in images that are the centre of a radial gradient like the one shown in the left picture below. Any ideas on how I could use a Hough transform or some other computer vision method?

Thanks

example search image:

• Great question! – Spacey Aug 1 '12 at 18:11
• Also, take a look at Roberts' Cross: (en.wikipedia.org/wiki/Roberts_Cross) as an example of a way to estimate gradients. – Spacey Aug 1 '12 at 19:25
• looks like a smaller sobel operator. I'm not sure how to use that to find a radial gradient though – waspinator Aug 1 '12 at 19:32
• @waspinator: well have you run a sobel operator on your image and looked at the output? It's like the 2D equivalent of taking the derivative of a 1D function, so it should cross 0 at the local minima or maxima? – endolith Aug 1 '12 at 19:33
• For a simple Hough-like approach that would probably work you could try this: for every pixel of the image, calculate the gradient direction and render a short line segment in the direction of the gradient starting at this pixel into an accumulator. The center points you are looking for should be the highest peaks in the accumulator (by a large margin). – koletenbert Aug 1 '12 at 21:18

I was working in opencv and trying to find the peak of a gradient generated by distance transform. I realised that using morphologic operations (erosion/dilatation) in grey-scsale images was very useful in this case. If you erode dilate a grey-scale image, any pixel is going to take the value of the lower/highest neighbour. You can therefore find peaks of intensity in gradients by subtracting grey-scale image from the same dilated/eroded image. Here is my result:

And a way to do it in OpenCV/Cpp:

#include "opencv2/imgproc/imgproc.hpp"
#include "opencv2/highgui/highgui.hpp"

int main( int argc, char** argv ){

cv::Mat objects, img ,peaks,BGR;
std::vector<std::vector<cv::Point> > contours;
/* Converts it to Grayscale*/
cv::cvtColor(BGR,img,CV_BGR2GRAY);
/* Devine where are the objects*/
cv::threshold(img,objects,0,255,cv::THRESH_BINARY);
/* In order to find the local maxima, "distance"
* is subtracted from the result of the dilatation of
* "distance". All the peaks keep the save value */
cv::dilate(img,peaks,cv::Mat(),cv::Point(-1,-1),3);
cv::dilate(objects,objects,cv::Mat(),cv::Point(-1,-1),3);

/* Now all the peaks should be exactely 0*/
peaks=peaks-img;

/* And the non-peaks 255*/
cv::threshold(peaks,peaks,0,255,cv::THRESH_BINARY);
peaks.convertTo(peaks,CV_8U);

/* Only the zero values of "peaks" that are non-zero
* in "objects" are the real peaks*/
cv::bitwise_xor(peaks,objects,peaks);

/* The peaks that are distant from less than
* 2 pixels are merged by dilatation */
cv::dilate(peaks,peaks,cv::Mat(),cv::Point(-1,-1),1);

/* In order to map the peaks, findContours() is used.
* The results are stored in "contours" */
cv::findContours(peaks, contours, CV_RETR_CCOMP, CV_CHAIN_APPROX_SIMPLE);
/* just draw them and save the image */
cv::drawContours(BGR,contours,-1,cv::Scalar(255,0,0),-1);
cv::imwrite("result.png",BGR);

return 1;
}


Here is what I have so far. The way I'm populating my Hough space is far from optimal. I'm pretty sure there's some vectorization I can do to make it faster. I'm using Matlab R2011a. Original Image

Suggestions are appreciated, Thank you.

clear all; clc; close all;

[rows, columns] = size(img);
[x y] = meshgrid(1:columns, 1:rows);
u = dx;
v = dy;
imshow(img);
hold on
quiver(x, y, u, v)

%% create Hough space and populate
hough_space = zeros(size(img));

for i = 1:columns
for j = 1:rows

X1 = i;
Y1 = j;
X2 = round(i + dx(j,i));
Y2 = round(j + dy(j,i));
increment = 1;

slope = (Y2 - Y1) / (X2 - X1);
y_intercept = Y1 - slope * X1;

X3 = X1 + 5;

if X3 < columns && X3 > 1
Y3 = slope * X3 + y_intercept;
if Y3 < rows && Y3 > 1
hough_space = func_Drawline(hough_space, Y1, X1, floor(Y3), floor(X3), increment);
end
end
end
end

imtool(hough_space)


I modified a draw line function I found on matlab central to increment by a pixel by a value instead of setting a pixel to a value

function Img = func_DrawLine(Img, X0, Y0, X1, Y1, nG)
% Connect two pixels in an image with the desired graylevel
%
% Command line
% ------------
% result = func_DrawLine(Img, X1, Y1, X2, Y2)
% input:    Img : the original image.
%           (X1, Y1), (X2, Y2) : points to connect.
%           nG : the gray level of the line.
% output:   result
%
% Note
% ----
%   Img can be anything
%   (X1, Y1), (X2, Y2) should be NOT be OUT of the Img
%
%   The computation cost of this program is around half as Cubas's [1]
%   [1] As for Cubas's code, please refer
%
% Example
% -------
% result = func_DrawLine(zeros(5, 10), 2, 1, 5, 10, 1)
% result =
%      0     0     0     0     0     0     0     0     0     0
%      1     1     1     0     0     0     0     0     0     0
%      0     0     0     1     1     1     0     0     0     0
%      0     0     0     0     0     0     1     1     1     0
%      0     0     0     0     0     0     0     0     0     1
%
%
% Jing Tian Oct. 31 2000
% scuteejtian@hotmail.com
% This program is written in Oct.2000 during my postgraduate in
% GuangZhou, P. R. China.
% Version 1.0

Img(X0, Y0) = Img(X0, Y0) + nG;
Img(X1, Y1) = Img(X1, Y1) + nG;
if abs(X1 - X0) <= abs(Y1 - Y0)
if Y1 < Y0
k = X1; X1 = X0; X0 = k;
k = Y1; Y1 = Y0; Y0 = k;
end
if (X1 >= X0) & (Y1 >= Y0)
dy = Y1-Y0; dx = X1-X0;
p = 2*dx; n = 2*dy - 2*dx; tn = dy;
while (Y0 < Y1)
if tn >= 0
tn = tn - p;
else
tn = tn + n; X0 = X0 + 1;
end
Y0 = Y0 + 1; Img(X0, Y0) = Img(X0, Y0) + nG;
end
else
dy = Y1 - Y0; dx = X1 - X0;
p = -2*dx; n = 2*dy + 2*dx; tn = dy;
while (Y0 <= Y1)
if tn >= 0
tn = tn - p;
else
tn = tn + n; X0 = X0 - 1;
end
Y0 = Y0 + 1; Img(X0, Y0) = Img(X0, Y0) + nG;
end
end
else if X1 < X0
k = X1; X1 = X0; X0 = k;
k = Y1; Y1 = Y0; Y0 = k;
end
if (X1 >= X0) & (Y1 >= Y0)
dy = Y1 - Y0; dx = X1 - X0;
p = 2*dy; n = 2*dx-2*dy; tn = dx;
while (X0 < X1)
if tn >= 0
tn = tn - p;
else
tn = tn + n; Y0 = Y0 + 1;
end
X0 = X0 + 1; Img(X0, Y0) = Img(X0, Y0) + nG;
end
else
dy = Y1 - Y0; dx = X1 - X0;
p = -2*dy; n = 2*dy + 2*dx; tn = dx;
while (X0 < X1)
if tn >= 0
tn = tn - p;
else
tn = tn + n; Y0 = Y0 - 1;
end
X0 = X0 + 1; Img(X0, Y0) = Img(X0, Y0) + nG;
end
end
end

• I think I will attribute the bounty to your answer, as nobody else bothered to contribute. It's not exactly what I want but it's the closest of the 3. Did you further improve this method? – Cape Code Jan 21 '14 at 20:46

Run a Histogram of Oriented Gradients over patches of the image - the peak in each of those histograms will give you the dominant direction of that patch (like the arrows you show).

Find where all those arrows intersect - if that point is inside the object it could be the centre of a radial gradient.