how can I create a 4D matrix and use it for interpolation

Hello everybody,

I have 3 vectors conditions:

first vector: Wpig=0.5:0.5:7;

second vector: Chla=0.5:0.01:0.6;

third vector: Starch=0:0.04:0.24;

Wavelength:linspace(300,700,31)

for each possible condition I calculate Ea using a function I named "proprad(Wpig,Chla,Starch)", so for exemple for proprad(0.5,0.5,0)=Ea and Ea is a vector composed by 31 values (it's for this I created a vector named Wavelength having the same lenght that Ea vector result). Summarizing, usign my function, when you put 3 scalars: Wpig,Chla,Starch, you get a vector result (Ea).

So, I want to create a 4D matrix in wich the 3 firsts dimensions will correspond to my 3 conditions (Wpig,Chla,Starch) and the 4th dimension will be the vectors estimated by my function.

My first idea is to create a zeros 4D matrix like this:

Mat_Ea=zeros(length(Wpig),length(Chla),length(Starch),length(Wavelength));

and then remplace each dimension by the value, giving something like this:

1 dimension : Wpig

2 dimension: Chla

3 dimension: Starch

4 dimension: Ea

finally, my goal is to use the 4D matrix, for interpolation. So for exemple, if I want to know what's the value of Ea for conditions like Wpig=0.93,Chla=0.57 and Starch=0.12, I cand interpolate in the 4D matrix to estimate it.

Could someone help me to solve this coding problem?

Thank you very much for your time and your help!

Answer 1

[Fabio Freschi]

% your vectors
Wpig=0.5:0.5:7;
Chla=0.5:0.01:0.6;
Starch=0:0.04:0.24;

Now you need to sample uniformly the space using your vectors. Use may use ndgrid for this

[X,Y,Z] = ndgrid(Wpig,Chla,Starch);

X, Y, and Z are now 3d matrices, where

• X: Wpig
• Y: Chla
• Z: Starch

Note that you have numel(Wpig)*numel(Chla)*numel(Starch) = 1078 elements in these matrices.

You must now evaluate your function. If it admits matrices as inputs, you can supply directly the 3d matrices. For example

Ea = cos(X)+sin(Y).*tan(Z);

If it admits vectors as inputs or you evaluate each value of Ea providing only scalar values, e.g. in a for loop, you can create Ea as a vector, then reshape it. For example

% this three lines are not necessary, but makes things more clear
X1 = X(:);
Y1 = Y(:);
Z1 = Z(:);
% preallocation
Ea = zeros(numel(X),1);
% evaluation
for i = numel(X1)
    Ea = cos(X1(i))+sin(Y1(i)).*tan(Z1(i));
end
% reshape the result vector to fit X,Y,Z dimensions
Ea = reshape(Ea,size(X));

Now you can use intepolation. interpn can interpolate n-d gridded data in WpigQ, ChlaQ, StarchQ query points

% query points
WpigQ = [5.5; 6.3];
ChlaQ = [0.53; 0.58];
StarchQ = [0.111; 0.134];

% interp
EaInterp = interpn(X,Y,Z,Ea,WpigQ,ChlaQ,StarchQ);

[Fernando Robert Ferrel Ballestas]

Hi Fabio!

Thanks in advance for answering my question. The grid you propose me to do it's a good idea to organize the conditions I want to evaluate in my function "proprad(Wpig, Chla, Starch)". However, as you said, I will have 1078 elements in this matrix and for each combination I will obtain a vector of 31 elements, so evaluating the 1078 elements in my function, I will have 1078 vectors "Ea" with length 31.

So, I don't know if knowing this information you suggest me to do anyway this?:

% preallocation

Ea = zeros(numel(X),1);

% evaluation

for i = numel(X1)

Ea = Ea = cos(X1(i))+sin(Y1(i)).*tan(Z1(i));

end

% reshape the result vector to fit X,Y,Z dimensions

Ea = reshape(Ea,size(X));

Moreover, X doesn't have the same dimension of Y , Z and Ea ...

Thank you a lot one more time!

[Fabio Freschi]

I think I have missed something.

1) in my code X,Y,Z have the same dimensions

2) if you start the code with

% your vectors
Wpig = 0.5:0.5:7;
Chla = 0.5:0.01:0.6;
Starch = 0:0.04:0.24;
Wavelength = linspace(300,700,31);

% X,Y,Z,W are 4d matrices (14*11*7*31)
[X,Y,Z,W] = ndgrid(Wpig,Chla,Starch,Wavelength);

Could you fill Ea so that it hase the same dimesnsions as X,Y,Z,W?

If so, the code can be easily modified for your case

[Fernando Robert Ferrel Ballestas]

No, i can't. because for Wpig, Chla and Starch scalar combinaison like for exemple:

proprad0.5,0.5,0) I obtain Ea vector of 31 elements. (31x1) .

I put the vector Wavelength because it has the same dimension of Ea and it could be used for the 4D-matrix construction but I'm not sure.

What do you recommend me to do? Thanks!