Distance between centers

 
/* sphdj - distance to center j */

void sphdj(d1,d2,w,s,j) double *d1, *d2, *w, s[][7]; int j; {
double xx, yy, zz, rr;
xx=w[0]-s[j][0]; xx*=xx;
yy=w[1]-s[j][1]; yy*=yy;
zz=w[2]-s[j][2]; zz*=zz;
rr=s[j][3]; rr*=rr;
*d1=xx+yy+zz-rr;
*d2=xx+zz-rr;
}