MATLAB: Is ScatteredInterpolant (linear) returning a value far outside of input data

buginterpolationMATLABscatteredinterpolant

So I'm using ScatteredInterpolant to interpolate the values of a measurement on the surface of a helicopter blade. I'm showing it graphically in the figure. The input data is shown with the faint lines at a 45 degree angle. The output is the coloring of the surface. Right in the center is a bright spot. The highest input data value is ~119. The output value at that location is about 2512.

I've attached the data which generated that figure. The function used to generate the Interpolation is:

Interpolant = scatteredInterpolant(inX, inY, inZ, inV,'linear','nearest');
outV = Interpolant(outX, outY, outZ);

For reference, the bad output data point is the 2713 entry in outV. How can I possibly get an output value that is more than an order of magnitude bigger than any input value, when using linear interpolation? Is this a bug?

A quick way to visualize the error:

figure;
plot3(inX, inY, inV, 'o');
hold on
plot3(outX, outY, outV, 'o');

Best Answer

You are making a major error here in what you are trying to do.

plot3(inX,inY,inZ,'bo')
box on
grid on

Now, rotated one way, we see this surface.

Now, rotated slightly differently, we see this:

So you have essentially a nice surface in 2-dimensions. But there is NO information about what happens when you move off that surface. Now, add in the points you are trying to "interpolate".

hold on
plot3(outX,outY,outZ,'r+')

Yep. You are trying to predict values that lie HUGELY away from the 2-manifold of data that you have. Massive extrapolation here.

So, regardless of what scatteredInterpolant did, your results will be complete, meaningless garbage. You don't have three independent variables, comprising a true point cloud. You have a 2-manifold of data.

So what will scatteredInterpolant try to do here? Linear interpolation for any point will involve a tessellation of that domain. It will not be some sort of triangulation of the surface, but a set of tetrahedra spanning, dissecting that domain. At each of those points, you have a value of V that you wish to interpolate.

So scatteredInterpolant uses a tessellation of that domain. Thinking that you have a real point cloud.

tess = delaunayn([inX,inY,inZ]);
whos tess
Name          Size              Bytes  Class     Attributes
tess      43112x4             1379584  double

So there are 43112 tetrahedra. Its kind of hard to visualize, since there are so many edges. But the volume looks like this:

Most of those tetrahedra are HUGE things. Spanning across that void where no data exists. scatteredInterpolant has no clue what you did. At any poiint in that domain, it finds the tretrahedron that contains your point. Then it does LINEAR interpolation across that tetrahedron. But remember that most of the tetrahedra are huge things, connecting one "wing" of that surface with the other. (No pun intended with use of the word wing here, on what is essentially airfoil data.)

ScatteredInterpolant just does what it is told, having no idea that when you try to interpolate some point in that volume, it is creating meaningless gibberish as a result. I'm sorry, but you simply cannot use scatteredInterpolant to produce a meaningful result from this data, as you are trying to do.

Do I really need to figure out why you got something mildly strange? Garbage in, garbage out. I guess I could go a bit more deeply into what you have, if you really want it. But it would be a pure exercise in trying to debug code that I don't have, because as I recall, scatteredInterpolant is supplied to us in compiled form.

Related Solutions

MATLAB: Scattered interpolation unexpected peaks in extrapolation results

You might appropriately rescale your x/y data before interpolation

A=load('PV_curve.txt');
x=A(:,1);
y=A(:,2);
z=A(:,3);
xscale = 1/10; % <= experiment with it
yscale = 1;
xn = x*(xscale/max(x(:)));
yn = y*(yscale/max(y(:)));
xi = linspace(0,xscale);
yi = linspace(0,yscale);
[X,Y] = meshgrid(xi,yi);
I = scatteredInterpolant(xn,yn,z,'linear','linear');
Z = I(X,Y);
xi = xi*max(x(:))/xscale;
yi = yi*max(y(:))/yscale;
close all
figure
surf(xi,yi,Z);
hold on
plot3(x,y,z,'.r');

MATLAB: How do you fit a curve through 3D points, using splines, constrained within a volume.

Ugh, enough to say about this that it seems worth an answer instead of just a comment on Sean's answer.

