calculate volume under of 3D isosurface plot

I have a 3D matrix data V(201,101,101). I want to calculate the isosurface volume for a particular isovalue.

How to calculate the volume.

I used following code to generate isosurface.

figure;
p = patch(isosurface(x,y,z,V,-2));
isonormals(x,y,z,V,p);
cdata = smooth3(rand(size(V)),'box',7);
isocolors(x,y,z,cdata,p)
p.Facecolor = 'interp';
p.Edgecolor = 'none';

[darova]

alphaShape?

[Abhishek Saini]

I didin't understand, what do you mean by alphashape?

[Akira Agata]

I believe darova-san is reffering the following functions:

https://jp.mathworks.com/help/matlab/ref/alphashape.html

https://jp.mathworks.com/help/matlab/ref/alphashape.volume.html

[Abhishek Saini]

How to get the values of x,y, and z.

I used

p = issurface(x,y,z,V,isovalue).

It returns

p.vertices and

p.faces.

What should I use to calculate the colum and how does it will give volume.

I have attached isosurface image.

Answer 1

[Constantino Carlos Reyes-Aldasoro]

What isosurface and isovalue are doing are generating a surface around those points that are above a certain value. All code you showed is good for visualising, but what you really want here is not a visual display but a number of the volume.

If you want the volume that is the same as how many points satisfy your condition of voxel intensity > value, so what you have to do is the following:

volumeOfInterest = (V >-2 );

That will generate a matrix of same dimensions as V with 1 for those points above the isovalue and 0 below, then simply sum those values

volume = sum (volumeOfInterest(:));

And then, if needed you will need to callibrate that value to the physical quantity of your interest (say mm^3) with a conversion between voxels and that quantity.

Hope that helps

[Abhishek Saini]

Thanks Constantino for the detailed answer. If I am not wronf the volume is the sum of number of voxels having values 1 (lets say 200). This means it is scalar quantity. So the actual volume will be

volumeOfInterest = (V >-2 );
volume = sum (volumeOfInterest(:));
actual_volume = volume * volume_of_1_voxel; % unit of voxel i.e. mm^3

Is this correct.

[Walter Roberson]

Abhishek Saini: please only use Flag to notify the moderators of problems with posts.

[Constantino Carlos Reyes-Aldasoro]

True, the code I wrote was only calculating the number of voxels, but if you want to callibrate that number into a certain physical quantity, like m^3 or mm^3 you will need to do that multiplication. Will edit my previous response so that it has the complete answer in case someone new arrives to this post.