CODEWITHSUNDEEP
matlabimage-processingbounding-box
user1384636
user1384636

Reputation: 501

extraction of minimal bounding box

I've got a volume (image) and a minimal bounding box which is wrapping a certain structure. I used the following algorithm to extract the minimal bounding box

http://uk.mathworks.com/matlabcentral/fileexchange/18264-minimal-bounding-box/content/minboundbox.m

So that I get a bounding box which looks like this

Minimal Bounding Box in the volume space

Now, I want to create a new volume from the MBB. In other words, I want to map every point within the MBB to a new box which is axis-parallel.

I can get the dimensions of the new box from the corner points

 dim = [0 0 0];
 x = cp(:,1);
 y = cp(:,2);
 z = cp(:,3);
 dim(3) = sqrt( (x(1)-x(2))^2 + (y(1)-y(2))^2 + (z(1)-z(2))^2 );
 dim(1) = sqrt( (x(1)-x(4))^2 + (y(1)-y(4))^2 + (z(1)-z(4))^2 );
 dim(2) = sqrt( (x(1)-x(5))^2 + (y(1)-y(5))^2 + (z(1)-z(5))^2 );

Now, I can apply the rotation matrix returned by the minboundbox algorithm, 

 A = zeros(4,4); A(1:3, 1:3) = R; A(4,4) = 1;
 tform = affine3d(A);
 N = numel(img);
 [X,Y,Z] = ind2sub(size(img), 1:N);
 V = img(1:N);
 [Xt, Yt, Zt] = transformPointsForward (tform,X,Y,Z);
 Xt = reshape(Xt, size(img));
 Yt = reshape(Yt, size(img));
 Zt = reshape(Zt, size(img));

And now I'm stuck. I need to:

any idea how can I do that?

UPDATE the main problem was not how to get Xt, Yt, Zt...

the problem, since Xt, Yt and Zt form an irregular grid of points, with certain values associated, how can I get a regular grid with the interpolated point back? I tried this but it went out of memory

Vq = griddata(Xt, Yt, Zt, double(V), 1:dim(1), 1:dim(2), 1:dim(3));

Upvotes: 1

Views: 369

Answers (1)

Rollen
Rollen

Reputation: 1210

It's easier to map the space of the MBB to the unit cube if you want to extract only the region of the MBB. Here's how I would do it (note this is mostly pseudocode, won't guarantee it works. I'll try it out on my matlab instance on my other computer soon...)

Translation = eye(4)
Translation(1:3, 4) = transpose(-cornerpoints(1,:))

BoxAxisX = transpose(cornerpoints(:,4) - cornerpoints(:,1))
BoxAxisX = BoxAxisX / norm(BoxAxisX)

BoxAxisY = transpose(cornerpoints(:,5) - cornerpoints(:,1))
BoxAxisY = BoxAxisY / norm(BoxAxisY)

BoxAxisZ = transpose(cornerpoints(:,2) - cornerpoints(:,1))
BoxAxisZ = BoxAxisZ / norm(BoxAxisZ)

Rotation = eye(4)
Rotation(1:3,1:3) = [ BoxAxisX, BoxAxisY, BoxAxisZ ]

Scale = diag( [ (1/dim(1)) (1/dim(2)) (1/dim(3)) 1 ] )

A = affine3d( Scale * Rotation * Translation )

[Xt, Yt, Zt] = transformPointsForward (A,X,Y,Z);

% Now all your points that are within the box should satisfy unit constraints
GoodIndices = Xt <= 1 && Xt >= 0 && Yt <= 1 && Yt >= 0 && Zt <= 1 && Zt >= 0
GoodX = Xt(GoodIndices)
...

The general idea is to turn the box back into a unit cube with one point on the origin. Here are the steps:

  • Choose a base point; this will become your origin.
  • Find vectors from the base point that form the bounds of the box; this is your current coordinate frame.
  • Construct the rotation matrix from those vectors; This is the matrix going from the standard basis vector (X,Y,Z) to your box rotation using the base point you've chosen. Invert this matrix.
  • Translate the box back such that the base point is now at the origin
  • Now use that inverted rotation matrix and apply that to your box.
  • Note your box is now axis aligned at the origin. Now all you need to do is scale your box appropriately.

The reason why I don't use the rotation matrix they give is that I don't really know which axis it will go towards and what point it rotates around (rotating around the origin! what!). So... It is just easier to construct the rotation matrix using the corner points you know allowing you to choose the base point.

As for the Out of Memory Exceptions...there are lot of posts on how to reduce memory usage. That image size is huge though..(~240MB)

Upvotes: 1

Related Questions