(x21,x22,x23,x41,x42,x43,x81)
; angdif.pam
; Returns:  x81 - angular distance

   x81=COS(x23)*SIN(x22)
   x82=SIN(x23)*SIN(x22)
   x83=COS(x22)

   x84=COS(x43)*SIN(x42)
   x85=SIN(x43)*SIN(x42)
   x86=COS(x42)

   x88=ABS(x81*x84+x82*x85+x83*x86)
   IF(x88.GE.1.0) x88=1.0
   x81=ACO(x88)

   RE