dp@cfa.harvard.edu (Dave Plummer) (04/26/91)
I am looking for an algorithm (C or Fortran code would be nice) that finds the area of intersection between a circle and a square of arbitrary position and size. After thinking about it for a few minutes I realized that this is more complicated than it sounds. The intersecting area will vary from zero (no overlap) to the full area of the circle or square (one enclosed by the other). Please post any suggestions, references or code or e-mail them to: dp@cfa.harvard.edu - Dave Plummer
mgobbi@cs.ubc.ca (Mike Gobbi) (04/26/91)
In article <1991Apr25.192400.8438@cfa203.harvard.edu> dp@cfa.harvard.edu (Dave Plummer) writes: >I am looking for an algorithm (C or Fortran code would be nice) that >finds the area of intersection between a circle and a square of arbitrary >position and size. After thinking about it for a few minutes I realized >that this is more complicated than it sounds. > >The intersecting area will vary from zero (no overlap) to the full area of >the circle or square (one enclosed by the other). > >Please post any suggestions, references or code or e-mail them to: >dp@cfa.harvard.edu > >- Dave Plummer Here's a thought (I presume that you just want a number -- not a polygon): Start with the area of the circle. Call it A. Take the "top" edge of the square and extend it into a line. This line will either (a) pass "below" the circle (tangency counts as a miss) (b) pass "above" the circle (ditto) (c) intersect the circle at exactly two points, P & Q In case (a) A := 0 (b) A := A (c) A := A - (pie + wedge) where conect P & Q to center of circle (point C) pie = area of pie-shape bounded by circumference, PC, QC and lying "above" the line. (this is an easily-computed fraction of the circle) wedge = area of triangle PCQ repeat in a similar fashion for the other three sides of the square, possibly reducing the value of A every time. Note that this algorithm should work on a circle & any convex polygon, not just squares. As I said earlier, it won't give you the closed curve defined, just the area (though it should be fairly simple to modify it to save the points of intersection and determine which path (curved or straight) is the bound). -- __ /..\ In quest of knowledge... --mm--mm-- Mike Gobbi
davis@bedlam.asd.sgi.com (Tom Davis) (04/26/91)
If you just need a numerical approximation, why not approximate the
circle with an n-sided regular polygon, and use standard clipping
procedures to clip that polygon against the square. Then if the (clipped)
vertices are (x[0], y[0]), (x[1], y[1]), ... (x[m-1], y[m-1]), the area is
given by:
1/2 * sum { x[i](y[i+1] - y[i]) - y[i](x[i+1] - x[i])}
The sum goes from i = 0 to i = m-1, and assume x[m] = x[0] and y[m] = y[0].
Pick a big enough n, and you ought to do OK.lambert@silver.cs.umanitoba.ca (Tim Lambert) (04/26/91)
>>>>> On 25 Apr 91 19:24:00 GMT, dp@cfa.harvard.edu (Dave Plummer) said: > I am looking for an algorithm that finds the area of intersection > between a circle and a square of arbitrary position and size. All right, I'll have a go at it. You could modify one of the standard polygon clipping algorithms to work on circles, to give you the perimeter of the resulting shape and compute the area from that, but a square and a circle seem simple enough that something just for this case would be nicer. First, let's rotate, scale and translate everything so that the circle has centre (0,0) and radius 1. Let the rectangle have corners (l,b) (left bottom) and (r,t) (right top). Now cut the rectangle up into four triangular pieces by connecting the four corners to the origin. (If the origin is not in the square, some of these pieces may be negative, but that's OK.) We'll add up the areas (possibly negative) of the intersection between these triangles and the circle. I'll describe how to work out the area of the intersection between the circle and the bottom triangle (corners (0,0) (l,b) (r,b)). The area of this triangle is -0.5*(r-l)*b. (Note that if b is positive this area is _negative_) Now let's look at where the line y=b intersects the circle. (Drawing a picture for each following case will probably help understand what is going on.) Case 1: b>1 or b<-1 -- no intersection between line and circle. Both (l,b) and (r,b) are outside the circle. The triangle-circle intersection is just a sector of the circle between angles atan2(b,l) and atan2(b,r). Area is atan2(b,r) - atan2(b,l) Note1: if b is positive, this is negative, which is correct. Note2: you could also calculate this area as acos(b^2+r*l/((sqrt(b^2+r^2)*sqrt(b^2+l^2)))) (and then fix the sign) Case 2: -1<=b<=1 -- intersection at (A,b) and (-A,b) where A=sqrt(1-b^2) Now -A<=A and l<=r so there are six possible orderings of -A,A,l,r Case 2a: -A<=l<=r<=A The triangle is entirely inside the circle. Area is -0.5*(r-l)*b Case 2b: A<=l<=r (l,b) and (r,b) are outside the circle. Intersection is a sector as in case 1. Case 2c: l<=r<=-A Same as 2b. Case 2d. -A<=l<=A<=r The intersection between the triangle and circle can be cut into two pieces: the triangle (0,0) (l,b) (A,b) and the sector between the angles atan2(b,A) and atan2(b,l). Area is -0.5*(A-l)*b + atan2(b,r) - atan2(b,A) Case 2e. l<=-A<=r<=A Similar to 2d Area is -0.5*(r+A)*b + atan2(b,-A) - atan2(b,l) Case 2f. l<=-A<=A<=r This is the complicated one. We have to cut it up into two sectors and a triangle. Area is atan2(b,-A) - atan2(b,l) -0.5*(A+A)*b + atan2(b,r) - atan2(b,A) Now code that all up into a function tri_circ(r,l,b) that computes the area between the bottom triangle and the circle. The area of the circle rectangle intersection is: tri_circ(r,l,b) + tri_circ(t,b,r) + tri_circ(l,r,t) + tri_circ(b,t,l) (I just rotate the picture through 90 degrees each time) Well, almost. I assumed r>=l in the above analysis so I'd better reflect about the y axis (which will multiply the area by -1) as well in the last two cases. So the true answer is: tri_circ(r,l,b) + tri_circ(t,b,r) - tri_circ(r,l,t) - tri_circ(t,b,l) If the intersection is empty all the positive and negative areas above will cancel out and the area will be zero (or very close, given floating point arithmetic). This is a good test that you got all the signs right. Tim P.S. Alert readers might be concerned that if (c,d) is in the third quadrant and (e,f) in the second, then atan2(d,c) - atan2(f,e) doesn't give the angle correctly (it's 2*pi too small). Fortunately this can't happen in the above algorithm, since we only compute angles between points on the same horizontal line.
lambert@silver.cs.umanitoba.ca (Tim Lambert) (04/26/91)
>>>>> On Thu, 25 Apr 91 20:39:16 GMT, mgobbi@cs.ubc.ca (Mike Gobbi) said: > In article <1991Apr25.192400.8438@cfa203.harvard.edu> dp@cfa.harvard.edu (Dave Plummer) writes: >>I am looking for an algorithm that finds the area of intersection >>between a circle and a square of arbitrary position and size. Let's call the piece of a circle formed by cutting a piece of along a line a "segment" (is there an official name for this?) I've edited down Mike's method. > Start with the area of the circle. Call it A. > Take the "top" edge of the square and extend it into a line. This line > will either > (c) intersect the circle at exactly two points, P & Q > In case (c) A := A - (segment) > repeat in a similar fashion for the other three sides of the square, > possibly reducing the value of A every time. The trouble with this is that after cutting off one piece what you have left is no longer a circle, and if on the next step you cut off a segment that intersects one that you already cut off, you will subtract the area of the intersection between the two segments twice. Hmmm.... Maybe you could fix this by adding the area of the intersection of two segents back on. You can easily detect this situation, because it only happens if you have a corner of the square inside the circle. Tim