double AzimuthAngle(double x1,double y1,double x2,double y2)
{
double dx,dy,angle;
dx=x2-x1;
dy=y2-y1;
if(x2 == x1)
{
angle=PI/2.0;
if(y2==y1) angle=0.0;
else if(y2<y1) angle=3.0*PI/2.0;
}
else if((x2>x1) && (y2>y1)) angle=atan(dx/dy);
else if((x2>x1) && (y2<y1)) angle=PI/2+atan(-dy/dx);
else if((x2<x1) && (y2<y1)) angle=PI+atan(dx/dy);
else if((x2<x1) && (y2>y1)) angle=3.0*PI/2.0+atan(dy/-dx);
return (angle*180/PI);
}