Edge orientation histograms in global and local features (0ctave/Matlab)

Posted on by | 40 Replies | About me

This is the last type of histograms I used in my project of training an Adaboost classifier to distinguish two artistic styles.

The basic idea in this step is to build a histogram with the directions of the gradients of the edges (borders or contours). It is possible to detect edges in an image but it in this we are interest in the detection of the angles. This is possible trough Sobel operators. The next five operators could give an idea of the strength of the gradient in five particular directions (Fig 1.).

Fig. 1 The sobel masks for 5 orientations: vertical, horizontal, diagonals and non-directional

The convolution against each of this mask produce a matrix of the same size of the original image indicating the gradient (strength) of the edge in any particular direction. It is possible to count the max gradient in the final 5 matrix and use that to complete a histogram (Fig 2.)

Fig 2. Edge Orientation Histogram

In terms of avoiding the amount of non important gradients that could potentially be introduced by this methodology, an option is to just take into account the edges detected by a very robust method as the canny edge detector. This detector returns a matrix of the same size of the image with a 1 if there is an edge and 0 if there is not and edge. Basically it returns the contours of the objects inside the image. If you just consider the 1's we are just counting the most pronounced gradients.

I am also interested in calculate global and local histograms (I have already talk about this in previous posts).  For example Fig 1, Fig 2 and Fig 3 presents the regions for three different type of region divisions: 1, 3, 8 respectively.

Fig 1. 1x1 region divisions
Fig 2. 3x3 region divisions
8x8 region divisions

 

 

 

 

 

 

 

I found this code but I had to do several modifications because of my particular requirements. The most importants are:

  • I need to work just with gray scale images
  • I took out an initial filter that seems to be unnecessary
  • I need to extract histograms of different regions
  • I need a linear response. Just a vector with the responses together

I am posting the code with all the modifications:

# parameters
# - the image
# - the number of vertical and horizontal divisions
function [data] = edgeOrientationHistogram(im, r)

% define the filters for the 5 types of edges
f2 = zeros(3,3,5);
f2(:,:,1) = [1 2 1;0 0 0;-1 -2 -1];
f2(:,:,2) = [-1 0 1;-2 0 2;-1 0 1];
f2(:,:,3) = [2 2 -1;2 -1 -1; -1 -1 -1];
f2(:,:,4) = [-1 2 2; -1 -1 2; -1 -1 -1];
f2(:,:,5) = [-1 0 1;0 0 0;1 0 -1];

% the size of the image
ys = size(im,1);
xs = size(im,2);

# The image has to be in gray scale (intensities)
if (isrgb(im))
    im = rgb2gray(im);
endif

# Build a new matrix of the same size of the image
# and 5 dimensions to save the gradients
im2 = zeros(ys,xs,5);

# iterate over the posible directions
for i = 1:5
    # apply the sobel mask
    im2(:,:,i) = filter2(f2(:,:,i), im);
end

# calculate the max sobel gradient
[mmax, maxp] = max(im2,[],3);
# save just the index (type) of the orientation
# and ignore the value of the gradient
im2 = maxp;

# detect the edges using the default Octave parameters
ime = edge(im, 'canny');

# multiply against the types of orientations detected
# by the Sobel masks
im2 = im2.*ime;

# produce a structur to save all the bins of the
# histogram of each region
eoh = zeros(r,r,6);
# for each region
for j = 1:r
    for i = 1:r
        # extract the subimage
        clip = im2(round((j-1)*ys/r+1):round(j*ys/r),round((i-1)*xs/r+1):round(i*xs/r));
        # calculate the histogram for the region
        eoh(j,i,:) = (hist(makelinear(clip), 0:5)*100)/numel(clip);
    end
end

# take out the zeros
eoh = eoh(:,:,2:6);

# represent all the histograms on one vector
data = zeros(1,numel(eoh));
data(:) = eoh(:);

 
The makelinear function doesn't exist in Octave. All it does is converting a matrix into a vector. If the function doesn't exist you can use the following function.

# makelinear.m
# converts any input matrix into a 1D vector (output)

function data = makelinear(im)
data = zeros(numel(im),1);
data(:) = im(:);

