Insert new control points over an image

Question

I'm trying to make an automatic non-rigid image registration using B-Splines as a deformation model, but I'm having some problems:

1. Initially, I put a few control points (5x5), and after a long time, the result is good, but the processing time is prohibitive.
2. I've read lots of articles in image registration (Rueckert, et.al, 1999; for example) and to improve the performance, the authors use multi-level B-Spline deformation, but the algorithms are not totally clear for me.
3. In this article, the author presents a set of formulas (section 4.2) to increase the number of control points without changing the shape of the B-Spline, but my implementation does not return what I expected (for example, I expected apply this new 10x10 control point transformation and have the same 5x5 transformation).

Here is my MATLAB implementation:

function ref_coeff = refine_coeff_grid(c, ft, delta)

%c -- Control points
%ft -- image to be transformed
%delta -- spacing between the control points

[X,Y]=ndgrid(-delta:delta:(size(ft, 1) +(delta*2)), -delta:delta:(size(ft, 2)+(delta*2)));

[r, co, s] = size(c);

i = 2:(r-2);
j = 2:(co-2);

new_c = zeros(2*(r-3)+3, 2*(co-3)+3, s);
new_c(:, :, 1) = X;
new_c(:, :, 2) = Y;

new_c(2*i, 2*j, :) = (c(i-1, j-1, :) + c(i-1, j+1, :) + c(i+1, j-1, :) + c(i+1, j+1, :) + ...
                      6*(c(i-1, j, :) + c(i, j-1, :) + c(i, j+1, :) + c(i+1, j, :)) + 36*c(i, j, :))/64;

new_c(2*i, 2*j+1, :) = (c(i-1, j, :) + c(i-1, j+1, :) + c(i+1, j, :) + c(i+1, j+1, :) + 6*c(i, j, :) + 6*c(i, j+1, :))/16;

new_c(2*i+1, 2*j, :) = (c(i, j-1, :) + c(i, j+1, :) + c(i+1, j-1, :) + c(i+1, j+1, :) + 6*c(i, j, :) + 6*c(i+1, j, :))/16;

new_c(2*i+1, 2*j+1, :) = (c(i, j, :) + c(i, j+1, :) + c(i+1, j, :) + c(i+1, j+1, :))/4;    

ref_coeff = new_c;
end

What am I doing wrong?

Answer 1

There is a good implementation of the Rueckert et al (1999) b-spline registration available on matlab central, written by Dirk-Jan Kroon.

I would advise you to make use of existing implementations like this since writing a classic registration algorithm from scratch means re-inventing the wheel.

Note that this b-spline registration implementation expects coarsely (e.g. rigid or affine) pre-registered data. It won't work if reference and target are e.g. completely differently scaled or strongly rotated etc. So pre-register first.