3

What is the algorithm used in OpenCV function convexityDefects() to calculate the convexity defects of a contour?

Please, describe and illustrate the high-level operation of the algorithm, along with its inputs and outputs.

Jeru Luke
  • 20,118
  • 13
  • 80
  • 87
Crazy Cat
  • 186
  • 3
  • 13
  • The inputs are described in the [documentation](https://docs.opencv.org/3.4.3/d3/dc0/group__imgproc__shape.html#gada4437098113fd8683c932e0567f47ba) -- two lists of coordinates, one defining the contour, the other defining the convex hull corresponding to that contour. The library is open source, so you can [read the implementation](https://github.com/opencv/opencv/blob/808ba552c532408bddd5fe51784cf4209296448a/modules/imgproc/src/convhull.cpp#L268) -- if there's some specific part of the algorithm you don't understand, feel free to ask a concrete question. – Dan Mašek Sep 24 '18 at 14:42
  • 1
    Im looking for the explanation for the implementation, i mean in words, not by program. – Crazy Cat Sep 24 '18 at 16:00
  • OK, hope that explains it. Nevertheless, it's valuable to read the code too, and try to understand it -- it's an excellent way to learn, an essential skill for a programmer. – Dan Mašek Sep 24 '18 at 17:44

1 Answers1

11

Based on the documentation, the input are two lists of coordinates:

  • contour defining the original contour (red on the image below)
  • convexhull defining the convex hull corresponding to that contour (blue on the image below)

Example contour and convex hull

The algorithm works in the following manner:

If the contour or the hull contain 3 or less points, then the contour is always convex, and no more processing is needed. The algorithm assures that both the contour and the hull are accessed in the same orientation.

N.B.: In further explanation I assume they are in the same orientation, and ignore the details regarding representation of the floating point depth as an integer.

Then for each pair of adjacent hull points (H[i], H[i+1]), defining one edge of the convex hull, calculate the distance from the edge for each point on the contour C[n] that lies between H[i] and H[i+1] (excluding C[n] == H[i+1]). If the distance is greater than zero, then a defect is present. When a defect is present, record i, i+1, the maximum distance and the index (n) of the contour point where the maximum located.

Distance is calculated in the following manner:

dx0 = H[i+1].x - H[i].x
dy0 = H[i+1].y - H[i].y

if (dx0 is 0) and (dy0 is 0) then
    scale = 0
else
    scale = 1 / sqrt(dx0 * dx0 + dy0 * dy0)

dx = C[n].x - H[i].x
dy = C[n].y - H[i].y

distance = abs(-dy0 * dx + dx0 * dy) * scale

It may be easier to visualize in terms of vectors:

  • C: defect vector from H[i] to C[n]
  • H: hull edge vector from H[i] to H[i+1]
  • H_rot: hull edge vector H rotated 90 degrees
  • U_rot: unit vector in direction of H_rot

H components are [dx0, dy0], so rotating 90 degrees gives [-dy0, dx0].

scale is used to find U_rot from H_rot, but because divisions are more computationally expensive than multiplications, the inverse is used as an optimization. It's also pre-calculated before the loop over C[n] to avoid recomputing each iteration.

|H| = sqrt(dx0 * dx0 + dy0 * dy0)

U_rot = H_rot / |H| = H_rot * scale

Then, a dot product between C and U_rot gives the perpendicular distance from the defect point to the hull edge, and abs() is used to get a positive magnitude in any orientation.

distance = abs(U_rot.C) = abs(-dy0 * dx + dx0 * dy) * scale


In the scenario depicted on the above image, in first iteration, the edge is defined by H[0] and H[1]. The contour points tho examine for this edge are C[0], C[1], and C[2] (since C[3] == H[1]).

First edge

There are defects at C[1] and C[2]. The defect at C[1] is the deepest, so the algorithm will record (0, 1, 1, 50).

The next edge is defined by H[1] and H[2], and corresponding contour point C[3]. No defect is present, so nothing is recorded.

The next edge is defined by H[2] and H[3], and corresponding contour point C[4]. No defect is present, so nothing is recorded.

Since C[5] == H[3], the last contour point can be ignored -- there can't be a defect there.

unlut
  • 3,525
  • 2
  • 14
  • 23
Dan Mašek
  • 17,852
  • 6
  • 57
  • 85
  • 2
    Nice illustrations! As additional verbal explanation for visualizing the `distance = abs(-dy0 * dx + dx0 * dy) * scale` calculation: it’s a dot product formed by projecting vectors from `H[i]` to each `C[n]` onto a unit vector in the hull segment direction from `H[i]` to `H[i+1]` that’s been rotated 90° (thus the `scale` and `dx0` <-> `dy0` swaps with the sign change.) – rob3c May 21 '19 at 16:55
  • @rob3c Thanks. If you want, you can [edit] the answer, and add that explanation there. – Dan Mašek May 21 '19 at 18:18