How to find the difference in elevation of a point above a mesh?

Hi,

I have 2 datasets of x,y,z (Lidar and Total Station) and what am trying to do is produce a trimesh of the Total Station data and find the difference of elevation for each Lidar point above the surface.

Since the 2 data sets do not overlap in x or y, I was hoping this could be a simple way to compare elevation differences.

BUT...the output I get from the dalaunay gives me the vertices and I don't know how to subtract the lidar xyz from the trimesh.

Any suggestions would be greatly appreciated.

My code so far:

% generate random xyz
ts = rand(100,3);

% offset lidar data in x,y,z
lidar = [ts(:,1)*1.1,ts(:,2)*1.1,ts(:,3)*2];

% compute delaunay
dt = delaunay(ts(:,1:2));

% compute convexhull to look only at 'in' points 
k = convhull(ts);
in = inpolygon(lidar(:,1),lidar(:,2),ts(k,1),ts(k,2));

% plot mesh
trimesh(dt,ts(:,1),ts(:,2),ts(:,3))
hold on

% plot lidar
plot3(lidar(in,1),lidar(in,2),lidar(in,3),'k*');

view(90,90)
xlabel('x')
ylabel('y')
zlabel('z')

Answer 1

[Norris]

So i figured it out!

I ended up using the interp2 function to interpolate the x,y on the surface.

Here is the code for anyone else to use.

% generate random xyz
ts = rand(100,3);
x_ts = ts(:,1);
y_ts = ts(:,2);
z_ts = ts(:,3);

% offset lidar data in x,y,z
lidar = [ts(:,1)*1.1,ts(:,2)*1.1,ts(:,3)*2];
x_lidar = lidar(:,1);
y_lidar = lidar(:,2);
z_lidar = lidar(:,3);

% compute delaunay
dt = delaunay(ts(:,1:2));

% compute convexhull to look only at 'in' points 
k = convhull(ts);
in = inpolygon(lidar(:,1),lidar(:,2),ts(k,1),ts(k,2));

% linearly space vectors for meshgrid
cellsize = 100;
xlin = linspace(min(x_ts),max(x_ts),cellsize);
ylin = linspace(min(y_ts),max(y_ts),cellsize);

% make meshgrid
[x_gr, y_gr] = meshgrid(xlin, ylin);

% interpolation
f = TriScatteredInterp(x_ts,y_ts,z_ts,'linear');
z_gr = f(x_gr,y_gr); % get z values for grid

% interpolate to get z on mesh
z_surf = interp2(x_gr,y_gr,z_gr,x_lidar,y_lidar);

figure
% plot mesh
mesh(x_gr, y_gr, z_gr)
hold on
plot3(x_lidar,y_lidar,z_surf,'k*') % z on mesh
plot3(x_lidar,y_lidar,z_lidar,'o') % orig z

view(90,90)
xlabel('x')
ylabel('y')
zlabel('z')