line-circle intersection test

//circle-line intersection
//nasarouf@cs.ubc.ca

// from tju1607.cpp http://acm.tju.edu.cn/toj/showp1607.html
// actually does a ray intersection test. but should work for lines too.
// finds the int pt closest to pos
// for ray intersection, if t<0, there is no intersection

#define X first
#define Y second
double hypot(double x, double y){   return sqrt(x*x + y*y); }
#define area2(x1,y1,x2,y2,x3,y3) ((x1)*(y2)+(x2)*(y3)+(x3)*(y1))-((y1)*(x2)+(y2)*(x3)+(y3)*(x1))
typedef pair<double,double> pdd;
pdd unit(pdd d){   double h=hypot(d.X,d.Y);   return make_pair(d.X/h,d.Y/h); }
pdd operator+(pdd a, pdd b){   return make_pair(a.X+b.X,a.Y+b.Y); }
pdd operator*(pdd a, double t){   return make_pair(a.X*t,a.Y*t); }
double dot(pdd a, pdd b){   return a.X*b.X + a.Y*b.Y; }

//dir has to be a unit vector
//returns: pair<yes=1/no=0, t(param, how far is the intersection pt)>
pair<int,double> line_intersects_circle(double x, double y, double r, pdd pos, pdd dir){
  pos.X-=x;   pos.Y-=y;
  x=y=0;
  if (dir.X<0){     dir.X=-dir.X;     pos.X=-pos.X;   }
  if (dir.Y<0){     dir.Y=-dir.Y;     pos.Y=-pos.Y;   }
  pdd pos1=pos+dir;
  double h=area2(pos.X,pos.Y,pos1.X,pos1.Y,x,y); //since base=|dir|=1
  if (r<=abs(h)) return make_pair<int,double>(0,0.0);
  pdd dir1=make_pair(dir.Y,-dir.X);
  pos1=dir1*h;//+dir*(-sqrt(r*r-h*h));
  pos1=pos1+pos*(-1);
  double d,d1=sqrt(r*r-h*h);
  if (abs(dir.X)>abs(dir.Y)) d=pos1.X/dir.X;
  else d=pos1.Y/dir.Y;
  if (d>0) return make_pair<int,double>(1,d-d1);
  else return make_pair<int,double>(1,d+d1);
}