Creating mesh from volumetric array

I am currently utilizing data taken from CT scans edited within MATLAB to output an obj mesh for potential augmented reality visualization.

For this, I have created an edited matrix which is a 276x276x8 uint8. I can utilize volumeViewer to show the data has been segmented and only showing the desired patient data, however I am struggling to take this array and input into any meshing format for output.

The current code looks like this:

folder2 = "C:\Users\Hidde\Documents\MATLAB\CT_Abdomen";

[V,sp] = dicomreadVolume(folder2);

V = squeeze(V);

V2 = im2uint8(V);

maxValue = max(V2,[],'all');

T = maxValue *(.80);

T2 = maxValue * (.4);

[a,b,c] = size(V2);

V3 = zeros(a,b,c);

for l = 1:a

for m = 1:b

for n = 1:c

if V2(l,m,n) <= T && V2(l,m,n) >= T2

V3(l,m,n) = V2(l,m,n);

else

V3(l,m,n) = 0;

end

end

end

end

V3 = im2uint8(V3);

V3 = V3(225:500,125:400,1:8);

volumeViewer(V3)

V3 is 276x276x8 uint8 and V4 is 276x276x8 double format.

Can someone help me convert this data into a format which works with mesh, meshgrid, or trimesh?

Shown is the volumeViewer of V3.

Answer 1

[Kevin Holly]

Here I made a meshgrid and then created a 3D scatter plot with the intensity displayed as a color.

% Randomly generate an volume
V = uint8(255*rand(276,276,8));

x=1:size(V,1);
y=1:size(V,2);
z=1:size(V,3);

[X,Y,Z]=meshgrid(x,y,z);

x=reshape(X,1,[]);
y=reshape(Y,1,[]);
z=reshape(Z,1,[]);

intensity = reshape(V,1,[]);
M = [x;y;z;intensity];

figure
scatter3(M(1,:),M(2,:),M(3,:),[],M(4,:),'filled')
colormap(bone)
colorbar

Answer 2

[Riccardo Scorretti]

Dear Alan,

perhaps for your purpose you could try to use the marching cube algorithm. A Matlab implementation is available here: https://fr.mathworks.com/matlabcentral/fileexchange/32506-marching-cubes

Given the kind of data you want to display, I sincerely think that this is your best option. An alternative could be to generate a mesh composed by tetrahedra (for instance by using iso2mesh: http://iso2mesh.sourceforge.net/cgi-bin/index.cgi), but this is a very complex procedure, and even if you manage to do it, I'm pretty sure it will not provide a better result.

Here you are a complete code which downloads a computational phantom from the web and render the blood. You have to download the function MarchingCube.m to run this code.

% Download a computational phantom (license GPL3)
buffer = webread('https://version.aalto.fi/gitlab/ilaakso/alvar/-/raw/master/Alvar_v16.mat?inline=false');
fid = fopen('alvar.mat', 'w');
fwrite(fid, buffer) ; clear buffer
fclose(fid);

% Load it and undersample it (the model is huge)
load alvar.mat
under = 4;
voxelData = voxelData(1:under:end, 1:under:end, 1:under:end);
x = x(1:under:end) ; y = y(1:under:end) ; z = z(1:under:end);



% Isolate the blood and artery (this is nasty, but I didn't manage to do otherwise)
buffer = voxelData;
buffer(buffer ~= 5  &  buffer ~= 8) = 0;

% Render it by using the marching cube algorithm
[x_, y_, z_] = meshgrid(single(y), single(x), single(z));   % *** mind the order of input args ***
[F, V, col] = MarchingCubes(x_, y_, z_, buffer, 1);
view(3) ; axis image ; camproj perspective

% The model can be also rendered in a separate window as:
figure
patch('vertices', V, 'faces', F, 'edgecolor', 'none', ...
  'facecolor', 'b');
view(3) ; axis image ; camproj perspective

On my PC this code produced the following figure:

Enjoy