Asked by bobby
on 15 Apr 2012

I have x, y axis values and corresponding z values. This allows me to plot the shape using surf function. Example:-

How to calculate surface area of this irregular shape?

Please i'm a novice learning this stuff, i have come along quite well since last week or two. Thank you for all your help. The other answers provided before for irregular shapes have me confused, since i don't know how to triangulate.

Additionally i also want to calculate surface area above a certain height. So i'll cut all the values below a certain height and then possibly calculate surface area above a cut-off height.

Edit:

Thanks to inputs by the community and basically whole programme written by Richard. Here's the code.

function sa(Z, cutoff)

%Credit: Richard Brown, MATLAB central forum (http://www.mathworks.com/matlabcentral/answers/)

dx=0.092; % x-axis calibration

dy=0.095; % y-axis calibration

[m, n] = size(Z);

areas = 0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(1:m-1,1:n-1))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(1:m-1,1:n-1))).^2) + ...

0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(2:m,2:n))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(2:m,2:n))).^2);

zMean = 0.25 * (Z(1:m-1,1:n-1) + Z(1:m-1,2:n) + Z(2:m,1:n-1) + Z(2:m,2:n));

areas(zMean <= cutoff) = 0;

surfaceArea = sum(areas(:));

sprintf('Total surface area is %2.4f\n', surfaceArea)

return

end

Answer by Richard Brown
on 16 Apr 2012

Edited by Richard Brown
on 15 Jul 2013

Accepted Answer

OK, well how about splitting each quadrilateral into two triangles, and just summing up the areas? I'm sorry there's no way I can make this look pleasant, but ...

[m, n] = size(Z);

areas = 0.5*sqrt((dx*dy)^2 + (dy*(Z(1:m-1,2:n) - Z(1:m-1,1:n-1))).^2 + ...

(dx*(Z(2:m,1:n-1 - Z(1:m-1,1:n-1)))).^2) + ...

0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(2:m,2:n))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(2:m,2:n))).^2)

surfaceArea = sum(areas(:))

edit: Jul 15, 2013 There was a mistake, and dx and dy were backwards in the first two lines of the code. The code has been corrected now.

bobby
on 16 Apr 2012

That is one awesome code Richard. I can't thank you enough. Here i was searching for triangulations and TINs since last thursday. It gives me result within 3% tolerance level, which is absolutely awesome. Just a small correction in the code, it is missing a bracket.

areas = 0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(1:m-1,1:n-1))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(1:m-1,1:n-1))).^2) + ...

0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(2:m,2:n))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(2:m,2:n))).^2);

surfaceArea = sum(areas(:))

bobby
on 16 Apr 2012

Richard Brown
on 16 Apr 2012

For each of the areas in the 'areas' array, you could work out the mean Z value

zMean = 0.25 * (Z(1:m-1,1:n-1) + Z(1:m-1,2:n) + Z(2:m,1:n-1) + Z(2:m,2:n))

and then

areas(zMean > cutoff) = 0;

surfaceArea = sum(areas(:));

It's a bit sloppy - you should probably do this for the individual triangles, but it will certainly give you a reasonable approximation

Sign in to comment.

Answer by Walter Roberson
on 15 Apr 2012

I do not know if you will be able to calculate the surface area of the regular shape interspersed with the irregular tops. When I look at the image with the irregular tops, it looks to me as if the tops are fairly irregular, possibly even fractal. I don't know if a meaningful surface area could be calculated: the surface area of a fractal is infinite.

For the first figure, I can think of a crude way to find the area, but I think there are simpler ways. I would need to think further about good ways to find the area. But to cross-check: are your x and y regularly spaced, or irregularly ?

bobby
on 16 Apr 2012

Sean de Wolski
on 16 Apr 2012

I wasn't commenting on Richard's script at all, but rather the fractal nature of surface area (and length) of irregular shapes. I apologize if that came across the wrong way.

http://www.google.com/#hl=en&gs_nf=1&cp=19&gs_id=1k&xhr=t&q=how+long+is+the+coast+of+britain&fp=46b06bf633a62cc4

bobby
on 16 Apr 2012

Sign in to comment.

Answer by Richard Brown
on 16 Apr 2012

If your data is smooth enough (assuming that's what you're after), then there is a really quick way to work out approximately the surface area, using the normals that Matlab computes when plotting the surface. To make it simple I'll assume you have a uniform grid with spacings dx and dy.

dA = dx * dy;

h = surf( ... )

N = get(h, 'VertexNormals');

N_z = squeeze(N(:, :, 3));

% Normalise it

N_z = N_z ./ sqrt(sum(N.^2, 3));

Area = dA * sum(1./N_z(:));

Richard Brown
on 16 Apr 2012

bobby
on 16 Apr 2012

Richard Brown
on 16 Apr 2012

Sign in to comment.

Answer by bobby
on 16 Apr 2012

function surfacearea(Z, cutoff)

dx=0.092;

dy=0.095;

[nrow, ncol] = size(Z);

for m=1:nrow

for n=1:ncol

zMean = mean(Z(1:m-1,1:n-1) + Z(1:m-1,2:n) + Z(2:m,1:n-1) + Z(2:m,2:n));

if(zMean >= cutoff)

areas = 0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(1:m-1,1:n-1))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(1:m-1,1:n-1))).^2) + ...

0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(2:m,2:n))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(2:m,2:n))).^2);

end

end

end

surfaceArea = sum(areas(:));

return

end

This is how the final code looks like with a cut-off limit included. Any inputs into speeding up this code? I have 1024*1280 matrix.

Richard Brown
on 17 Apr 2012

I think you may have slightly misunderstood some of my earlier answers - you don't need to have any loops at all. All you should need to do is:

areas = 0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(1:m-1,1:n-1))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(1:m-1,1:n-1))).^2) + ...

0.5*sqrt((dx*dy)^2 + (dx*(Z(1:m-1,2:n) - Z(2:m,2:n))).^2 + ...

(dy*(Z(2:m,1:n-1) - Z(2:m,2:n))).^2);

zMean = 0.25 * (Z(1:m-1,1:n-1) + Z(1:m-1,2:n) + Z(2:m,1:n-1) + Z(2:m,2:n))

areas(zMean > cutoff) = 0;

surfaceArea = sum(areas(:));

And that should run pretty fast

bobby
on 17 Apr 2012

Oh yes. This does not require loops. Let me run it again. Thanks Richard.

bobby
on 17 Apr 2012

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 1 Comment

## Richard Brown (view profile)

## Direct link to this comment

https://www.mathworks.com/matlabcentral/answers/35501-surface-area-from-a-z-matrix#comment_159835

Sign in to comment.