Thread Subject: Intersect circle and irregular shape

Hi, I am looking for ways to calculate the percentage of the
circumference of a circle enclosed within an intersecting irregular
polygon. For example, if I have a circle centered within a triangle, I
want to know how much of the circumference of the circle is contained
within the triangle. This could be anywhere between 0 and 100. I know
I can do this in a topological GIS, but I really need Matlab for this.
I want to use this to build an edge correction program. I couldn't
find anything in the Matlab file central or anywhere else.

Thank you for any help you may provide,

GL

glaroc <glaroc@gmail.com> wrote in message <66c9ae27-c839-41e2-
a9a3-d5dcc7f311a7@i12g2000prf.googlegroups.com>...
> Hi, I am looking for ways to calculate the percentage of the
> circumference of a circle enclosed within an intersecting irregular
> polygon. For example, if I have a circle centered within a triangle, I
> want to know how much of the circumference of the circle is contained
> within the triangle. This could be anywhere between 0 and 100. I know
> I can do this in a topological GIS, but I really need Matlab for this.
> I want to use this to build an edge correction program. I couldn't
> find anything in the Matlab file central or anywhere else.
>
> Thank you for any help you may provide,
>
> GL
--------
  Here's an outline of how I would tackle your problem, GL. First, find all the
points of intersection between the circle and the various line segments of the
polygon. Then sort such intersection points in accordance with the angles
they make with the circle's center. These points now divide up the circle into
non-overlapping arcs. Find the arc "midpoint" on the circle between each
pair of successive intersection points. Finally, use the function 'inpolygon' to
determine which of these midpoints lie inside the polygon. Each of the entire
circular arcs containing such midpoints must also lie within the polygon.
Then the sum of these "inside" arcs' angles divided by 2*pi times 100 gives
you the percentage you desire. I leave you to deal with the special case when
there are no intersection points.

  Here is a procedure for finding the intersection between a circle having
radius r and center (a,b) with a line segment whose endpoints are (c,d) and
(e,f). Any point on the infinite line through these points has the parametric
representation

 x = (c+e)/2+(c-e)/2*t
 y = (d+f)/2+(d-f)/2*t

with parameter t. If t lies between -1 and +1, then such a point is actually
within the line segment. To also lie on the circle, we must have

 (x-a)^2+(y-b)^2 = r^2

By substituting in the parametric expressions in t for x and y in this last
equality, we get a quadratic equation in t which we can use matlab's 'roots'
function to solve. This done as follows:

 t = roots([(c-e)^2+(d-f)^2, ...
            2*((c-e)*(c+e-2*a)+(d-f)*(d+f-2*b)), ...
            (c+e-2*a)^2+(d+f-2*b)^2-4*r^2]);

where the three expressions show the derived coefficients of the above
quadratic. There will be either two, one, or no distinct real roots to this,
according as the circle intersects that many points of the infinite line. As
mentioned above, those that have a real t value between -1 and +1 yield the
points of intersection with the line segment, and their coordinates are
determined by substituting t into the above formulas for x and y.

Roger Stafford

Thanks a lot for your detailed response. However, tell me if I am
wrong, this type of approach would be pretty difficult to apply in a
situation like this:

http://bayimg.com/cAimkAAbn

In an ideal world, I would need something that could handle even more
complicated polygons.

GL