It sounds like you have a facet of a tetrahedron, and a set of points that lie inside that facet, and are in the plane of the facet. Now, you want to find a spline that passes through the points, and is ensured to lie entirely inside the facet, as well as staying entirely in-plane.

The obvious answer is you can project the points into a 2-dimensional subspace, one that is defined by the facet. Then one could do any fits needed in the planar subspace, and the problem is completely a 2-d problem. At the end, you can reverse the projection, returning the interpolated points back into 3-d.

In fact, I think this will not be necessary, since a spline interpolant will result in a linear combination of the points. If the points are all in-plane, then any linear combination of those points will also be in-plane too. (I'll think about this claim, but I'm pretty confident it will hold.)

The interpolation itself can be done using tools like pchip, my own interparc, etc. A problem is how to enforce that the points lie inside a polygonal domain as defined by that facet. The suggestion of vert2con will fail here however. Why?

vert2con will produce a set of linear constraints that could in theory be applied to the spline. But the spline is a nonlinear function, even though it is linear in the spline parameters. Those linear inequalities cannot be applied to the spline, as a FUNCTION. You could sample the spline at a gazillion (the technical term for a LOT) points, then apply those inequalities to each point. Essentially this will be a NECESSARY constraint on the spline, that it lies inside the volume, but it will not be a sufficient condition. This is the same as how monotonicity constraints can be applied to a spline. If you constraint the spline to be monotone at a finite set of points, then you will have a necessary, but not sufficient condition that the spline be monotonic as a function. So there could be a small region (or regions) where the spline actually fails to satisfy the constraints.

And of course, you need to be careful. For example, suppose the points happen to lie at the vertices of the triangular facet? For example, suppose you choose the points as rows of xy below:

xy = [1 0;1 1;0 1];

They define a triangular domain. Can you now find a spline that lies entirely inside that triangular region, AND also interpolate those points? Of course not, if the spline interpolant will be at least a C1 function. No matter how you define a smooth interpolant through those points, if it is differentiable at the vertices of the triangle, then it must go outside of the triangle at that second vertex.

xyt = [0 0;0 1;1 0];
xy = rand(100,2);
ind = sum(xy,2)<1;
xy = xy(ind,:);
plot(xyt(:,1),xyt(:,2),'ob-')
hold on
plot(xy(:,1),xy(:,2),'g+',xyint(:,1),xyint(:,2),'r:')

As you can see, the interpolated curve falls entirely inside the triangular domain, mainly because I used a pchip interpolant for the purpose. That helps, but it won't ensure the curve always does. I'm quite confident that with only a tiny amount of effort, I can come up with an example where that fails. For example...

xy = [0.1 0.1;0.3 0.69;0.69 0.3;0.2 0.2];
xyt = [0 0;0 1;1 0];
xyint = interparc(1000,xy(:,1),xy(:,2),'pchip');
plot(xyt(:,1),xyt(:,2),'ob-')
hold on
plot(xy(:,1),xy(:,2),'g+',xyint(:,1),xyint(:,2),'r:')

I could easily enough show you how to do the interpolation using pchip directly, without need to use interparc here. Regardless, the point is that ensuring the curve falls entirely inside a polygonal domain is difficult, UNLESS your spline is a piecewise linear one. Then it is entirely trivial. But few people ever say spline when they mean a piecewise linear (connect-the-dots) function.

xyint = interparc(1000,xy(:,1),xy(:,2),'linear');
plot(xyt(:,1),xyt(:,2),'ob-')
hold on
plot(xy(:,1),xy(:,2),'g+',xyint(:,1),xyint(:,2),'r:')

As you can see, a piecewise linear interpolant will never have a problem, but since it is just connect-the-dots, I'm not sure you would be happy with that result. Any smooth interpolant will potentially have problems with some set of points.

Again, the fact that I've done my examples in 2-d is not relevant, since one can always use projective means to turn a 3-d problem that lies entirely in a plane into a 2-d one. And given my claim that sine a spline interpolant will yield a linear combination of the points, then you don't really need to do the projection at all. Any linear combination of a set of points that lie in a plane will result in a point that also lies in the same plane.

If I've misunderstood the question, please explain and I'll see what I can do.

Best Answer

Related Solutions

MATLAB: Scattered interpolation unexpected peaks in extrapolation results

MATLAB: How do you fit a curve through 3D points, using splines, constrained within a volume.

Related Question