I have one 3D dataset and one 2D dataset which is a slice through the first volume. They are at different scales, resolutions and in a different coordinate system, but of both I know the affine transformation to world coordinates. I then think I know how to apply those, but how do I get an image back from those transformed coordinates again using sinc interpolation? I would like to learn how to perform this/how this works. The first comments below already pointed me at existing functions inside matlab that perform linear interpolation, but I would also like to know how to do this myself so I can use sinc interpolation (and others).
I can round the coordinates and get values for those, which would be nearest neighbor interpolation. I would like to lose as little information as possible irrespective of computation time, I think I should use sinc interpolation then. When I have the transformed coordinates, how do I make (for example sinc) an interpolation algorithm?
for example:
%%get data
A = rand(200,250,250); % I want to get the slice from this that corresponds to B
B = rand(1200,1200); % I want to get the data from A that corresponds to the same pixels
%%get coordinates for A
siza = size(A); clear xx;clear yy;clear zz;
[xx,yy,zz] = meshgrid(1:siza(1), 1:siza(2), 1:siza(3));
coora = [xx(:)';yy(:)';zz(:)'; ones(size(zz(:)))'];
%%get coordinates for B
sizb = size(B); clear xx;clear yy;clear zz;
[xx,yy] = meshgrid(1:sizb(1),1:sizb(2)); zz = zeros(size(xx));
coorb = [xx(:)';yy(:)';zz(:)'; ones(size(zz(:)))'];
%%define affine transformation matrices
T3d = [-0.02 0.02 1 -88 ;
-0.98 0 -0.02 130;
0 0.98 -0.02 -110;
0 0 0 1 ];
T2d = [-0.2 0 0 126;
0 0.2 -0.2 -131;
0 0 2 43 ;
0 0 0 1 ];
%%transform A coordinates to world coordinates and world coordinates to B coordinates
cooraInBref = T3d*inv(T2d)*coora;
aslice = zeros(size(B));
%% then nearest neighbor would go something like this (dont exactly know how to do this either):
cooraInBround = round(cooraInBref);
for idx = 1:length(coorb);
if cooraInBround(3,idx) == 0
aslice(cooraInBround(1,idx),cooraInBround(2,idx)) = ...;% dont know how to do this
end
end
%% how would I implement sinc interpolation instead of rounding the transformed coordinates
Related questions that dont quite help me further:
Nearest-neighbor interpolation algorithm in MATLAB
"Nearest neighbor"-like interpolation in MATLAB
How to apply an affine transformation (4x4 matrix) to ndgrid/meshgrid results?
Python/PIL affine transformation
How do I rotate a 3D matrix by 90 degrees counterclockwise?
Resizing a 3D image (and resampling)
Ready to use functions from matlab as pointed out in comments by chappjc anonsubmitter85 and cape code:
http://mathworks.com/help/matlab/math/interpolating-scattered-data.html
http://mathworks.com/help/matlab/ref/griddata.html?refresh=true
My other SE question I use to get the affine matrix
Using scatteredInterpolant as suggested in the comments and the linked question took almost ten minutes and resulted in a slice with all NaN`s, am now trying other methods.
I now start the other way around: get the 2d points to sample and then transform those to the 3d volume and then use interp3 to calculate the values. Its pretty quick and works well. This code only works for getting a single slice, but i think you can easily adapt it to get a whole transformed volume. I still dont know how to do sinc interpolation though.
Still not sure if I should have used zeros or ones for Z.