Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Polybool gives subtraction with wrong vertices

Subject: Polybool gives subtraction with wrong vertices

From: Kenan

Date: 1 May, 2013 22:57:09

Message: 1 of 4

Hello all,

I am using the function polybool from the Mapping Toolbox and my version of Matlab is R2013a.

I have problems when I perform a simple subtraction of two polygonal regions (a rectangle described by [x1 y1] and a triangle described by [x0 y0]). See attached code and figures.

px = .7;
py = .75;
lx = (5*2^.5-6)/20;
ly = (10*3^.5-9)/60;

[x1,y1] = poly2cw([px; px+lx; px+lx; px; px],[py; py; py+ly; py+ly; py]); % rectangle
[x0,y0] = poly2cw([15;15;14;15]/19,[15;16;16;15]/19); % triangle
[x2,y2] = polybool('intersection',x1,y1,x0,y0); % intersection

[x3,y3] = polybool('subtraction',x0,y0,x2,y2); % soustraction

figure; plot(x0,y0,'b--',x1,y1,'k:',x3,y3,'r+-'); hold on
for i=1:length(x3)
    text(x3(i)*1.003,y3(i)*1.003,num2str(i))
end

As you can see, there is an extra vertex that should not be there (vertex 2-3). So instead of being 1-2,2-3,3-4,4-5,5-6; the subtraction should be 1-3,3-4,4-5,5-6. (6 being equal to 1).


This wouldn't be such a bad issue if it was all, but it causes the rest of the algorithm to behave abnormally (after performing a constrained Delaunay triangulation, the function isInterior gives an incorrect result)

First I check if the subtraction [x3 y3] is a connected polygon with ~isShapeMultipart(x3,y3): it is.

If I perform a non constrained Delaunay triangulation, the result contains triangles outside of [x3 y3] (in the general case) and there is no direct way of finding which ones. Hence it is easier to perform a constrained Delaunay triangulation. If I understand the word constrain correctly, Matlab performs the triangulation making sure the vertices in the constrain are vertices of the triangulation. However, it does not mean that these vertices will be the only vertices of the triangulation. I need someone to confirm this statement, or explain to me what is a constrain exactly. So I construct the constraint matrix for the Delaunay triangulation and perform the triangulation:

C3 = zeros(length(x3)-1,2);
    for j=1:2
        for i2=1:length(x3)-1
            C3(i2,j) = i2+j-1;
        end
    end
C3(end,end)=1;
DT31 = delaunayTriangulation(x3,y3,C3);

Finally, the statement inside31 = isInterior(DT31) gives 0 0 0 meaning that all 3 triangles in DT31 are outside of the constraint C3.

The code for constructing the constrain matrix C3 could be adapted to the weird result of polybool but I would lose robustness in my code. Plus, I would rather correct the weird result the earliest possible (ie just after [x3,y3] = polybool('subtraction',x0,y0,x2,y2); % soustraction) rather than making corrections afterwards.

Thank you for your precious help and/or remarks.

See here for illustrations: http://stackoverflow.com/questions/16326810/polybool-gives-subtraction-with-wrong-vertices

Subject: Polybool gives subtraction with wrong vertices

From: Kenan

Date: 7 May, 2013 18:57:09

Message: 2 of 4

NOTE: I mistook 'vertex' with 'edge'.

I received an answer from the technical support: "the algorithms are kind of sensitive to these rounding errors. I do not know whether this can be improved, but I will inform development, maybe they can make polybool more robust against these problems."

Subject: Polybool gives subtraction with wrong vertices

From: Amy Haskins

Date: 7 May, 2013 22:13:10

Message: 3 of 4

Hi Kenan,

In regards to your question, it looks like you're running into some rounding errors here. If you subtract the rectangle from the large triangle directly, you will get the result you would expect:

[x3,y3] = polybool('subtraction',x0,y0,x1,y1); % subtraction

When you do this as two steps, there is a very small amount of rounding error in the calculation which causes the extra line in your results. If we shift the bottom vertex down even slightly, the correct answer is given:

y2([1,4]) = y2([1,4]) - .000001;
[x3,y3] = polybool('subtraction',x0,y0,x2,y2); % subtraction

In a world of infinite percision, the two-step approach would work just fine. However, since we are dealing with floating point, there is about a 50/50 chance that each edges of the smaller triangle will be ever so slightly inside the larger one.

This is just a natural limitation of doing computation geometry with floating point numbers, the results from polybool are not surprising.

Amy

