MATLAB Answers

Estimate area inside point cloud

6 views (last 30 days)
Hi, I'm trying to determine the energy consumed during a hysteresis cycle. That should be the area locked inside a point cloud such as the one in the figure, where two 600 element vectors act as border.
I have tried trapz(H,B), what results in a negative solution, so I was wondering if there was a way of maybe defining the point cloud as a curve and integrating it.
All data comes from experimental results, so I don't have a function to use in integral2.

Accepted Answer

Star Strider
Star Strider on 1 May 2021
Edited: Star Strider on 1 May 2021
I have tried trapz(H,B), what results in a negative solution
The data are still reliable. The negative solution is the result of the direction that trapz integrates the data.
Example —
t = linspace(0, 2*pi, 1000);
x = 2*cos(t);
y = 2*sin(t);
A = trapz(x, y)
A = -12.5663
The area is still , so:
r = 2;
A = pi*r^2
A = 12.5664
And the sign here is irrelevant.
EDIT — (1 May 2021 at 18:00)
Another illustration —
t = linspace(2*pi, 0, 1000); % Reversed Direction Produces Positive Result
x = 2*cos(t);
y = 2*sin(t);
A = trapz(x, y)
A = 12.5663
.

More Answers (1)

John D'Errico
John D'Errico on 1 May 2021
As long as the curve is a sequence of points defining the perimeter and they are distinct, you can just compute the area using a tool like polyarea. Or you can create a polyshape, and then compute the area of that curve. Be careful though. You cannot apply any such tool (including trapz) to a randomly sequenced point cloud. They MUST be in sequence around the perimeter.
For example:
n = 100;
theta = linspace(0,2*pi,n);
x = cos(theta);
y = sin(theta);
plot(x,y,'o')
axis equal
Now compute the area.
polyarea(x,y)
ans = 3.1395
The area of a circle of unit radius would be pi, so roughly 3.14. This will be a slight underestimate, becuase the polygonal region lies entirely inside the circle. So 3.1395 is decent.
The nice thing is polyarea does not care if the sequnce traverses the circle in a counter-clockwise or clockwise manner. That is, the above sequence went clockwise around the circle. So this next area, which is traversed in a clockwise manner, is also positive.
polyarea(cos(flip(theta)),sin(flip(theta)))
ans = 3.1395
And we could have used a polyshape too.
ps = polyshape(x(1:end-1),y(1:end-1));
I dropped the last point from each of x and y there so polyshape would not get upset, since that point is essentially also the first point.
plot(ps)
axis equal
area(ps)
ans = 3.1395

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!