How to obtain the chain code
In the previous post I discussed simple techniques to estimate the boundary length of a binarized object. These techniques are based on the chain code. This post will detail how to obtain such a chain code. The algorithm is quite simple, but might not be trivial to understand. Future posts will discuss other measures that can be derived from such a chain code.
In short, the chain code is a way to represent a binary object by encoding only its boundary. The chain code is composed of a sequence of numbers between 0 and 7. Each number represents the transition between two consecutive boundary pixels, 0 being a step to the right, 1 a step diagonally right/up, 2 a step up, etc. In the post Measuring boundary length, I gave a little more detail about the chain code. Worth repeating here from that post is the figure containing the directions associated to each code:
The chain code thus has as many elements as there are boundary pixels. Note that the position of the object is lost, the chain code encodes the shape of the object, not its location. But we only need to remember the coordinates of the first pixel in the chain to solve that. Also note, the chain code encodes a single, solid object. If the object has two disjoint parts, or has a hole, the chain code will not be able to describe the full object.
Image to chain code
Let’s assume we have a binary image with a single object in it. Using DIPimage, we could generate a suitable test image thus:
img = rr<30;
The algorithm also needs an array with the definition of each of the chain codes. Call this a “cipher” if you will. Indexing in this array with one code from the chain gives the change in coordinates as you go from one pixel to the next:
directions = [ 1, 0 1,-1 0,-1 -1,-1 -1, 0 -1, 1 0, 1 1, 1];
cc is one code from the chain, yields the
[x,y] increment for the coordinates.
+1 is necessary because MATLAB indexing starts at 1, whereas the first code is 0.
The algorithm starts at any pixel that is on the object’s boundary. We will use MATLAB’s function
find to find the
first “on” pixel in the image:
indx = find(dip_array(img),1)-1;
We subtract one from the index to obtain 0-based indices, which are much easier to use than MATLAB’s 1-based indices. Next we convert this index into image coordinates:
sz = imsize(img); start = [floor(indx/sz(2)),0]; start(2) = indx-(start(1)*sz(2));
MATLAB is column-major, meaning that indices increase downward first. The function
find returned the first object
pixel using this column-major indexing, meaning that none of the columns to the left have any object pixels, and in this
column, none of the pixels above have any object pixels. The boundary possibly continues downward on one side, and
up/right on the other side. That is, from this point the possible steps are chain codes 0, 1, 6 and 7. Codes 3, 4 and 5
all point to the column to the left, which we know is empty, and code 2 points up, where we also know the object cannot
be. We choose to follow the boundary clockwise (though it is equally valid to define the chain code counter-clockwise if
one so wishes), and thus we need to try the possible steps in this order: 1, 0, 7, 6. If the step 1 yields an object
pixel, this is the pixel that continues the boundary in the clockwise direction. If this pixel is empty, we try to see
if we find the object using step 0, etc. Once the next pixel is found, this step is repeated again and again until we
have travelled all the way around the object and return to the initial pixel at coordinates
start. However, each next
step has different limits for which directions are possible. The first direction to try (the “most clockwise” neighbour,
if you will) is the one that is two codes up from the previous step. That is, the boundary makes a 90° left turn. It
cannot make a 135° left turn, because that pixel is also a neighbour to the previous pixel, and if the object is there,
the previous step would have pointed to it instead. As before, we first try the “most clockwise” step first, then
continue down the codes until we find an object pixel. The following code implements this algorithm:
cc = ; % The chain code coord = start; % Coordinates of the current pixel dir = 1; % The starting direction while 1 newcoord = coord + directions(dir+1,:); if all(newcoord>=0) && all(newcoord<sz) ... && img(newcoord(1),newcoord(2)) cc = [cc,dir]; coord = newcoord; dir = mod(dir+2,8); else dir = mod(dir-1,8); end if all(coord==start) && dir==1 % back to starting situation break; end end
The loop above combines the loop over all the neighbours of a pixel (to find the next boundary pixel) and the loop over all the boundary pixels. This makes the code a lot simpler. We have a starting coordinate and a starting direction. Inside the loop, we test the pixel pointed at by this direction. If it contains an object pixel, we add this direction as a code in the chain, set the current pixel to the newly found boundary pixel, initialize the starting direction to the 90° left turn, and continue on with the loop. If the pixel pointed at contains background (or is outside the image domain), we decrease the direction by one and continue on with the loop. The loop continues until we are back at the starting position (that is, at the same coordinates and direction).
As I said before, this algorithm is very simple, though it might be difficult to understand. I recommend you take a paper and a pen and draw out the steps of the algorithm, until you are satisfied that it always works as intended. Also note, for example, the special case of an object with a single pixel: The algorithm looks in all directions once, then meets the finishing requirements and quits. Thus, a single pixel object has 0 elements in the chain code.
Chain code to image
The inverse algorithm, walking over the codes in the chain and drawing the boundary back in the image, is much more
trivial. We first create an empty, binary image, then start at coordinates
start, and modify the coordinates by
directions(cc(ii)+1,:), in a loop over
border = newim(sz,'bin'); coords = start; for ii=1:length(cc) border(coords(1),coords(2)) = 1; coords = coords + directions(cc(ii)+1,:); end joinchannels('rgb',img,border')*255
The final command,
joinchannels, overlays the two binary images in colour, the original one in red and the newly drawn
boundary in green. Where the two overlap, the pixels appear yellow. As you can see, the whole boundary is yellow,
indicating that the pixels represented by the chain code are the boundary pixels of the object.
Just for the MATLAB vectorization aficionados out there, I’ve put together this second version of the last algorithm. See if you can figure out how it works!
stride = [sz(2);1]; step = directions*stride; indx = cumsum([start*stride;step(cc+1)]); border = newim(sz,'bin'); border(indx) = 1;