Roger Stafford wrote:
> glaroc <glaroc@gmail.com> wrote in message <66c9ae27-c839-41e2-
> a9a3-d5dcc7f311a7@i12g2000prf.googlegroups.com>...
> > Hi, I am looking for ways to calculate the percentage of the
> > circumference of a circle enclosed within an intersecting irregular
> > polygon. For example, if I have a circle centered within a triangle, I
> > want to know how much of the circumference of the circle is contained
> > within the triangle. This could be anywhere between 0 and 100. I know
> > I can do this in a topological GIS, but I really need Matlab for this.
> > I want to use this to build an edge correction program. I couldn't
> > find anything in the Matlab file central or anywhere else.
> >
> > Thank you for any help you may provide,
> >
> > GL
> --------
> Here's an outline of how I would tackle your problem, GL. First, find all the
> points of intersection between the circle and the various line segments of the
> polygon. Then sort such intersection points in accordance with the angles
> they make with the circle's center. These points now divide up the circle into
> non-overlapping arcs. Find the arc "midpoint" on the circle between each
> pair of successive intersection points. Finally, use the function 'inpolygon' to
> determine which of these midpoints lie inside the polygon. Each of the entire
> circular arcs containing such midpoints must also lie within the polygon.
> Then the sum of these "inside" arcs' angles divided by 2*pi times 100 gives
> you the percentage you desire. I leave you to deal with the special case when
> there are no intersection points.
>
> Here is a procedure for finding the intersection between a circle having
> radius r and center (a,b) with a line segment whose endpoints are (c,d) and
> (e,f). Any point on the infinite line through these points has the parametric
> representation
>
> x = (c+e)/2+(c-e)/2*t
> y = (d+f)/2+(d-f)/2*t
>
> with parameter t. If t lies between -1 and +1, then such a point is actually
> within the line segment. To also lie on the circle, we must have
>
> (x-a)^2+(y-b)^2 = r^2
>
> By substituting in the parametric expressions in t for x and y in this last
> equality, we get a quadratic equation in t which we can use matlab's 'roots'
> function to solve. This done as follows:
>
> t = roots([(c-e)^2+(d-f)^2, ...
> 2*((c-e)*(c+e-2*a)+(d-f)*(d+f-2*b)), ...
> (c+e-2*a)^2+(d+f-2*b)^2-4*r^2]);
>
> where the three expressions show the derived coefficients of the above
> quadratic. There will be either two, one, or no distinct real roots to this,
> according as the circle intersects that many points of the infinite line. As
> mentioned above, those that have a real t value between -1 and +1 yield the
> points of intersection with the line segment, and their coordinates are
> determined by substituting t into the above formulas for x and y.
>
> Roger Stafford

glaroc <glaroc@gmail.com> wrote in message <cfe96c89-937d-4cd8-be63-
ed7ceb3607fa@k39g2000hsf.googlegroups.com>...
> Thanks a lot for your detailed response. However, tell me if I am
> wrong, this type of approach would be pretty difficult to apply in a
> situation like this:
>
> http://bayimg.com/cAimkAAbn
>
> In an ideal world, I would need something that could handle even more
> complicated polygons.
>
> GL
-------
  No, that polygon/circle situation should come out fine with the method I
described. The bottom polygon segment should produce two real t values
lying between -1 and +1. The next segment on the southeast side would
have two real t's but only one in the above range. The next segment beyond
that would have a similar result. The remaining three segments would have
only complex-valued t's. So that leaves you with four actual points of
intersection. The midpoint of the arc betweent the first two will lie outside
the polygon, the next arc midpoint is inside, the third arc midpoint is outside,
and finally the arc midpoint between the fourth point and the first is again
inside the polygon. So you add the arc lengths of the second and fourth of
these arcs and divide by 2*pi to get your percentage.

Roger Stafford

glaroc <glaroc@gmail.com> wrote in message <66c9ae27-c839-
41e2-a9a3-d5dcc7f311a7@i12g2000prf.googlegroups.com>...
> Hi, I am looking for ways to calculate the percentage of
the
> circumference of a circle enclosed within an intersecting
irregular
> polygon. For example, if I have a circle centered within
a triangle, I
> want to know how much of the circumference of the circle
is contained
> within the triangle. This could be anywhere between 0 and
100. I know
> I can do this in a topological GIS, but I really need
Matlab for this.
> I want to use this to build an edge correction program. I
couldn't
> find anything in the Matlab file central or anywhere else.
>
> Thank you for any help you may provide,
>
> GL

GL,

I would take 100 or 1000 or 10 000 or 10^k evenly spaced
points on the circle and calculate the percentage that lie
inside the figure. If the figure is quite irregular, then
perhaps choosing 10^k randomly distibuted points and
running the code many times would be best.

This is only an approximation, of course, but it should be
easy to code, fast and give good results.

Charles

