Interpolation of a 2D image series at specific slice locations

I have imported a series of 2D DICOM files containing 14 MR images. The resolution of the image data matrix is 512x512x14. These MR slices are located in the axial plane at real world slice locations -40.063mm to 11.937mm in 4mm intervals. I want to reslice this series along the same axis to have slice locations at xx.466mm in 1mm intervals such that I end up with a new data matrix with resolution 512x512x52 (slice locations from -39.466mm to 11.534mm). Should trilinear interpolation be applied in this case? Could someone help me with this interpolation problem? Also, the voxel resolution is (0.3516 x 0.3516 x 4)mm/pixel.

Many thanks in advance!

Answer 1

[Matt J]

You should use interp1, and a little bit of reshaping:

 zold=linspace(-40.063,  +11.937, 14).';
 znew=linspace(-39.466 ,  11.534, 52).';

array=reshape(MRvolume,[],14).';

new=interp1(zold,array,znew).';

result=reshape(new, [512,512,52] );

[Bualookkaew Sakchatchawan]

May I ask for my question? If I have a series of 2D axial CT with data matrix as 512x512x78, but I want to show the real 3D volume. How can I perform volume with correct the pixel spacing (0.7813) and slice thickness at 2 mm?

Sine

[Matt J]

@Bualookkaew Sakchatchawan

What about the solution already proposed? Why not just use that?

[Sebastian Joseph]

What if the case of paired images and intetpolation..??

[Matt J]

@Sebastian Joseph, I'm not sure what you mean by a "paired image". If you simply mena you have two images, you can either apply the interpolation steps above to each image individually.

[Sebastian Joseph]

Paired image I meant source images and if I want to extract middle image (between left and right) from interpolation hw I should proceed..if I'm applying to each image individually then not getting the middle one. On this I'm getting result either similar to left image or right image..

So confused why it happens..?

[Sebastian Joseph]

%% clc; clear all; close all;

% Initialization of Input Images

A1 = im2double(rgb2gray(imread('cuphex 0degree.png'))); A5 = im2double(rgb2gray(imread('cuphex 50degree.png')));

% Fourier Transorm of Input Image

A1 = fftshift(fft2(ifftshift(1.*A1))); A5 = fftshift(fft2(ifftshift(1.*A5)));

[M,N] = size(A1); Nx = N; Mx = 2*N;

height = size(A1,1); width = size(A1,2);

%Taking Vectors from input images

for ii = 1:height

    input_fftA1 = A1(ii,:,:);
    input_fftA1(:,1:2:2*N) = A1(ii,:,:);

    [x y] = size(input_fftA1);
    [x0 y0] = meshgrid(1:y, 1:x);

    input_fftA5 = A5(ii,:,:);
    input_fftA5(:,2:2:2*N) = A5(ii,:,:);

    [x1 y1] = size(input_fftA5);
    [x01 y01] = meshgrid(1:y1, 1:x1);

% figure(1), imagesc(abs(input_fftA1)); % figure(3), imagesc(abs(input_fftA5));

end

% Interpolation Plane mj1 = 1.5:2:2*Mx+1; mj2 = 1.5:2:2*M+1; mj3 = mj2'+mj1;

% figure(), plot(mj1); % figure(), imagesc(mj1); % hold all % figure(), plot(mj2); % figure(), imagesc(mj2); % figure(), plot(mj3); % figure(), imagesc(mj3);

% Interpolation

 Out_img = interp2(x0,y0,mj3,x01,y01, 'linear'); %spline , 'extrap'

 imagesc(Out_img);

 out = Out_img;
 out(isnan(out))=0;%replace Nan
 plot(out);
out1=reshape(out,length(N),[]);

plot(out1); imagesc(out1);

% IFFT of the result % Fourier Transorm of Input Image

A3 = fftshift(ifft2(ifftshift(out1))); imagesc(A3);

[Sebastian Joseph]

Something wrong in the paste command or I did wrong..the link is

[Sebastian Joseph]

https://www.mathworks.com/matlabcentral/answers/803986-how-to-do-interpolation-for-paired-images#comment_1475700