Interpolation/Align two nonuniform grid in same time and space and convert it to a uniform Grid where both values are known at the same point.

Hello everyone,

I am not the professinal with Matlab and I have a problem to solve. I got a commercial program, which gives me the values of E-field and W/kg for a phantom.

When I export the values so I can use it in Matlab, the two coordinate systems doesn't fit to each other like one to one. I have made an illustration in 2D but it is a 3D problem tho.

It is two non-uniform grids, but it is in the same space

In my case I got two coordinate system, which you can see down below: The distance from each pixel/voxel, just like the E-field grid is between 1.9mm-0.9mm. I need to make the grid uniform and change the distance to 2mm for each pixel/voxel. And I have to know the E and Kg/W at the midpoint because I need them for further calculations. The E-field value is a complex number, which include real and imaginary numbers but we can forget about that for now.

So when I got my data, I got the x,y,z coordinates which are vectors and an E-value and a S-value, which fits with the coordinates as a n*m*l (matrix) since we need 3 positions to place the point.

if true
  % code
load('E1.mat')
Ex=x_axis*1000; % 195x1
Ey=y_axis*1000; % 190x1
Ez=z_axis*1000; % 181x1
E %195x190x191
load('SAR1.mat')
KgWx=x_axis*1000; % 194x1
KgWy=y_axis*1000; % 189x1
KgWz=z_axis*1000; % 180x1
KgW %194x189x180
% Creating a mesh
[Exx, Eyy, Ezz] = meshgrid(Ey,Ex,Ez);
[KgWxx, KgWyy, KgWzz] = meshgrid(KgWy,KgWx,KgWz);
% Making a meshgrid with a spacing of 2 making more points
[Exx0, Eyy0, Ezz0] = meshgrid(min(Ex):2:max(Ex), min(Ey):2:max(Ey), min(Ez):2:max(Ez));
[KgWxx0, KgWyy0, KgWzz0] = meshgrid(min(KgWx):2:max(KgWx), min(KgWy):2:max(KgWy), min(KgWz):2:max(KgWz));
% As I said I also got my E-value and kgW-value corresponding to the x,y & z, so I though I could regard E & kgW as *v* for griddata(x,y,z, *v*,xq,yq,zq);
% The griddata function doesn't work somehow...
%vq = griddata(Exx,Eyy,Ezz,E,Exx0,Eyy0,Ezz0);
end

But I couldn't get further than that.... and I don't know to to combine the grids so they overlap each other.

It will be a great help if someone got an idea how to solve it, and thank you for your time! :)

Answer 1

[Matt J]

Since the original samples of both E and KgW are gridded, it sounds like you could just use do gridded interpolation twice. The only thing you need to do is specify different initial grid vectors {ne,me} for E and {nkgw,mgw} for KgW:

    Fe=griddedInterpolant({ne,me},E)
    Fkgw=griddedInterpolant({nkgw,mgw},KgW);

    E_resamp=Fe({nCommon,mCommon},E);
    KgW_resamp= Fkgw({nCommon,mCommon},E);

and a common set of grid vectors {nCommon,mCommon} for the target sampling.

[Nick Yin Larsen]

Hello Matt and thank you for your fast reply :)

I will try to do as you sugested tmr, since I was dumb and forgot to export the KgW-values to my laptop...

So the dimensions for my E-field is 195x190x181 and since the E-field include a real and imaginary value for each axis I have to code something like this when I use :F = griddedInterpolant(X,Y,Z,V);

EX= eRealX+1i*eImgX; EY= eRealY+1i*eImgY; EZ= eRealZ+1i*eImgZ;

E = EX+EY+EZ; F = griddedInterpolant(195,190,191,E);

This part: E_resamp=Fe({nCommon,mCommon},E); KgW_resamp= Fkgw({nCommon,mCommon},E);

What should I place inside the nCommon, and do you have a suggestion for the equal space of 2mm for each pixel/voxel? :)

[Matt J]

I have to code something like this when I use

It's hard to read. You need to mark-up your code using the "{} Code toolbar" button.

What should I place inside the nCommon, and do you have a suggestion for the equal space of 2mm for each pixel/voxel?

The Grid Vectors nCommon, mCommon, etc... are the coordinates of the locations that you want to interpolate at. If 2 mm is the sample spacing that you want, that would be the place to put it.

