The title is probably the three-point coordinates (ax, ay), (bx, by), (cx, cy) of the triangle and the radius r of the circle center (ox, oy), find the intersection area of the triangle and the circle
This question can obviously obtain a certain area, but there are many special situations to consider. It is not as simple as I thought at the beginning. The code was written a long time ago, and it is quite brain-intensive. I will post it as a souvenir. If anyone has a simpler algorithm, please let me know
get_S
function res= get_S(ax,ay,bx,by,cx,cy,ox,oy,r)
ans=0.0;
% Move the center of the circle to the origin
ax=ax-ox;
ay=ay-oy;
bx=bx-ox;
by=by-oy;
cx=cx-ox;
cy=cy-oy;
ox=0;
oy=0;
ans=get_area(ax,ay,bx,by,ox,oy,r) + get_area(bx,by,cx,cy,ox,oy,r) + get_area(cx,cy,ax,ay,ox,oy, r);
res = abs(ans);
end
get_area
function res = get_area(ax,ay,bx,by,ox,oy,r)
sign=get_sign(ax,ay,bx,by);% judge the sign of the directed area
ans=0;
oa=dis(ox,oy,ax,ay),ob=dis(ox,oy,bx,by),ab=dis(ax,ay,bx,by);
l=dis_line(ax,ay,bx,by);
if oa==0 || ob==0 % If one side has a length of 0, then the area is 0
res = 0;
return;
end
% In the first case, both points are inside the circle
if oa<=r & &ob<=r
ans=abs(multi(ax,ay,bx,by))/2.0;% The intersection area is the area of triangle oab
res = sign*ans;
return
% In the second case, both points are outside the circle and the distance from the center of the circle to the straight line ab is greater than the radius
elseif (oa>=r & amp; & amp; ob>=r & amp; & amp; l>=r)
[t1x,t1y]=get_point(ax,ay,r); the intersection point of %oa and the circle
[t2x,t2y]=get_point(bx,by,r); the intersection point of %ob and the circle
d=dis(t1x,t1y,t2x,t2y);
ang=acos(get_angle(r,r,d));% angle t1ot2
ans=abs(ang*r*r/2.0);%The intersection area is the area of fan-shaped ot1t2
res= sign*ans;
return
% In the third case, both points are outside the circle and the distance from the center of the circle to the straight line is less than the radius, but the triangle oab has a base angle that is an obtuse angle, that is, the straight line ab does not intersect the circle
elsif (oa>=r & amp; & amp; ob>=r & amp; & amp;l<=r & amp; & amp;(get_angle(ab,oa,ob)<=0||get_angle(ab, ob,oa)<=0))
[t1x,t1y]=get_point(ax,ay,r); the intersection point of %oa and the circle
[t2x,t2y]=get_point(bx,by,r); the intersection point of %ob and the circle
dist=dis(t1x,t1y,t2x,t2y);%Length of line segment t1t2
ang=acos(get_angle(r,r,dist));% angle t1ot2
ans=abs(ang*r*r/2.0);%The intersection area is the area of fan-shaped ot1t2
res= sign*ans;
return;
%In the fourth case, both points are outside the circle and the distance from the origin to the line ab is less than the radius, and the base angles of the triangle oab are both acute angles, that is, there are two intersection points between the line ab and the circle
elseif(oa>=r & amp; & amp;ob>=r & amp; & amp;l<=r & amp; & amp;get_angle(ab,oa,ob)>0 & amp; & amp;get_angle( ab,ob,oa)>0)
%c, d are the two intersection points of the straight line ab and the circle
if(ax~=bx)% the slope of the straight line ab exists
k=(ay-by)/(ax-bx);% slope of straight line ab
h=ay-k*ax;% line ab intercept
% Solve the quadratic equation in one variable to find the intersection point of the circle and the straight line ab
a0=1.0 + k*k;
b0=2.0*k*h;
c0=h*h-r*r;
cx=(-b0 + sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
cy=k*c.x + h;
dx=(-b0-sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
dy=k*d.x + h;
else% The slope of the straight line ab does not exist
cx=ax;
dx=ax;
cy=sqrt(r*r-ax*ax);
dy=-sqrt(r*r-ax*ax);
end
[t1x,t1y]=get_point(ax,ay,r); the intersection point of %oa and the circle
[t2x,t2y]=get_point(bx,by,r); the intersection point of %ob and the circle
d1=dis(cx,xy,dx,dy);% cd length of line segment
d2=dis(t1x,t1y,t2x,t2y);%Length of line segment t1t2
ang1=acos(get_angle(r,r,d1));% angle cod
ang2=acos(get_angle(r,r,d2));% angle t1ot2
s1=abs(ang1*r*r/2.0);% small sector ocd area
s2=abs(ang2*r*r/2.0);% small and large sector ot1t2 area
s3=abs(multi(cx,xy,dx,dy))/2.0;% triangle ocd area
ans=s2 + s3-s1;%The intersecting area is the area of the triangle ocd plus the difference between the areas of the two sectors
res= sign*ans;
return
% In the fifth case, point a is outside the circle and point b is inside the circle
elseif(oa>=r & amp; & amp; ob<=r)
The two points of %c and d are the two intersection points of the straight line ab and the circle, and the point e is the intersection point between the two points ab
if(ax~=bx)% the slope of the straight line ab exists
k=(ay-by)/(ax-bx);% slope of straight line ab
h=ay-k*ax;% line ab intercept
% Solve the quadratic equation in one variable to find the intersection point of the circle and the straight line ab
a0=1.0 + k*k;
b0=2.0*k*h;
c0=h*h-r*r;
cx=(-b0 + sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
cy=k*cx + h;
dx=(-b0-sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
dy=k*dx + h;
%The intersection point should be the one between the two points c and d whose abscissa coordinates are between the two points (a.x,b.x)
if(ax<=cx & amp; & amp;cx<=bx||ax>=cx & amp; & amp;cx>=bx)
ex=cx;
ey=cy;
else
ex=dx;
ey=dy;
end
else% The slope of the straight line ab does not exist
dx=ax;
cx=dx;
cy=-sqrt(r*r-a.x*a.x);
dy=sqrt(r*r-a.x*a.x);
% The intersection point should be the one between the two points c and d whose ordinate is between the two points (a.y, b.y)
if(a.y<=c.y & amp; & amp;c.y<=b.y||a.y>=c.y & amp; & amp;c.y>=b.y)
ex=cx;
ey=cy;
else
ex=dx;
ey=dy;
end
end
[ t1x,t1y]=get_point(ax,ay,r); the intersection point of %oa and the circle
dist=dis(t1x,t1y,ex,ey);%Length of line segment t1e
ang=acos(get_angle(r,r,dist));% angle eot1
s1=abs(ang*r*r/2.0);% sector oet1 area
s2=abs(multi(bx,by,ex,ey))/2.0;% triangle obe area
ans=s1 + s2;%The intersecting area is the area of the triangle obe plus the area of the sector oet1
res= sign*ans;
return;
% In the fifth case, point b is outside the circle and point a is inside the circle
elseif(ob>=r & &oa<=r)
The two points of %c and d are the two intersection points of the straight line ab and the circle, and the point e is the intersection point between the two points ab
if(ax~=bx)% the slope of the straight line ab exists
k=(ay-by)/(ax-bx);% slope of straight line ab
h=ay-k*ax;% line ab intercept
% Solve the quadratic equation in one variable to find the intersection point of the circle and the straight line ab
a0=1.0 + k*k;
b0=2.0*k*h;
c0=h*h-r*r;
cx=(-b0 + sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
cy=k*cx + h;
dx=(-b0-sqrt(b0*b0-4.0*a0*c0))/(2.0*a0);
dy=k*dx + h;
%The intersection point should be the one between the two points c and d whose abscissa coordinates are between the two points (a.x,b.x)
if(ax<=cx & amp; & amp;cx<=bx||ax>=cx & amp; & amp;cx>=bx)
ex=cx;
ey=cy;
else
ex=dx;
ey=dy;
end
else% The slope of the straight line ab does not exist
dx=ax;
cx=dx;
cy=-sqrt(r*r-a.x*c.x);
dy=sqrt(r*r-a.x*a.x);
% The intersection point should be the one between the two points c and d whose ordinate is between the two points (a.y, b.y)
if(ay<=cy & amp; & amp;cy<=by||ay>=cy & amp; & amp;cy>=by)
ex=cx;
ey=cy;
else
ex=dx;
ey=dy;
end
end
[t1x,t1y]=get_point(bx,by,r); the intersection point of %oa and the circle
dist=dis(t1x,t1y,ex,ey);%Length of line segment t1e
ang=acos(get_angle(r,r,dist));% angle eot1
s1=abs(ang*r*r/2.0);% sector oet1 area
s2=abs(multi(ax,ay,ex,ey))/2.0;% area of triangle oae
ans=s1 + s2;%The intersection area is the area of the triangle oae plus the area of the sector oet1
res= sign*ans;
return
end
Here are the helper functions:
function [x,y] = get_point(ax,ay,r)
if ax~=0 % slope exists
k=ay/ax;
x=abs(r)/sqrt(1.0 + k*k);
if(ax<0)
x=-x;% judge the sign of ans.x
end
y=k*x;
else% slope does not exist
x=0;
if(ay>0)
y=r;% judge the sign of ans.y
else
y=-r;
end
end
end
function res = get_sign(ax,ay,bx,by)
if multi(ax,ay,bx,by)>0
res = 1;
end
res= -1;
end
function res = dis(ax,ay,bx,by) res=sqrt((ax-bx)*(ax-bx) + (ay-by)*(ay-by)); end
function res = get_angle(a,b,c) res=(a*a + b*b-c*c)/(2.0*a*b); end
function res=dis_line(ax,ay,bx,by) res=abs(multi(ax,ay,bx,by))/dis(ax,ay,bx,by);
function [theta1,theta2] = get_theta(x0,y0,r,x1,y1)
l2=(y1-y0)^2 + (x1-x0)^2;
syms x y
[x,y]=solve((x-x0)^2 + (y-y0)^2-r^2,(x-x1)^2 + (y-y1)^2-l2); % get intersection point analytical solution of
v1=[x(1)-x1,y(1)-y1];
v2=[x(2)-x1,y(2)-y1];
theta1=eval(acos(dot(v1,[1,0])/norm(v1)));
if eval(v1(2))<0
theta1=2*pi-theta1;
end
theta2=eval(acos(dot(v2,[1,0])/norm(v2)));
if eval(v2(2))<0
theta2=2*pi-theta2;
end
if theta1>theta2
t=theta1;
theta1=theta2;
theta2 = t;
end