Interpolation across only 1 dimension for a multi-dimensional array

I have an 3D-array V that is 2x1000x5 and contains the values of two functions (dimension 1) each of the same two inputs (dimensions 2 and 3). For resource reasons, I cannot increase the the number of points in each input grid before the calculation of V, but I would still like to see estimates of the function along the intermediate points, hence the need for interpolation. My goal is to interpolate accross the first input (dimension 2), while leaving the other two dimensions unchanged.

The code I have is the following:

n_a=1000; %original size of dimension 2
int_n_a=n_a*10; %interpolated size of dimension 2

a_grid = (linspace(0,1,n_a).^3)'
int_a_grid = (linspace(0,1,int_n_a).^3)'

n_z=5; %size of dimension 3
z_grid = [0.53;0.71;0.96;1.30;1.76]

int_V=interpn(1:2,a_grid,z_grid,V,1:2,int_a_grid,z_grid)

However, this does not seem to do what I intended. The code does not fail, but the values within the object int_V are clearly not the result of the intented interpolation. What is the correct specification for interpolating a multi-dimensional array along only 1 dimension?

Answer 1

[Walter Roberson]

Use interp1(). The first input must be a vector; the second input can be an array of any size; permute() it so that you are interpolating over the first dimension.

[Anthony Swaminathan]

This works perfectly, but I run into a follow-up issue of the array indices no longer being positive integers, presumably because there are new indices in between the original integer indices due to interpolation. Is there a way for me to re-index to return to having only positive integers?

[Walter Roberson]

You wouldn't normally need the original indices afterwards. You would normally create a query vector of points (typically equi-distant, not always) and that is typically all you would need.

Every once in a while, it is useful to calculate something along the lines of

RI = interp1(MonoTonicValues, 1:length(MonoTonicValues), QueryVector)

where 1:length is effectively indices, and what you get out of the interp1 is "relative position" of where the query vector values would be in the list, expressed as a fraction -- so for example 3.81 indicating 0.81 of the way between 3rd and 4th point. Those relative values can then be used to calculate based on other related arrays, at reduced cost,

flRI = floor(RI);
off1 = RI - flRI;
off0 = RI - off1;
Interpolated = off0 * SomeVector(flRI) + off1 * SomeVector(flRI+1);

The need to do this is reduced as of R2023b, in the case that scatteredInterpolant is being used:

R2023b release notes:

The scatteredInterpolant object can now interpolate multiple data sets at the same query points. Specify the Values property as a matrix, where the number of rows is the same as the number is sample points and each column in Values represents the values of a different function at the sample points. For example, if the sample points are column vectors with 10 elements, you can specify Values as a 10-by-4 matrix to interpolate using four different sets of values.