Perimeter of a implicit curve

I tried to make a program that returns the perimeter of a plane curve defined implicitly. The difficulty is that this type of curves can have more connected regions.

As an example we can take

f=@(x,y) sin(x).*sin(y)-0.5;
fp=fimplicit(f,[-8,8,-8,8]);

I observed that the points used by the fimplicit command are in such an order that makes the connected regions separately, and I detect the indices where we pass from one connected region to another. My idea was that if the segment connected two consecutive points from Points is greater than a small enough number (but not too small) tol than this is the step where a connected region is finished and we go to the next one. I'm not sure that this is indeed the way fimplicit command works. I see this picture that I make using the plot command:

plot(real(Points),imag(Points)) %see below what is the vector Points

Here is my code:

tol=0.01;
fp=fimplicit(f,[-8,8,-8,8],'meshdensity', 1000,'visible', 'off');
Points=fp.XData+1i.*fp.YData; %these are the x and y coordinates of the points used by fimplicit to represent the curve f(x,y)=0.
delete(fp);
Points=Points(isfinite(Pointsnew)==1);
d=diff(Points);
P=sum(hypot(real(d),imag(d)).*(hypot(real(d),imag(d))<tol)); %this is the perimeter

It returns very good values (3 or 4 decimals are exact). But it depends on how I choose the tolerance tol. If I take it too low the perimeter will be less than it is indeed. Moreover if I hace a part of a connected region on a corner then I want to take into account the part of the perimeter of the rectangle that bounds that region. With my program that is impossible.

Is there a way to write that program without the use of tol and to take into account the corners if needed?

[Star Strider]

First, ‘Pointsnew’ is nowhere to be found.

Second, the visible fimplicit plot has no connected lines.

Third, what perimeter do you want to calculate? Is it the perimeter of the entire region (in which case the boundary function would be appropriate to define it), or the perimeters of each or all of the closed contours (essentially trivial)?

[Maxim Bogdan]

I made a typo with Pointsnew;

If I delete NaN points with the command Points=Points(isfinite(Pointsnew)==1) and use plot afterwards I will get to the image that I posted above.

I want to calculate the perimeter indicated by black lines in the picture inserted by Matt J in his answer (so with the corners taken into account).

Answer 1

[Matt J]

Perhaps as follows (see also this thread)

plotRange=[-1,+1, -1,+1]*8;
plotCorners=[-1 1;1 1;1 -1; -1 -1]*8;

f=@(x,y) (sin(x).*sin(y)-0.5);
fp=fimplicit(f,2*plotRange,'MeshDensity',3000); 
XY=[fp.XData;fp.YData].'; close

shpRegions=polyshape(XY);

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

shpBounds=polyshape(plotCorners);

shpRegions=intersect(shpBounds,shpRegions);

Perimeter=perimeter(shpRegions)

Perimeter = 109.2290

plot(shpRegions,'FaceColor','r')

[Maxim Bogdan]

What is with that warning? There seem not to be duplicated points...

[Maxim Bogdan]

I do not know from what reason, when I put your code into a function it returns the result 107.8392 and gives me the following plot which do not include the corners:

[Maxim Bogdan]

Here is the function code:

function P=perimetru(f,a1,b1,a2,b2)
plotRange=[a1,b1,a2,b2];
plotCorners=[a1 b2;b1 b2;b1 a2; a1 a2];
fp=fimplicit(f,plotRange,'MeshDensity',3000); 
Points=[fp.XData;fp.YData].'; 
close
shpRegions=polyshape(Points);
shpBounds=polyshape(plotCorners);
shpRegions=intersect(shpBounds,shpRegions);
P=perimeter(shpRegions);
plot(shpRegions,'FaceColor','k')
end

[Matt J]

There are duplicate points:

size(XY)

ans = 1×2

39379 2

size(unique(XY,'rows','stable'))

ans = 1×2

39329 2

There are also NaNs:

any(isnan(XY(:)))

ans = logical
   1

[Maxim Bogdan]

If I remove the NaN points then all is crashed. I get the following figure:

[Maxim Bogdan]

Do you know a reason why the corners are not taken into account for the above writen function, although if I put in the command window your code it makes the figure with the good corners?

[Matt J]

Well, you don't seem to have included samples of the boundary rectangle in the list of points in any way.