"Kenan" wrote in message <kls6k5$osr$1@newscl01ah.mathworks.com>...
>
> I am using the function polybool from the Mapping Toolbox and my version of Matlab is R2013a.
>
> I have problems when I perform a simple subtraction of two polygonal regions (a rectangle described by [x1 y1] and a triangle described by [x0 y0]). See attached code and figures.
>
> px = .7;
> py = .75;
> lx = (5*2^.5-6)/20;
> ly = (10*3^.5-9)/60;
>
> [x1,y1] = poly2cw([px; px+lx; px+lx; px; px],[py; py; py+ly; py+ly; py]); % rectangle
> [x0,y0] = poly2cw([15;15;14;15]/19,[15;16;16;15]/19); % triangle
> [x2,y2] = polybool('intersection',x1,y1,x0,y0); % intersection
>
> [x3,y3] = polybool('subtraction',x0,y0,x2,y2); % soustraction
>
> figure; plot(x0,y0,'b--',x1,y1,'k:',x3,y3,'r+-'); hold on
> for i=1:length(x3)
> text(x3(i)*1.003,y3(i)*1.003,num2str(i))
> end
>
> As you can see, there is an extra vertex that should not be there (vertex 2-3). So instead of being 1-2,2-3,3-4,4-5,5-6; the subtraction should be 1-3,3-4,4-5,5-6. (6 being equal to 1).
>
>

Subject: Polybool gives subtraction with wrong vertices

From: Kenan

Date: 8 May, 2013 18:54:08

Message: 4 of 4

Thank you Amy. Yes, I had noticed I was obtaining the correct result when using "triangle MINUS rectangle" instead of "triangle MINUS (rectangle INTER triangle)".
However, I was curious as to why there was a difference between the two methods (that are on paper the same). I understand it comes from rounding errors. However, correct me if I'm wrong, but there is no way of knowing *for sure* that the direct method will never fail either, right?

Talking about computational geometry and rounding errors, do you think I will be luckier if I use functions from the PDE Toolbox (initmesh instead of delaunayTriangulation for instance), or the rounding errors would be the same?

Thank you again.

Kenan

"Amy Haskins" <amy.haskins.nospam@mathworks.com> wrote in message <kmbu9m$q70$1@newscl01ah.mathworks.com>...
> Hi Kenan,
>
> In regards to your question, it looks like you're running into some rounding errors here. If you subtract the rectangle from the large triangle directly, you will get the result you would expect:
>
> [x3,y3] = polybool('subtraction',x0,y0,x1,y1); % subtraction
>
> When you do this as two steps, there is a very small amount of rounding error in the calculation which causes the extra line in your results. If we shift the bottom vertex down even slightly, the correct answer is given:
>
> y2([1,4]) = y2([1,4]) - .000001;
> [x3,y3] = polybool('subtraction',x0,y0,x2,y2); % subtraction
>
> In a world of infinite percision, the two-step approach would work just fine. However, since we are dealing with floating point, there is about a 50/50 chance that each edges of the smaller triangle will be ever so slightly inside the larger one.
>
> This is just a natural limitation of doing computation geometry with floating point numbers, the results from polybool are not surprising.
>
> Amy
>
> "Kenan" wrote in message <kls6k5$osr$1@newscl01ah.mathworks.com>...
> >
> > I am using the function polybool from the Mapping Toolbox and my version of Matlab is R2013a.
> >
> > I have problems when I perform a simple subtraction of two polygonal regions (a rectangle described by [x1 y1] and a triangle described by [x0 y0]). See attached code and figures.
> >
> > px = .7;
> > py = .75;
> > lx = (5*2^.5-6)/20;
> > ly = (10*3^.5-9)/60;
> >
> > [x1,y1] = poly2cw([px; px+lx; px+lx; px; px],[py; py; py+ly; py+ly; py]); % rectangle
> > [x0,y0] = poly2cw([15;15;14;15]/19,[15;16;16;15]/19); % triangle
> > [x2,y2] = polybool('intersection',x1,y1,x0,y0); % intersection
> >
> > [x3,y3] = polybool('subtraction',x0,y0,x2,y2); % soustraction
> >
> > figure; plot(x0,y0,'b--',x1,y1,'k:',x3,y3,'r+-'); hold on
> > for i=1:length(x3)
> > text(x3(i)*1.003,y3(i)*1.003,num2str(i))
> > end
> >
> > As you can see, there is an extra vertex that should not be there (vertex 2-3). So instead of being 1-2,2-3,3-4,4-5,5-6; the subtraction should be 1-3,3-4,4-5,5-6. (6 being equal to 1).
> >
> >

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us