isec_circlecircle

Help Contents geomland Intersection isec_circlecircle
isec_lineline isec_linecircle

geomland.isec_circlecircle :a :b


Circle-Circle intersection

Definition:

to geomland.isec_circlecircle :a :b
	;01        (a.x-X)^2 + (a.y-Y)^2 = a.r^2
	;02        (b.x-X)^2 + (b.y-Y)^2 = b.r^2
	;03=01-02  ((a.x-X)-(b.x-X))*((a.x-X)+(b.x-X)) + ((a.y-Y)-(b.y-Y))*((a.y-Y)+(b.y-Y)) = a.r^2-b.r^2
	;04        (a.x-X-b.x+X)*(a.x-X+b.x-X) + (a.y-Y-b.y+Y)*(a.y-Y+b.y-Y) = a.r^2-b.r^2
	;05        (a.x-b.x)*(a.x+b.x-2*X) + (a.y-b.y)*(a.y+b.y-2*Y) = a.r^2-b.r^2
	;06        (a.x-b.x)*(a.x+b.x)-2*(a.x-b.x)*X + (a.y-b.y)*(a.y+b.y)-2*(a.y-b.y)*Y = a.r^2-b.r^2
	;06        -2*(a.x-b.x)*X -2*(a.y-b.y)*Y + a.x^2-b.x^2 + a.y^2-b.y^2 = a.r^2-b.r^2
	;08        A=-2*(a.x-b.x)  B=-2*(a.y-b.y)  C=+a.r^2 -b.r^2 -a.x^2 +b.x^2 -a.y^2 +b.y^2
	make local "A -2 * (:a.center.x-:b.center.x)
	make local "B -2 * (:a.center.y-:b.center.y)
	make local "C :a.radius*:a.radius - :b.radius*:b.radius - :a.center.x*:a.center.x + :b.center.x*:b.center.x - :a.center.y*:a.center.y + :b.center.y*:b.center.y

	if :A<>0
	[
		;09        A*X+B*Y=C
		;10        X=(C-B*Y)/A= C/A-(B/A)*Y   (if A not 0)
		;11(10,01) (a.x-C/A+(B/A)*Y)^2 + (a.y-Y)^2 = a.r^2
		;12        (a.x-C/A)^2 +2(a.x-C/A)*(B/A)*Y + ((B/A)*Y)^2 + (a.y^2-2*a.y*Y+Y^2 = a.r^2
		;13        ((B/A)^2+1)*Y^2+(2(a.x-C/A)*(B/A)-2*a.y)*Y+((a.x-C/A)^2+a.y^2-a.r^2)=0
		;14        P=((B/A)^2+1)    Q=((a.x-C/A)*(B/A)-a.y)   R=((a.x-C/A)^2+a.y^2-a.r^2)
		;          (multiply back with A^2
		;14a       P=B^2+A^2    Q=(a.x*A-C)*B-a.y*A^2   R=(a.x*A-C)^2+(a.y^2-a.r^2)*A^2
		make local "P :B*:B+:A*:A
		make local "Q (:A*:a.center.x-:C) * :B-:a.center.y*:A*:A
		make local "R (:A*:a.center.x-:C) * (:A*:a.center.x-:C) +:A*:A * (:a.center.y-:a.radius) * (:a.center.y+:a.radius)
		;15        P*Y^2 + 2*Q*Y +R = 0
		;16        D=sqrt(Q^2-P*R)  (if Q^2>=P*R)
		make local "D :Q*:Q-:P*:R
		; D<0 => no intersection
		if :D<0 [output set]

		; D=0 => 1 point
		if :D=0
		[ 
			make local "Y1 (-:Q/:P)
			make local "X1 (:C-:B*:Y1)/:A
			output (set point :X1 :Y1)
		]

		;16a       D=sqrt(D)
		make "D sqrt :D

		;17        Y1 = (-Q+D)/P    X1 = (C-B*Y1)/A
		;18        Y2 = (-Q-D)/P    X2 = (C-B*Y2)/A
		make local "Y1 (-:Q+:D)/:P
		make local "X1 (:C-:B*:Y1)/:A
		make local "Y2 (-:Q-:D)/:P
		make local "X2 (:C-:B*:Y2)/:A

		output (set point :X1 :Y1 point :X2 :Y2)
	]
	[
		;09        A*X+B*Y=C   (A=0)
		;10        B*Y=C
		if :B=0 
		[ 
			if :C=0
			[ output :a ]
			[ output set ]
		]

		make local "Y :C/:B
		;11        (a.x-X)^2 = a.r^2-(a.y-Y)^2
		;12        (a.x-X)^2 = (a.r-a.y+Y)*(a.r+a.y-Y)
		;13        (a.x-X)   = +- sqrt[(a.r-a.y+Y)*(a.r+a.y-Y)]
		;14        X = a.x -+ sqrt[(a.r-a.y+Y)*(a.r+a.y-Y)]
		make local "X1 :a.center.x-(sqrt (:a.radius-:a.center.y+:Y)*(:a.radius+:a.center.y-:Y)) 
		make local "X2 2*:a.center.x-:X1
		output (set point :X1 :Y point :X2 :Y)
	]  
end

isec_lineline isec_linecircle