Intersection of two circles

[Pete]

I should clarify, by Mr. Bill's invention, I was referencing Bill Gates, as the inventor of the computer. You probably got that, but just to be sure. Yeah, yeah, yeah, it's a stretch, he only bought an operating system, but Wosniack's first name isn't Bill. Anyway, SNL 1970's skits aside, so what are the practical implications? Generally, I live by two rules and a saying. The rules are, Rule #1 don't do anything without having an idea. Rule #2, makes sure the idea is a good one, and the saying: "Just because you can doesn't mean you should." So if this is a boat load of un-fun work, I hope it proves useful. I mean a lot of math is impressive onto itself, but when it can be used for practical purposes, it's priceless.

Pete

[SierraKen]

No pressure Static, whatever you want to do of course. I was just playing around with it.

[bplus]

Any method that can be reduced or duplicated by printing the two circles on two pages and holding the pages up to the light to look for intersections is... cute, but I hope we don't need to explain the fallibility of these methods. At this point I have a choice: (i) either spruce up the sample code so it works nicely in all quadrants, or (ii) write up the derivation so then anyone can do it properly. Choices choices... where's my bong?

Just use bplus tried and true code taken from some C algo thing or prove it fails at some point.  It is only Analytic Geometry, pre-Calculus stuff for goodness sake. Wait... I should have notes somewhere from July 3....,

update: oh that's the translation from SmallBASIC, dang this is going to take some digging...

[bplus]

From 2015-09-12 2 pages of SmallBASIC code notes say Stackoverflow:

Yes this looks like my notes: https://stackoverflow.com/questions/3349125/circle-circle-intersection-points

First answer (of 5) written by Paul Bourke.

[Qwerkey]

So, I've had an independent go at it (even though, I'm sure, bplus has got the maths completely covered).

In the diagram (intersecting circles.jpg), we need to find the positions of C and D.

So, we need to consider triangles ABC and ABD.  I'm using just the centres outside each other condition in this case.

Firstly, the triangle ABD is just a reflection of ABC about the line AB, so only do ABC.  Secondly, the whole thing can be treated as if AB is the x- axis, and then do the (trivial) rotational translation conversion - STxAxTIC & bplus will know this and can easily be found in Wiki.

So consider the triangle ABC at the bottom of the diagram.  The distance between A and B = x1 + x2, and the distance from line AB to C is y.  We need to calculate x1 (and x2) and y.  We know AB, R1 and R2.

So, there is a formula to get triangle angles if all three sides are known.  See the (triangle angles.jpg) document, from Wiki (see the URL).

Thus, in our case alpha = arccos((R1^2+(x1+x2)^2-R2^2)/2*R1*(x1+x2))

When we know alpha, we know x1 and y, and therefore the position of C.

Using the angle theta and the rotational translation we can then get the actual position of C and similarly D.

Ah, trigonometry!  How wonderful.


Quod erat demonstrandum.  I claim my Librarian's Prize.  First prize is a 5 minute lecture on US politics given by Pete.  The booby prize is his hour-long lecture.

https://en.wikipedia.org/wiki/Solution_of_triangles#Three_sides_given_(spherical_SSS)

[Pete]

...I claim my Librarian's Prize.  First prize is a 5 minute lecture on US politics given by Pete.  The booby prize is his hour-long lecture.
https://en.wikipedia.org/wiki/Solution_of_triangles#Three_sides_given_(spherical_SSS)

Oh come on Richard, the contest isn't over yet, and that prize tiering isn't at all fair! I mean sure, if Bill wins, he'd of course enjoy my 5-min lecture on U.S. politics, but we all know Bill likes boobies a whole lot more! Way to rig the contest!!!

Pete :D

[SierraKen]

