In reply to Seth Ness.

IF p=(x1, y1) and q=(x2, y2) then the dot product of those two vectors is computed as
p.q = |p| |q| cos(theta), where theta is the angle between them.
Therefore
pq = x1*x2 + y1*y2;
pNorm = sqrt(x1*x1 + y1*y1);
qNorm = sqrt(x2*x2 + y2*y2);
theta = arcos( pq / (pNorm*qNorm) );

In the SAS/IML language, the computation is easier since you can use vectors:

proc iml;
p = {2, 1};
q = {1, 3};
theta = arcos( p`*q / (norm(p)*norm(q)) ); /* angle in radians */
degree = theta*180/constant('pi'); /* angle in degrees */
print theta degree;