MATLAB: Rough interpolation result using Griddata

griddata

I'm using meshgrid and griddata to interpolate scattered data in 2d. The input data seems to be very smooth, but i cannot retrieve such a smooth result with griddata. Does anyone has a suggestion how to to this properly? The figure depicting the results is attached.
Code:
clearvars; close all; clc;
load bulk.mat
X = bulk(:,1);
Y = bulk(:,2);
Z = bulk(:,3);
% now with meshgrid and griddata
xspace = linspace(0,5000,200);
yspace = linspace(0,2e6,160);
[xq,yq] = meshgrid(xspace,yspace);
zq = griddata(X,Y,Z,xq,yq);

Best Answer

It does not appear that the data you plot is the data you provide in the .mat file. If not, then how can we rationally help you?
The pictures that you show are of a function that is not apparently single valued, in the form z(x,y). The plots you show have a function that rolls over on top of itself. A single valued relation is crucial for griddata to work. (The same requirement would apply for my own gridfit. So don't try gridfit and then get upset when it fails to work well either.)
Note that the data you do provide appears to have a singularity at x=y=0. That will tend to make it difficult to produce a well behaved surface, due to the extremely high slopes.
Lets look at the data you have given us, plotted just as (x,y) pairs.
plot(x,y,'.')
MATLAB does not magically know how to connect these points. All you have given us is a list of those points as scattered pairs. So while your eye may be able to connect them into a surface, all that MATLAB can do is to form a delaunay triangulation in (x,y). It has no idea how the points are connected otherwise.
Here is a small portion of the triangulation that would have been formed:
These long, thin triangles are a problem in interpolation, as well as the choice of which triangles were formed. When you have a relationship with huge gradients, linear interpolation using long thin triangles is never going to work well. That, by the way, is what griddata uses here by default. Yes, it is true that the other methods of griddata will also work poorly on this data. Again, singularities produce problems. If singularities were easy to deal with, they would not have a special name for them, would they?
As well, the fact that x is scaled from 0 to 5e3 and y varies from 0 to 5e6 is also a problem, since it makes linear interpolation that is based on a delaunay triangulation difficult.
Worse, your data appears to be slightly concave along that top edge of the plot in (x,y). The result are some extremely long triangles along that edge. You can see them in the plot above. That will cause major problems.
One solution to the problem might be to rotate the entire data set in 3 dimensions, avoiding any issue of singularities, or single valued functions. That transformation might differ for each problem though, and you would need to understand how to do that.
A better solution would be to build the triangulation more intelligently. Then you could use tools to interpolate that functional form. But that would require you to learn to use such tools.
Related Question