OK because well, I put myself up to it, I decided to try this whole thing from scratch myself using the POINT command and my circle code. After working on it all morning and some afternoon, I think I got it down. It might not be perfect, but it's as good as I can do it. What it does is, first it draws the 2 circles you want at the beginning. Then it draws the first one again only this time it looks for the a slight variance in white which is the second one. Using POINT it then displays the coordinate and keeps going until the whole first circle is finished. Because this code breaks the coordinates up more and more the bigger the circle is, I added some IF/THEN statements to make the width lesser the bigger the circles are. Please try it and tell me what you think. I'm sure it's not 100% perfect, but it passes all the tests so far. :) Of course, mathematically 2 circles are only supposed to have up to 2 points of intersection, but like SMcNeill did, it just displays the coords it overlaps. So it can be anywhere from no intersection (which that is detected as well) to maybe 20 intersections or something like that. I draw little red circles around where the intersections are and list the intersection coordinates.

Edit: I added a couple lines of code so the user can't make the 2 circles too similar so only around up to 27 intersections can happen.

Code: QB64: [Select]

1. SCREEN _NEWIMAGE(800, 600, 32)
2. _TITLE "Circles Intersections Detector by Ken G."
3. begin:
4. _LIMIT 100
5. t = 0
6. i = 0
7. CLS
8. begin2:
9. INPUT "X1 (-300 to 300) = ", x1
10. IF x1 < -300 OR x1 > 300 THEN PRINT "Out of range, try again.": GOTO begin2:
11. INPUT "Y1 (-300 to 300) = ", y1
12. IF x1 < -300 OR x1 > 300 THEN PRINT "Out of range, try again.": GOTO begin2:
13. INPUT "Radius 1 (Up to 300) = ", r1
14. IF r1 = 0 OR r1 < 0 THEN PRINT "Radius cannot be 0 or below.": GOTO begin2:
15. INPUT "X2 (-300 to 300) = ", x2
16. IF x2 < -300 OR x2 > 300 THEN PRINT "Out of range, try again.": GOTO begin2:
17. INPUT "Y2 (-300 to 300) = ", y2
18. IF x2 < -300 OR x2 > 300 THEN PRINT "Out of range, try again.": GOTO begin2:
19. INPUT "Radius 2 (Up to 300) = ", r2
20. IF r2 = 0 OR r2 < 0 THEN PRINT "Radius cannot be 0 or below.": GOTO begin2:
21. IF r1 > 300 OR r2 > 300 THEN PRINT "Radius cannot be larger than 300.": GOTO begin2:
22. IF x1 > x2 - 8 AND x1 < x2 + 8 AND y1 > y2 - 8 AND y1 < y2 + 8 AND r1 > r2 - 8 AND r1 < r2 + 8 THEN PRINT "Your 2 circles are too similar, try again.": GOTO begin2:
23.  
24.  
25. x1 = x1 + 400
26. y1 = y1 * -1 + 300
27. x2 = x2 + 400
28. y2 = y2 * -1 + 300
29.  
30. 'Set width between circle points.
31. IF r1 < 100 THEN sc1 = 1
32. IF r1 < 200 AND r1 > 99 THEN sc1 = .5
33. IF r1 < 301 AND r1 > 199 THEN sc1 = .25
34. IF r2 < 100 THEN sc2 = 1
35. IF r2 < 200 AND r2 > 99 THEN sc2 = .5
36. IF r2 < 301 AND r2 > 199 THEN sc2 = .25
37.  
38.  
39. CLS
40. FOR yy = 0 TO 600 STEP 10
41.     LINE (400, yy)-(400, yy + 5), _RGB32(255, 255, 255)
42. NEXT yy
43. FOR xx = 0 TO 800 STEP 10
44.     LINE (xx, 300)-(xx + 5, 300), _RGB32(255, 255, 255)
45. NEXT xx
46.  
47. 'Make first circle and on second GOTO when t=1, check POINT color.
48. one:
49. _LIMIT 100
50. seconds = seconds + sc1
51. s = (60 - seconds) * 6 + 180
52. x = INT(SIN(s / 180 * 3.141592) * r1) + x1
53. y = INT(COS(s / 180 * 3.141592) * r1) + y1
54. FOR z = -1 TO 1
55.     IF t = 1 AND POINT(x + z, y + z) = _RGB32(255, 254, 255) THEN
56.         FOR sz = 1 TO 5 STEP .5
57.             CIRCLE (x, y), sz, _RGB32(255, 0, 0)
58.         NEXT sz
59.         PRINT "("; x - 400; ","; y * -1 + 300; ")"
60.         i = 1
61.     END IF
62. NEXT z
63. FOR sz = .25 TO 2 STEP .25
64.     CIRCLE (x, y), sz, _RGB32(255, 255, 255)
65. NEXT sz
66. IF seconds > 60 THEN
67.     seconds = 0
68.     GOTO two:
69. END IF
70. GOTO one:
71. 'Make second circle.
72. two:
73. _LIMIT 100
74. IF t = 1 THEN GOTO nex:
75. seconds = seconds + sc2
76. s2 = (60 - seconds) * 6 + 180
77. x3 = INT(SIN(s2 / 180 * 3.141592) * r2) + x2
78. y3 = INT(COS(s2 / 180 * 3.141592) * r2) + y2
79. FOR sz = .25 TO 2 STEP .25
80.     CIRCLE (x3, y3), sz, _RGB32(255, 254, 255)
81. NEXT sz
82. IF seconds > 60 THEN
83.     seconds = 0
84.     GOTO three:
85. END IF
86. GOTO two:
87. three:
88. t = t + 1
89. IF t = 2 THEN GOTO nex:
90. GOTO one:
91.  
92. nex:
93. IF i = 0 THEN PRINT "No Intersections": PRINT
94. PRINT "Again (Y/N)?"
95. again:
96. ag$ = INKEY$
97. IF ag$ = "y" OR ag$ = "Y" THEN GOTO begin:
98. IF ag$ = "n" OR ag$ = "N" THEN END
99. GOTO again:
100.  