[Nick Yin Larsen]

Okay, I didn't see that function for the code. Here is my code after what I understand: (I replace Wkg with S and SAR)

clc; clear all;
load('E1.mat')
Ex=x_axis*1000;
Ey=y_axis*1000;
Ez=z_axis*1000;
load('SAR1.mat')
Sx=x_axis*1000;
Sy=y_axis*1000;
Sz=z_axis*1000;
eRealX = real(Real_Part_of_X_Comp_0s);
eRealY = real(Real_Part_of_Y_Comp_0s);
eRealZ = real(Real_Part_of_Z_Comp_0s);
eImgX = imag(Imaginary_Part_of_X_Comp_0s);
eImgY = imag(Imaginary_Part_of_Y_Comp_0s);
eImgZ = imag(Imaginary_Part_of_Z_Comp_0s);
EX= eRealX+1i*eImgX;
EY= eRealY+1i*eImgY;
EZ= eRealZ+1i*eImgZ;
E = EX+EY+EZ;
SAR0 = Value_of_SAR_0s;
Fe=griddedInterpolant({Ex,Ey,Ez},E);
FSAR=griddedInterpolant({Sx,Sy,Sz},SAR0);
% Last part doesn't work tho :(
% E-field part
d1 = min(Ex):2:max(Ex);
d2 = min(Ey):2:max(Ey);
[Xq,Yq] = ndgrid(d0,d1);
Vq = F(Xq,Yq);
mesh(Xq,Yq,Vq);
% SAR part
F = griddedInterpolant(X1,Y1,SAR0(:,:,100));
dim1 = min(Sx):2:max(Sx);
dim2 = min(Sy):2:max(Sy);
[Xq,Yq] = ndgrid(dim1,dim2);
Vq = F(Xq,Yq);
mesh(Xq,Yq,Vq);

Well it kind of works but it is only for 2D, and I don't think there is any interaction between the two results?

Is it possible if I can get your email so I can share a folder with you on Dropbox?

My email is : [e-mail]

Ones again Thank you for your time! :)

[Matt J]

griddedInterpolant works fine in 3D, or in any dimension, as the following simple example shows. I think the solution is for you to work with such examples until you have a full command of how griddedInterpolant works.

    >> A=reshape(1:8,2,2,2)

    A(:,:,1) =

         1     3
         2     4

    A(:,:,2) =

         5     7
         6     8

    >> ax=[-1000,1000]; F=griddedInterpolant({ax,ax,ax}, A);

    >> samps=[-1000, 0, 1000]; B=F({samps,samps,samps})

    B(:,:,1) =

        1.0000    2.0000    3.0000
        1.5000    2.5000    3.5000
        2.0000    3.0000    4.0000

    B(:,:,2) =

        3.0000    4.0000    5.0000
        3.5000    4.5000    5.5000
        4.0000    5.0000    6.0000

    B(:,:,3) =

        5.0000    6.0000    7.0000
        5.5000    6.5000    7.5000
        6.0000    7.0000    8.0000

[Nick Yin Larsen]

Okay,

So I should code something like this:

KgWx=x_axis; % 194x1
KgWy=y_axis; % 189x1
KgWz=z_axis; % 180x1
KgW % 194x189x180

F=griddedInterpolant({KgWx,KgWy,KgWz}, KgW);

Ex=x_axis*1000; % 195x1
Ey=y_axis*1000; % 190x1
Ez=z_axis*1000; % 181x1
% Create a new grid
KgW_new=F({Ex,Ey,Ez}) % 195x190x181

So What I have done is creating a new grid for KgW, so it fits with the E-grid if I am correct right? So now actually the grid's value is at the same point? :) Thank you very much Matt if my assuption is correct!

So right now it is interpolated but it is still non-uniform. So I can just make it uniform after this right?

[Matt J]

So right now it is interpolated but it is still non-uniform. So I can just make it uniform after this right?

Yes, but I wouldn't do it by interpolating twice. In other words, there is no need to first resample KgW onto a common grid with E and then resample again to make the grids uniform. That will compound interpolation/resampling error. Just resample each array separately onto the final uniform grid that you want, like I suggested in my original answer.

[Nick Yin Larsen]

Great! - everything is working now! Thank you for your time Matt!

And with the help from you it turned out to be:

Maybe it is hard to see here but now it is a uniform grid :D