# Finding zebra-like pattern in image (Detection of structured-light fringe centerline from photo)

I'm working in a project where fringes are projected against a subject, and a photo is taken. The task is to find the centerlines of the fringes, which represent, mathematically, the 3D curve of intersection between the fringe plane and the subject surface.

The photo is a PNG (RGB), and former attempts used grayscaling then difference thresholding to get a black-and-white, "zebra-like" photography, from which it was easy to find the midpoint of each pixel column of each fringe. The problem is that, by thresholding and also by taking the mean height of a discrete pixel column, we're having some precision loss and quantization, which is not desired at all.

My impression, by looking at the images, is that the centerlines could be more continuous (more points) and smoother (not quantized) if they were detected directly from the non-thresholded image (either RGB or grayscale), by some statistical sweeping method (some flooding / iterative convolution, whatever).

Below is an actual sample image:

Any suggestion would be much appreciated!

• it is very interesting. But by the way, I am doing some research using color stripe to detect 3d object. Because using color stripe, it is easy to find each stripe's correspondence from the projector.so using trigonometry the 3d information can be computed. How do you find the correspondence if the color is the same ? I guess your project is also about 3d reconstruction?
– user4628
May 22 '13 at 2:08
– Peter K.
May 22 '13 at 14:03
• Sorry for one more question instead of giving answer, In phase shift method we calculate phase at each pixel in projected image, but here why we need to find out center line of fringe, might be my question is to much silly but i don't no that, so please tel me the exact reason. U can delete my question after giving answer
– user11661
Nov 18 '14 at 12:30
• These are different methods. I am modeling a series of geometric planes by projecting a series of white stripes (each one forming a "plane" in 3D space). Thus, I need to find centerline of fringes, because planes have no thickness. Sure I could perform phase-shift analysis, but there is one problem: my projection is binary (black and white stripes alternating), the intensity does not vary sinusoidally, and so I cannot perform phase shifting (and don't need to, currently). Nov 18 '14 at 12:38

I suggest the following steps:

1. Find a threshold to separate the foreground from the background.
2. For each blob in the binary image (one zebra stripe), for each x, find the weighted center (by pixel intensity) in y direction.
3. Possibly, smooth the y values, to remove noise.
4. Connect the (x,y) points by fitting some kind of curve. This article might help you. You can also fit a high-level polynomial, though it is worse, in my opinion.

Here is a Matlab code that shows steps 1,2 and 4. I skipped the automatic threshold selection. Instead I chose manual th=40:

These are the curves that are found by finding the weighted average per column:

These are the curves after fitting a polynomial:

Here is the code:

function Zebra()
im = uint8(mean(im,3));

th = 40;
imBinary = im>th;
imBinary = imclose(imBinary,strel('disk',2));
% figure;imshow(imBinary);
labels = logical(imBinary);
props =regionprops(labels,im,'Image','Area','BoundingBox');

figure(1);imshow(im .* uint8(imBinary));
figure(2);imshow(im .* uint8(imBinary));

for i=1:numel(props)
%Ignore small ones
if props(i).Area < 10
continue
end
%Find weighted centroids
boundingBox = props(i).BoundingBox;
ul = boundingBox(1:2)+0.5;
wh = boundingBox(3:4);
clipped = im( ul(2): (ul(2)+wh(2)-1), ul(1): (ul(1)+wh(1)-1) );
imClip = double(props(i).Image) .* double(clipped);
rows = transpose( 1:size(imClip,1) );
%Weighted calculation
weightedRows  = sum(bsxfun(@times, imClip, rows),1) ./ sum(imClip,1);
%Calculate x,y
x = ( 1:numel(weightedRows) ) + ul(1) - 1;
y = ( weightedRows ) + ul(2) - 1;
figure(1);
hold on;plot(x,y,'b','LineWidth',2);
try %#ok<TRYNC>
figure(2);
[xo,yo] = FitCurveByPolynom(x,y);
hold on;plot(xo,yo,'g','LineWidth',2);
end
end
end

function [xo,yo] = FitCurveByPolynom(x,y)
p = polyfit(x,y,15);
yo = polyval(p,x);
xo = x;
end

• I found this very interesting. I use Python, but anyway I'm gonna have to study the rationale of all this. As an independent comment, I tend not to perform classic image processing (directly upon quantized image containers such as uint8 arrays), but instead load everything to memory as float arrays before applying the operations. Also, I'm surprised with the results from the lower half of your image, blue lines are not running along the expected fringe midlines... (?). Thanks for now, I'm gonna bring some feedback as soon as I get some result! Oct 15 '12 at 18:04
• @heltonbiker, please check the updated answer. You are right about floating point, I used it when I converted to double. About the results in the lower half, I need to check, it might be a software bug Oct 15 '12 at 18:07
• @heltonbiker, done. It was indeed a bug related to 1 based indexing. Oct 15 '12 at 18:26
• Excelent! Amazing, indeed. With this technique, and for my purposes, the smoothing not only won't even be needed, but also would be harmful. Thanks a lot for your interest! Oct 15 '12 at 19:15

I wouldn't use the RGB image. Color images are typically made by putting a "Bayer Filter" on the camera sensor, which usually reduces the resolution you can achieve.

If you use the grayscale image, I think the steps you described (binarize "zebra" image, find midline) are a good start. As a final step, I would

• Take each point in the midline you found
• take the grayvalues of the pixels in the "zebra" line above and below
• fit a parabola to these grayvalues using least mean squares
• the apex of this parabola is an improved estimate of the midline position
• Nice thoughts. I plan to use some sort of parabola or spline along the peak values of each pixel column, but I'm still wondering if I should examine a pixel column or instead a pixel "region" along the line... Gonna wait some more for more answers. Thanks for now! Oct 15 '12 at 16:18
• @heltonbiker - as a quick test use only the green channel. There are normally 2x as many green pixels on a colour sensor and it is less interpoalted than red and blue Oct 21 '12 at 4:03
• @MartinBeckett Thanks for your interest, I have already analyzed each channel, and indeed the green one seems to be much more resolved than, say, the red one. Plotting intensity values of vertical cross sections for each channel, though, the "stripe-pattern" does not seem to change so much between channels, and I am currently mixing them equally upon conversion to grayscale. Even though, I still plan to study the best linear combination between channels to get the best-contrast result, OR to acquire images already in grayscale. Thanks again! Oct 21 '12 at 22:55

Here is yet an alternative solution to your problem by modelling your question as a 'path optimization problem'. Though it is more complicated than the simple binarization-and-then-curvefitting solution, it is more robust in practice.

From the very high level, we should consider this image as a graph, where

1. each image pixel is a node on this graph

2. each node is connected to some other nodes, known as neighbours, and this connection definition is often referred as to the topology of this graph.

3. each node has a weight (feature, cost, energy, or whatever you want to call it ), reflecting the likelihood that this node is in an optimal central-line we are looking for.

As long as we can model this likelihood, then your problem of finding 'the centerlines of the fringes' becomes to the problem to find local optimal paths on the graph, which can be effectively solved by dynamic programming, e.g. Viterbi algorithm.

Here are some pros of adopting this approach:

1. all your results will be continuous ( unlike the threshold method that might break one center line into pieces )

2. a lot of freedoms to construct such a graph, you can select different features, and graph topology.

3. your results are optimal in the sense of path optimizations

4. your solution will be more robust against noise, because as long as noise is equally distributed among all pixels, those optimal paths remain stable.

Here is a short demonstration of the above idea. Since I donot use any prior knowledge to specify what are possible starting and ending nodes, I simply decode w.r.t every possible starting node.

For the fuzzy endings, it is caused by the fact that we are looking for optimal paths for every possible ending nodes. As a result, though for some nodes located in dark areas, the highlighted path is still its local optimal one.

For the fuzzy path, you could either smooth it after you find it or use some smoothed features instead of raw intensity.

It is possible to restore partial paths by changing starting and ending nodes.

It will not be difficult to prune these undesired local optimal paths. Because we have the likelihoods of all paths after viterbi decoding, and you may use various prior knowledge ( e.g. we see it is true that we only need one optimal path for those sharing the same source. )

For more details, you may refer to the paper.

 Wu, Y.; Zha, S.; Cao, H.; Liu, D., & Natarajan, P.  (2014, February). A Markov Chain Line Segmentation Method for Text Recognition. In IS&T/SPIE 26th Annual Symposium on Electronic Imaging (DRR), pp. 90210C-90210C.


Here is a short piece of python code using to make the above graph.


import cv2
import numpy as np
from matplotlib import pyplot
image_path = ;
img = cv2.imread( image_path, 0 );
rgb = cv2.imread( image_path, -1 );

# some feature to reflect how likely a node is in an optimal path
img = cv2.equalizeHist( img ); # equalization
img = img - img.mean(); # substract DC
img_pmax = img.max(); # get brightest intensity
img_nmin = img.min(); # get darkest intensity
# express our preknowledge
img[ img > 0 ] *= +1.0  / img_pmax;
img[ img = 1 :
prev_idx = vt_path[ -1 ].astype('int');
vt_path.append( path_buffer[ prev_idx, time ] );
time -= 1;
vt_path.reverse();
vt_path = np.asarray( vt_path ).T;

# plot found optimal paths for every 7 of them
pyplot.imshow( rgb, 'jet' ),
for row in range( 0, h, 7 ) :
pyplot.hold(True), pyplot.plot( vt_path[row,:], c=np.random.rand(3,1), lw = 2 );
pyplot.xlim( ( 0, w ) );
pyplot.ylim( ( h, 0 ) );

• This is a very interesting approach. I confess the topic of "graphs" has been obscure to me until recently when (on this same project) I only could solve another problem using graphs. After I "got it", I realized how powerful these shortest paths algorithms can be. Your idea is very interesting and it is not impossible that I would reimplement for this one if I have the need / oportunity. Thank you very much. Nov 11 '14 at 11:26
• As for your current results, from my experience it would be probably better to smooth the image first with gaussian and/or median filter, before building the graph. This would give much smoother (and more correct) lines. Also, one possible trick is to expand the neighbourhood to allow for "direct jumping" over two or more pixels (up to a given limit, say, 8 or 10 pixels). Of course a suitable cost function should be chosen, but I think it's easy to tune. Nov 11 '14 at 11:29
• Oh yes. I simply picked something at hand, you can definitely use other topology and energy functions. Actually, this framework is also trainable. In particular, you start with the raw intensity, decode for optimal paths, only pick up those optimal nodes with high confidences, and in this way you get 'labelled data'. With this small part of automatically labelled data you can learn many kinds of useful things. Nov 11 '14 at 19:04

Thought I should post my answer as it is bit different from other approaches. I tried this in Matlab.

• sum all channels and create an image, so all channels are weighted equally
• perform morphological closing and Gaussian filtering on this image
• for each column of the resulting image, find the local maxima and construct an image
• find the connected components of this image

One disadvantage I see here is that this approach will not perform well for some orientations of the stripes. In that case we have to correct its orientation and apply this procedure.

Here's the Matlab code:

im = imread('m0sy7.png');
imsum = sum(im, 3); % sum all channels
h = fspecial('gaussian', 3);
im2 = imclose(imsum, ones(3)); % close
im2 = imfilter(im2, h); % smooth
% for each column, find regional max
mx = zeros(size(im2));
for c = 1:size(im2, 2)
mx(:, c) = imregionalmax(im2(:, c));
end
% find connected components
ccomp = bwlabel(mx);


For example, if you take the mid column of the image, its profile should look like this: (in blue is the profile. in green are the local maxima)

And the image containing the local maxima for all columns looks like this:

Here are the connected components (though some stripes are broken, most of them get a continuous region):

• This is actually what we are doing now, with the only difference being how to find local maxima for each pixel column: we use a parabolic interpolation to find the exact vertex of the parabola passing through the pixel with maximum value and its upper and lower neighbors. This allow s for the result to be "between" pixels, which better represents the subtle smoothness of the lines. Thanks for your answer! Nov 12 '14 at 11:33