[STxAxTIC]

Hey all, finally gonna be able to give this more time tonight. Soon...

@bplus
Looks like you nailed this problem long ago, as I expect this should have been well-solved by *someone* out there. Leave it to stackoverflow, right? The equations in your sub versus my code look very similar (as they should in the end) so I was originally curious how you got it, thanks for following up. The solution you point to ends up being a tad simpler than mine, but then again I like to apply methods that don't need to be limited to circles, so the derivation was a longer walk.

@Querkey
Nice picture, nice explanation - I kindof followed it whilst on the road, I'll give it a second look. What I'm seeing is that most everyone solves this the same way, practically with the same picture.

As for the so-called prize, I'm not sure what to do. bplus certainly had the code sooner than anyone, explored this the earliest, in relative terms. Querkey on the other hand, posted his derivation right in the forums... but nobody derived the equations in my code right before my eyes yet...

[bplus]

Ha I think Qwerkey should get some kind of prize should he finish what he started, so far he has gotten us the "a" value and then kind of skipped over the rest. To get a prize one should write some code that produces some results that can be verified or laughed at :D

(Ken and Steve can pick up their prizes right now, Pete your on...)

Can he lay out the a distance along the two origins of circles then draw 2 perpendiculars the length of l in his diagram? (Wait has he solved for l as well as a?) Might not be as simple as Paul Bourke's but pretty easy to follow.

Here is a hint: the angle of the 2nd origin to the first is

angle = _ATAN2(oy2 - oy1, ox2-ox1)

Then point 3 between origins where the perpendiculars hang out the l's is
p3x = ox1 + a*cos(angle)
p3y = oy1 + a*sin(angle)

and so the intersect points should be:
i1x = p3x + l*cos(angle - pi/2)
i1y = p3y + l*sin(angle - pi/2)

i2x = p3x + l*cos(angle + pi/2)
i2y = p3y + l*sin(angle + pi/2)

