Matlab Delaunay Triangulation of Point Cloud - Color Matrix

space-dementia picture space-dementia · Sep 13, 2012 · Viewed 12k times · Source

I would like to create a plot of the 3D surface that spans over all points from an [X,Y,Z] point cloud. For example this is a scatter plot of my point cloud:

scatter3(X,Y,Z,5,C)

Scatter plot

As you can see each data point has an intensity value C.

I now carry out the triangulation

dt      = DelaunayTri(X,Y,Z); 
[tri Xb]= freeBoundary(dt); 

And I get the triangulated surface

figure 
trisurf(tri,Xb(:,1),Xb(:,2),Xb(:,3), 'FaceColor', 'cyan', 'faceAlpha', 0.8);

Surface

However, when I try to set the colour of the surface using

trisurf(tri,Xb(:,1),Xb(:,2),Xb(:,3),C,'EdgeAlpha',0,'FaceColor','interp')

I get the error message: "Warning: Color Data is not set for Interpolated shading", which is due to the fact that the size of C does not match Xb or tri.

How can I make sure I get the correct interpolated surface colour?

Answer

angainor picture angainor · Sep 13, 2012

You have changed the number of points in the plotted triangulation by calling freeBoundary: only the surface points are left, the inner points do not belong to the surface. Therefore, you have to extract C values that correspond to those points. You can use 'intersect(..., 'rows')' to map the surface points Xb onto the original point set XYZ. Based on this map you extract the needed values from C. The code below does this.

clear all;

XYZ = rand(100,3);
X=XYZ(:,1);
Y=XYZ(:,2);
Z=XYZ(:,3);
C=rand(size(X));

scatter3(X, Y, Z, 5,C);

dt = DelaunayTri(X, Y, Z);
[tri Xb]=freeBoundary(dt);

% map Xb onto XYZ
[~,IA,IB]=intersect(XYZ, Xb, 'rows');

% extract the needed colors using the IA map
Cn      = C(IA);

% permute the surface triangulation points using IB map
Xbn     = Xb(IB,:);

% map the point numbers used in triangle definitions
% NOTE: for that you need inverse map
iIB(IB) = 1:length(IB);
trin    = iIB(tri);

trisurf(trin,Xbn(:,1),Xbn(:,2),Xbn(:,3),Cn,'EdgeAlpha',0,'FaceColor','interp');