Path: utzoo!utgpu!watmath!att!tut.cis.ohio-state.edu!rutgers!mit-eddie!ll-xn!tj
From: tj@XN.LL.MIT.EDU (Thomas E. Jones)
Newsgroups: comp.graphics
Subject: Re: Tangents to Three Circles
Summary: Anyone else interested?
Message-ID: <1516@xn.LL.MIT.EDU>
Date: 8 Aug 89 21:33:47 GMT
References: <859@mrsvr.UUCP>
Reply-To: tj@xn.UUCP (Thomas E. Jones)
Organization: MIT Lincoln Laboratory, Lexington, MA
Lines: 123

With the help of Macsyma, I worked out the equations a few months ago.

Here's the program.  The few equations should be figurable from this.

It works fairly well.  Any suggestions?

------- CUT HERE --------
#include<stdio.h>
#include<math.h>

/* Find Radius and center of a circle that goes through
 * 3 points (x1,y1) (x2,y2) (x3,y3)
 * Copyright 1988, Thomas E. Jones MIT Lincoln Laboratory tj@xn.ll.mit.edu
 * May be freely distributed and used for any purpose.
 */

double PI;

double circ_x(x1,y1,x2,y2,x3,y3)
double x1,y1,x2,y2,x3,y3;
{ double y32,y22,y12,x12,x32,x22,num,den;
  y32=y3*y3;
  y22=y2*y2;
  y12=y1*y1;
  x12=x1*x1;
  x32=x3*x3;
  x22=x2*x2;
  num= -( (y2-y1)*y32 + (-y22+y12-x22+x12)*y3 + y1*y22 +
		(-y12+x32-x12)*y2+(x22-x32)*y1);
  den= 2*((x2-x1)*y3 + (x1-x3)*y2 + (x3-x2)*y1);
  if (den==0.0){
    return(9999.0);
  }
  else
    return(num/den);
}

double circ_y(x1,y1,x2,y2,x)
double x1,y1,x2,y2,x;
{  double result,num,den;
   if((y2==y1))
     result=9999.0;
   else
     result=((x1-x2)*x/(y2-y1)+(y1+y2)/2+(x2-x1)*(x1+x2)/(2*(y2-y1)));
   return(result);
}

double angle(x,y)
double x,y;
{ double pi;
  pi=3.14159265;
  if (x==0){
    if (y>=0)
      return(pi/2.0);
    else
      return(3.0*pi/2.0);
  }
  else {
    if (x>=0){
      if(y>=0)
        return(atan(y/x));
      else
        return(2*pi+atan(y/x));
    }
    else {
      if(y>=0)
        return(pi+atan(y/x));
      else
 	return(pi+atan(y/x));
    }
  }
}

int circle_conv(x1,y1,x2,y2,x3,y3,x,y,sa,ea,r)
int x1,y1,x2,y2,x3,y3,*x,*y,*r;
double *sa,*ea;
{
  double fx1,fy1,fx2,fy2,fx3,fy3,fx,fy,fr;
  fx1=x1*1.0;
  fy1=y1*1.0;
  fx2=x2*1.0;
  fy2=y2*1.0;
  fx3=x3*1.0;
  fy3=y3*1.0;
  fx=circ_x(fx1,fy1,fx2,fy2,fx3,fy3);
  fy=circ_y(fx1,fy1,fx2,fy2,fx);
  fr=sqrt((fx1- fx)*(fx1- fx)+(fy1- fy)*(fy1- fy));
  *x=fx;
  *y=fy;
  *r=fr;
  *sa=angle(fx3-fx,fy3-fy);
  *ea=angle(fx1-fx,fy1-fy);
  if ((fx==9999.0)||(fy==9999.0))
    return(0);
  else
    return(1);
}

/*
main(){
  int i,x1,y1,x2,y2,x3,y3,x,y,r;
  double sa,ea;
  for(i=1;i<=10;i++){
    printf("enter x and y\n");
    scanf("%lf %lf",&sa,&ea);
    printf("%f %f angle=%f\n",sa,ea,angle(sa,ea)*180.0/3.14159);
  }
  PI=M_PI;
  printf("Enter x1,y1,x2,y2,x3,y3\n");
  scanf("%d,%d,%d,%d,%d,%d",&x1,&y1,&x2,&y2,&x3,&y3);
  circle_conv(x1,y1,x2,y2,x3,y3,&x,&y,&sa,&ea,&r);
  printf("center X is %d\n",x);
  printf("center Y is %d\n",y);
  printf("radius = %d\n",r);
  printf("start angle = %f degrees\n",sa*180.0/PI);
  printf("end angle = %f degrees\n",ea*180.0/PI);
}

------ END OF PROGRAM -----

-- 
tj@xn.ll.mit.edu or tj@ll-xn.arpa          (one of these should work)
Thomas E. Jones, home (617) 279-0767 work (617) 981-5093