interpolating regular to irregular data locations using griddata

Hello geostats aficionados,

I'm using griddata to interpolate exhaustive, regular raster values to sparsely located sampling locations. My aim is the extract raster values at real world sampling locations. I have 78 sampling points, and the raster grid resolution is 1m. My out is a 78 x 1 array. Two concerns:

First MATLAB is throwing an error when viewing interpolation results because ZI is not a matrix. How can I view the interpolated points in the x,y,z coordinate space of the raster? And why aren't ZI values in a referenced matrix?

Second all extracted ZI values are NaNs - which are present only outside the extent of my (irregularly shaped) raster. Since XI,YI (target locations) and x,y (source locations) share the same projection and coords system, this is confusing.

Thanks!

% % EXTRACT SPATIAL INFO FROM REFERENCED RASTER [rast, rast_REFMAT, s1m_BBOX] = geotiffread('sedgwick_1m_utm10n_wgs84_updated.tif');

%%extract x,y vectors of raster nodes
X_size = rast_REFMAT(2,1); Y_size = rast_REFMAT(2,1);

%%create x,y meshgrid
[y, x] = meshgrid(s1m_BBOX(3):Y_size:s1m_BBOX(4)-Y_size, s1m_BBOX(1):X_size:s1m_BBOX(2)-X_size);
y = rot90(y); x = rot90(x); % rotate grids

x = double(x);
y = double(y);
z = double(rast);

%%LOAD CONTROL POINT (TARGET) COORDS
load xyvar
XI = xyvar(:,1);
YI = xyvar(:,2);

%%MESH TWO DATASETS USING GRIDDATA
ZI = griddata(x,y,z,XI,YI,'nearest');

figure
mesh(XI,YI,-ZI);

Error using mesh (line 76)
Z must be a matrix, not a scalar or vector.

Answer 1

[Walter Roberson]

mesh() cannot be used for scattered points. mesh() would not be able to derive the interconnections between the points in order to decide which points to connect to which.

If you do not wish to connect the locations, use scatter3() instead of mesh()

If you do wish to connect the locations, you could supply the known connection information and the data and create a patch(). If the connection information is not already known, you could create a Delaunay Triangulation and the use triplot() to plot it.

[Sam]

Thank you, Walter.

Scatter3() produces an empty plot with axes ranges that don't reflect the euclidean axes ranges of the source raster. Although the error is with mesh(), it seems the issue is with my griddata output? Or no...?

[Walter Roberson]

Did you try increasing the point-size for scatter? Did you assign the output of scatter3() to a handle and then get() on the handle to be sure that it has some actual data ? Check size(XI), size(YI), size(ZI)

The form of griddata() you used is defined to return a vector when XI and YI are a vector. mesh() cannot work on vectors.