40 thoughts on “Edge orientation histograms in global and local features (0ctave/Matlab)”

    • My bad, this blog was full of spamming comments and until today I had time to check
      and configure the spam detector.

      Makelinear is a simple function that take a matrix and convert it into an array. I had to search into my files but if you still need help, I would dig into them and try to get you the code.

      Sorry about the ultra late reply

  2. Hello [Your Name :)],

    I have two questions regarding this code of yours at
    http://robertour.com/2012/01/26/edge-orientation-histograms-in-global-and-local-features/

    1. Why the outcome of the code is in float, I think the histogram information must be in whole numbers; that is, how many pixels are there for which gray intensity.
    2. The number of vertical and horizontal divisions are same (‘r’), this may be modified, but when use it with say r=4 and image size of 144×176 the output data is of size 80 ( i am not able to understand why? and what this 80 is).

    Could you please help me find this out.

    Regards
    Vinay

    • Sorry about this late reply. I had a spam problem and I couldn’t fix it until today.
      I struggle to remember but I think I got the answers in case some body else is wondering the same.
      1. They are floats because I need to normalize. I do so to compare to other images of other sizes. I normalize dividing by numel(clip). numel(clip) is the number of pixel of the region.
      2. There are 4×4=16 regions, and the histogram for each region is of size 5: horizontal, vertical, diagonal left-right , diagonal right-left / and undirectional. Sometimes it is just impossible to know the direction of the edge. So, 5×16=80. All is represented in just one vector.

      Regards,
      My name = Roberto

    • I don’t understand what are u asking for because the sobel masks are in the code. They are just matrices.

      % define the filters for the 5 types of edges
      f2 = zeros(3,3,5);
      f2(:,:,1) = [1 2 1;0 0 0;-1 -2 -1];
      f2(:,:,2) = [-1 0 1;-2 0 2;-1 0 1];
      f2(:,:,3) = [2 2 -1;2 -1 -1; -1 -1 -1];
      f2(:,:,4) = [-1 2 2; -1 -1 2; -1 -1 -1];
      f2(:,:,5) = [-1 0 1;0 0 0;1 0 -1];

  4. Dear Roberto

    I find your code very helpful for my thesis, thanks for your efforts and help, but i want to know that… The max value is chosen as it is considered to be the representative of of edge type for a particular block. Here I want to know if two or more values are same and max…
    like [0, 0, -1, 0, 0]. In this case, the block as four edges, what has to be done. I mean which one has to be chosen…

    Regards

  5. That is going to depend on the function

    # calculate the max sobel gradient
    [mmax, maxp] = max(im2,[],3);

    The Octace documentation doesn’t specify, what is the strategy in that case:
    http://octave.sourceforge.net/octave/function/max.html

    It depends on the internal algorithm that octave uses. I could be the first it found. That will give a small bias in favour of the first direction of the array weight to array. I would adventure is meaningless for a few reasons:

    (1) The possibility of a tie on the max value is quite low because we are looking for a gradient.
    (2) If the tie exists, it would most probably be on 2 consecutive gradient directions. Maybe is possible to mathematically show that this is the only possible scenario for a tie (except maybe the last column, see (5))
    (3) If the octave algorithm is consistent (always return the first max), then we are being “fair” to everybody and shouldn’t matter at all because the bias is for all of them. It is a bias in the clockwise direction.
    (4) In any case if the there is a tie in two very different directions, then canny detector will almost certainly throw a 0, because it is not an edge.
    (5) If there is a tie in the max value against the last column (no-directionality), then, canny detector would throw a 0.

    If you are still unsure, I seriously recommend to ask in stackoverflow.com about the max function in Octave. I dig a bit in Google without luck. You can event try to ask the developers. Or do some tests to try to deduce the behaviour or max.

    For the reasons I give, I wouldn’t be so worried, specially if you are urged.

  7. They are consider one of the most informative feature, but the disadvantage is that there are expensive to calculate. You can afford them if you don’t mind a few minutes/hours/days depending on your task. However, they became useless for real time tasks, such as face/people/object detection on cameras.

  8. Hi
    i am working on image thresholding and i want to use the gradient histogram
    i am searching for a method to smooth this kind of histogram so that noisiy pixels and edge pixels are isolated
    thanks

  9. I haven’t work with noise reduction directly. This post was related to a machine learning algorithm (adaboost). Specific noise reduction didn’t make much sense for this project because I had a big set of images and I wouldn’t be sure in which of them the noise reduction was successful (unless I would check them manually)

    That said, the canny edge detector has a generic noise reduction stage. It is basically a Gaussian, it takes the surrounding (say 3 neighbours) pixels of each pixel in the image, and multiply by a Gaussian (represented by a matrix).

    I really cannot help you very much. I don’t know many others methods of noise reduction, but I can point you to the wikipedia: http://en.wikipedia.org/wiki/Noise_reduction#Removal

    Hope I helped a bit

    • It is actually very simple. Imagine you have any image and you use an algorithms to detect borders such as this .

      Now, you want to build an histogram that shows the distribution of the direction of the edges (border). Basically, how many vertical, horizontal or diagonals you have? (see the Fig 2 in this post). The algorithm count such borders. You can have as many directions as you want.

    • Probably. A histogram represents the distribution of data. Some of the algorithms around have patents so it might refer to a particular one.

      There are many ways of building an edge histogram. In this example, I just took 4 angles but you could potentially take as many as you want.

      In practice, 4 should be enough for most of the problems, but as everything else in computer vision, the missing clue for a problem could be anywhere.

      I will try to answer your other questions.

    • Hi Mabelle,

      The zeros indicate that no maximum direction (up, down,left, right or none) was found. In other words there was a tie, so they don’t really carry a lot of information, but noise… although I might be wrong and it could help.

      Cheers!

  14. Hello,

    I’m developing a software with OpenCv Library , my software is based on a paper about the multispectral stereo odometry. After edge extraction, from the image , i have to compute for each extracted feature the edge histogram descriptor.

    From the paper:

    “The spatial component of the descriptor is obtained as follows:

    • Select a region of P × P pixels centered at the keypoint of interest from the edge map.

    • Divide the region into 16 (4 × 4) subregions.

    • Compute local edge histograms for each subregion where five bins are used to categorize an edge: horizontal, vertical, 45◦ diagonal, 135◦ diagonal, and isotropic (no orientation).

    The resulting histogram vector formed of 80 bins (4 × 4 × 5) is then normalized. ”

    My question is:
    if the keypoint is on the border of the edgeMap, my window (PxP) is out of bound for every P>=3 . Some paper reccomended to use P=100. The same problem occurs also if the keypoint is near to the border.

    How should i proceed?

  15. Hi Pasquale,

    Thanks for your question. I am not really sure if I would be able to help you with this, but I will try.

    The decision is, in principle, arbitrary; and the answer depends on you intuition and knowledge of what you are doing. I can tell you a few general things:

    1. You can use as parameter the pixels of the region r (r = P/4) instead of P, so you garantee that you are working with P that can be divided by 4.

    2. However, in general if a region has 2 or 3 pixel less is not really a problem in terms of information. It could have been a problem in terms of the size of the resulting vector but since you are always using 80 bins, you are good to go.

    3. The border problem is slightly more complicated. I would try to run away from it by not accepting any keypoint that is less than P/2 of any border. If this doesn’t make sense, then

    It really depends what are you using this histograms for (I know you said multispectral stereo odometry but it seems it would take me a bit more time to understand what it is about). There is not really that much you can do either, are you sure those keypoints are important? Here, some hints:

    4. If you are detecting objects and assuming that the objects are in the center of the PxP area then you should do your best to keep things consistent, i.e. the 80 features (bins) should have a meaning. For example, the subregions that are outside or have less than say 50% of the pixels are simply discarded. In the histogram you indicate the missing values with -1. Hopefully, the learning algorithm is good enough to understand that those images that have -1s are “special” and give them a different treatment. if the learning algorithm doesn’t accept negatives, then try using 2, since the rest are normalized (from 0 to 1), then it should be ok. Actually, since you have the isometric bar, you could also use 0s on all the 5 bars, which is a rare case (that just noise conditions might produce).

    5. The other option is scaling P according to the position, so there are enough pixels for the division. Since you are normalizing the histograms, it is comparable. However, scaling is problematic, and I wouldn’t recommended unless you are very sure scaling makes sense in the problem domain. I cannot think any case, but you never know. Basically, you loose resolution and you are not comparing the same thing.

    I hope this makes sense, don’t hesitate in asking more. Also, this blog is part of a series, just in case it might be useful.

    Cheers

  19. I am doing my project on cbir .. but I am unable to find feature vector from sobel operator. . means I got binary image from sobel operator but not fearture vectors plz plz plz help me out .

  25. Sorry, I lost the your site address.
    I was looking for and found today
    I’ll test here

    Thank you for your help
    Thank you for your help
    Thank you for your help
    Thank you for your help

  27. I am confused between edge distribution histogram features and Edge orientation histograms. Both are same or different? In many papers 75-dimensional edge distribution histogram features are calculated.

    • Take this answer more as an opinion based on what I have read, bottom line is: the names are confusing and I don’t think there is consensus. My suggestion is to not think about the names so strictly, and cite the paper that more or less correspond to your implementation.

      The definitions are vague, or non existent. I think a lot of authors collide with similar names at more or less the same time – and, well, there are quite a few ways of estimating the edge orientations, and also of binning things.

      For the features part in “edge distribution histogram features”, it would refer to each bin of the histogram, it would be similar to say “edge orientation histograms features”. If you have 5 bins, well, you have five features once you input the histogram into a machine learning algorithm.

      It is also possible to calculate as many dimensional edges as directions there are in a 2d plane (I mean, infinite). However, you sobel mask needs to be created accordingly; you would need sobel masks bigger than 3×3 matrices to create the rotations. It gets ugly and computationally expensive. There are alternatives to Sobel masks to calculate rotations, e.g. using gaussians that might be more efficient.

      Here is also attempts to answers from other people.
      https://stackoverflow.com/questions/11505234/histogram-of-oriented-gradients-vs-edge-orientation-histograms#12745508

Leave a reply