if there is an intersection which might be good to check before doing all that math for side a ;-))
(PS assuming Qwerkey's side or length "a" is correct.)

Oh wait Qwerkey only has alpha solved, still need x1 (what I was calling side a, from the diagram I am familiar with.)

OK so if his alpha is right we can get x1 and l or y with alpha and R1 as he says algebra on sin and cos definitions.

 cos(alpha) = x1/R1 so x1 = R1*cos(alpha)
 sin(apha) = y/R1 so y or l = R1*sin(alpha)

OK now to check all this with real code.

[bplus]

Hey what do you know!

Code: QB64: [Select]

1. SCREEN 12 'intersect2circs test Qwerkey's and b+ derivation 2019-12-12
2.  
3. C1x = 150
4. C1y = 200
5. C2x = 300
6. C2y = 200
7. r1 = 130
8. r2 = 100
9.  
10. intersect2circs C1x, C1y, r1, C2x, C2y, r2, int1x, int1y, int2x, int2y, p3x, p3y
11. CIRCLE (C1x, C1y), r1, 15
12. CIRCLE (C2x, C2y), r2, 15
13. CIRCLE (int1x, int1y), 2, 14
14. CIRCLE (int2x, int2y), 2, 14
15. CIRCLE (p3x, p3y), 2, 14
16. LINE (int1x, int1y)-(int2x, int2y), 9
17.  
18. SUB intersect2circs (ox1, oy1, r1, ox2, oy2, r2, ix1, iy1, ix2, iy2, p3x, p3y)
19.  
20.     'if there is an intersection which might be good to check before doing all that math for side a ;-))
21.  
22.     d = ((ox1 - ox2) ^ 2 + (oy1 - oy2) ^ 2) ^ .5 'distance between the 2 origins = x1 + x2 in Qwerkey's diagram
23.  
24.     'alpha = _ACOS((r1 ^ 2 + (x1 + x2) ^ 2 - r2 ^ 2) / 2 * r1 * (x1 + x2))
25.     alpha = _ACOS((r1 ^ 2 + d ^ 2 - r2 ^ 2) / (2 * r1 * d)) '<<< need () around whole denominator
26.     '(PS assuming Qwerkey's side or length "a" is correct.)
27.     'It does look right now
28.  
29.     'solve  x1 and l of Qwerkey's diagram
30.     'cos(alpha) = x1/R1 so
31.     x1 = r1 * COS(alpha)
32.     'sin(apha) = y/R1 so y or
33.     l = r1 * SIN(alpha)
34.  
35.     'this is the angle of the 2nd circle origin to the first where alpha resides, on this line is point 3 from which perpendiculars go to intersect points
36.     angle = _ATAN2(oy2 - oy1, ox2 - ox1)
37.  
38.     'Then point 3 between origins where the perpendiculars hang out the l's is
39.     p3x = ox1 + x1 * COS(angle)
40.     p3y = oy1 + x1 * SIN(angle)
41.  
42.     'and so the intersect points should be:
43.     ix1 = p3x + l * COS(angle - _PI / 2)
44.     iy1 = p3y + l * SIN(angle - _PI / 2)
45.  
46.     ix2 = p3x + l * COS(angle + _PI / 2)
47.     iy2 = p3y + l * SIN(angle + _PI / 2)
48.  
49.  
50.  
51. END SUB
52.  
53.  

[STxAxTIC]

Bravo, bplus and Qwerkey.

I absolutely neglected acos and atan in my solution. Yours will inevitably be the best and simplest for two circles. I should make it tougher next time.

How about two ellipses?

[SMcNeill]

Bravo, bplus and Qwerkey.

I absolutely neglected acos and atan in my solution. Yours will inevitably be the best and simplest for two circles. I should make it tougher next time.

How about two ellipses?

My solution will still work; even if you want to detect intersections of polygons, lines, curves, or scribbles.  ;D

[STxAxTIC]

Veeeeery funny Steve! I mean, it's not wrong for shapes that render to the screen. But off-screen? Or even fractals... fascinating... This question has too many fun directions.

[bplus]

Yeah Steve's method will work with anything even solids overlaps, maybe paint an image in the overlap ;)

[STxAxTIC]

Alright, here is what my notes look like as a problem at the end of my Vectors document. This notation is a little dense, but needs no geometry, no pictures. We let vectors do the talking here.