"Charles Cuell" <cuell@math.usask.ca> wrote in message <fm33v1$hkn
$1@fred.mathworks.com>...
> I would take 100 or 1000 or 10 000 or 10^k evenly spaced
> points on the circle and calculate the percentage that lie
> inside the figure. If the figure is quite irregular, then
> perhaps choosing 10^k randomly distibuted points and
> running the code many times would be best.
>
> This is only an approximation, of course, but it should be
> easy to code, fast and give good results.
>
> Charles
-------
  Charles, I grant you this would be easier to code, but if the 10^k value
Glaroc choses to select is very large so as to be accurate, this would mean
that same number of calls on inpolygon. It might not execute as fast as he
would like, as compared with only one call per arc midpoint using the
intersection method.

Roger Stafford

In fact, I already tried this approximation method. The problem is
that it is quite slow since I have to do this for thousands of circle.
This involves a lot of inpolygon calculations. Roger's method seems
longer to code, but will likely be a lot faster.

I will give it a try when I have a minute.

GL
>
> $...@fred.mathworks.com>...> I would take 100 or 1000 or 10 000 or 10^k evenly spaced
> > points on the circle and calculate the percentage that lie
> > inside the figure. If the figure is quite irregular, then
> > perhaps choosing 10^k randomly distibuted points and
> > running the code many times would be best.
>
> > This is only an approximation, of course, but it should be
> > easy to code, fast and give good results.
>
> > Charles
>
> -------
> Charles, I grant you this would be easier to code, but if the 10^k value
> Glaroc choses to select is very large so as to be accurate, this would mean
> that same number of calls on inpolygon. It might not execute as fast as he
> would like, as compared with only one call per arc midpoint using the
> intersection method.
>
> Roger Stafford

glaroc <glaroc@gmail.com> wrote in message <1c0f9a46-130e-
4bd6-9918-8238171a5977@k2g2000hse.googlegroups.com>...
> In fact, I already tried this approximation method. The
problem is
> that it is quite slow since I have to do this for
thousands of circle.
> This involves a lot of inpolygon calculations. Roger's
method seems
> longer to code, but will likely be a lot faster.
>
> I will give it a try when I have a minute.
>
> GL
> >
> > $...@fred.mathworks.com>...> I would take 100 or 1000
or 10 000 or 10^k evenly spaced
> > > points on the circle and calculate the percentage
that lie
> > > inside the figure. If the figure is quite irregular,
then
> > > perhaps choosing 10^k randomly distibuted points and
> > > running the code many times would be best.
> >
> > > This is only an approximation, of course, but it
should be
> > > easy to code, fast and give good results.
> >
> > > Charles
> >
> > -------
> > Charles, I grant you this would be easier to code,
but if the 10^k value
> > Glaroc choses to select is very large so as to be
accurate, this would mean
> > that same number of calls on inpolygon. It might not
execute as fast as he
> > would like, as compared with only one call per arc
midpoint using the
> > intersection method.
> >
> > Roger Stafford
>

I would just go for coffee, or lunch, or a vacation :)
Point taken, anyway.

Charles

"Charles Cuell" <cuell@math.usask.ca> wrote in message <fm36sq$2m2
$1@fred.mathworks.com>...
> I would just go for coffee, or lunch, or a vacation :)
> Point taken, anyway.
> Charles
--------
  It would be a very long coffee break, Charles, if GL were to attempt to
extract all the accuracy out of his computations that matlab's double
precision is capable of, using the 10^k method. It would require something
of the order of 10^16 divisions of the full circle and therefore a like number
of calls on inpolygon to achieve that. That way I figure he might have his
computer occupied for several hundred years on each circle. With the web
example GL showed, there would only be four calls on inpolygon and six calls
on roots with the intersection method. It seems to me a better practice, both
from a computation efficiency point of view, and from an esthetic one too, to
find all the intersections between the circle and the polygon precisely and
work from there. The computations should all be accurate to the full one part
in 10^16 that way.

Roger Stafford

> I couldn't find anything in the Matlab file central or
anywhere else.

I think this is 99.99[inf]% solution to your problem.

http://www.mathworks.com/matlabcentral/fileexchange/
loadFile.do?objectId=18173&objectType=file

J. Luis

