Algorithm for shape registration : Iterative Closest Point

Just a few words about a work I am doing for G’MIC.

I’ve started to implement the Iterative Closest Point algorithm in G’MIC. This algorithm can be used to register two point clouds together (so basically two shapes as well) :

What I’ve done so far is just testing stuffs, so this is only the beginning. But my first tests are promising:

  • First experiment : Register a star shape (in red) and a rotated/de-zoomed star shape (in green).

(red shape is the target, green shape is the source and blue lines correspond to the closest matchs found at each iteration). The source / target shapes do not have the same number of points, of course.

  • Second experiment : Register a heart shape with a noisy heart shape. Even with the fact the source shape (green) has noisy coordinates, the algorithm converges into something correct.


The transformation here has 6 parameters and is a mix of rotation/translation/shearing.

As I said, this is only a prototype for now, but I’d like to make it available as a stdlib command in the future. I’m pretty sure we could find interesting uses of this
(@patdavid , you know what I’m talking about! :wink: ).

Feel free to play with the current prototype code if you wish :

foo :
  shape_heart 300 +erode. 3 -[-2,-1]
  1,1,1,2 eval.. "i?da_push([x,y]); end(resize(#-1,1,da_size(),1,2,0))" rm..

  shape_heart 200 rotate. 35 b. 1 ge. 50% +erode. 3 -[-2,-1] r. 300,300,1,1,0,0,0.5,0.5
  1,1,1,2 eval.. "i?da_push([x,y]); end(resize(#-1,1,da_size(),1,2,0))" rm..

  # Start iteration loop.
  repeat inf,f

  # Find closest points
  1,1,1,2 eval.. ">
    X = [i0,i1];
    dmin = inf;
    Pmin = [ 0,0 ];
    repeat (h#0,p,
      P = I[#0,p];
      d = norm(P - X);
      d<dmin?(dmin = d; Pmin = P);

  eval[0] ">I(#-1,i0,i1) = [ 255,0,0 ]"
  eval[1] ">P = [ i0,i1]; Q = I[#2,y]; polygon(#-1,2,P,Q,0.5,0,0,255)"
  eval[1] ">I(#-1,i0,i1) = [ 0,255,0 ]"
  w. 600,600

  # Solve least-square problem.
  +permute[2] cyzx y. # Y
  +permute[1] cyzx 1,100%,1,1,1 a[-2,-1] x
  l. r. 100%,200%,1,1,1 s y repeat $! r[$>] 6,100%,1,1,0,0,{$>%2} done a y endl # A
  solve.. . rm.

  # Apply estimated transformation on source coordinates.
  f[1] "
      A = i[#-1,0];
      B = i[#-1,1];
      C = i[#-1,2];
      D = i[#-1,3];
      E = i[#-1,4];
      F = i[#-1,5];
    X = i0;
    Y = i1;
    [ A*X + B*Y + C,
      D*X + E*Y + F ]"


I can’t imagine the use of this yet, but maybe some of this could help on a new filter. I was thinking of connecting points challenge which I have yet to start.

Oh man this is awesome!
I’ll fire up some tests for it soon! :heavy_heart_exclamation:

1 Like

Another great addition! It is certainly very useful for registration of point clouds and coarsely dithered images. Is this already available on 3.1.5?

No, I’ve actually put this project on hold.
I may return back to it later, especially if I can collaborate with @patdavid on a common project.

1 Like