"Roger Stafford"
<ellieandrogerxyzzy@mindspring.com.invalid> wrote in
message <fm3b02$giu$1@fred.mathworks.com>...
> "Charles Cuell" <cuell@math.usask.ca> wrote in message
<fm36sq$2m2
> $1@fred.mathworks.com>...
> > I would just go for coffee, or lunch, or a vacation :)
> > Point taken, anyway.
> > Charles
> --------
> It would be a very long coffee break, Charles, if GL
were to attempt to
> extract all the accuracy out of his computations that
matlab's double
> precision is capable of, using the 10^k method. It would
require something
> of the order of 10^16 divisions of the full circle and
therefore a like number
> of calls on inpolygon to achieve that. That way I figure
he might have his
> computer occupied for several hundred years on each
circle. With the web
> example GL showed, there would only be four calls on
inpolygon and six calls
> on roots with the intersection method. It seems to me a
better practice, both
> from a computation efficiency point of view, and from an
esthetic one too, to
> find all the intersections between the circle and the
polygon precisely and
> work from there. The computations should all be accurate
to the full one part
> in 10^16 that way.
>
> Roger Stafford
>
>

Ouch. Even with a government job, that's a long coffee
break. Thanks for the info.

Charles.

> > I couldn't find anything in the Matlab file central or
>
> anywhere else.
>
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis

Thanks! I was looking for this just before the new year, so I guess
its normal I didn't see that one...

> > I couldn't find anything in the Matlab file central or
>
> anywhere else.
>
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis

Thanks! I was looking for this just before the new year, so I guess
its normal I didn't see that one...

> > I couldn't find anything in the Matlab file central or
>
> anywhere else.
>
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis

Thanks! I was looking for this just before the new year, so I guess
its normal I didn't see that one...

> > I couldn't find anything in the Matlab file central or
>
> anywhere else.
>
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis

Thanks! I was looking for this just before the new year, so I guess
its normal I didn't see that one...

"Joaquim Luis" <jluis@--ualg--.pt> wrote in message <fm3cmg$r8q
$1@fred.mathworks.com>...
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis
-----------
  Joaquim, if I interpret Guillaume Jacquenot's m-file for the
"Polygons_intersection" function on the Matlab Central File Exchange
correctly, he has to approximate curved shapes like circles with discrete
polygons. His example has the following two lines:

 S(1).P(1).x = [1.5 1.5 0.5 1.5+cos(pi:pi/30:3*pi/2)];
 S(1).P(1).y = [0.5 1.5 1.5 1.5+sin(pi:pi/30:3*pi/2)];

which would contain a quarter circle as approximated by 16 polygonal
segments.

  The problem of finding the intersection between a line segment and a circle
is very straightforward and just involves the solution to a simple quadratic
equation, so I don't see why GL can't do it that way directly, rather than mess
around with such approximations.

Roger Stafford

> > I couldn't find anything in the Matlab file central or
>
> anywhere else.
>
> I think this is 99.99[inf]% solution to your problem.
>
> http://www.mathworks.com/matlabcentral/fileexchange/
> loadFile.do?objectId=18173&objectType=file
>
> J. Luis

Thanks! I was looking for this just before the new year, so I guess
its normal I didn't see that one...

"Roger
Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote
in message <fm3igl$kmn$1@fred.mathworks.com>...

> -----------
> Joaquim, if I interpret Guillaume Jacquenot's m-file
for the
> "Polygons_intersection" function on the Matlab Central
File Exchange
> correctly, he has to approximate curved shapes like
circles with discrete
> polygons. His example has the following two lines:
>
> S(1).P(1).x = [1.5 1.5 0.5 1.5+cos(pi:pi/30:3*pi/2)];
> S(1).P(1).y = [0.5 1.5 1.5 1.5+sin(pi:pi/30:3*pi/2)];
>
> which would contain a quarter circle as approximated by
16 polygonal
> segments.
>
> The problem of finding the intersection between a line
segment and a circle
> is very straightforward and just involves the solution to
a simple quadratic
> equation, so I don't see why GL can't do it that way
directly, rather than mess
> around with such approximations.

Roger,
Yes, that's right but in his first message he mentioned

"I am looking for ways to calculate the percentage of the
circumference of a circle enclosed within an intersecting
irregular polygon."

So I guess that the clipping polygons solution is, at
least, easier to use.

Joaquim Luis