subroutine sphere(xq,yq,zq,a,rho,xp,yp,zp,gx,gy,gz)
c
c  SPHERE calculates the three components of gravitational attraction
c  at a single point due to a uniform sphere. 
c  
c  Input parameters:
c    Observation point is (xp,yp,zp), and center of sphere is at 
c    (xq,yq,zq).  Radius of sphere is a and density is rho.  Density
c    in units of kg/(m**3). All distance parameters in units of km.
c
c  Ouput parameters:
c    Gravitational components (gx,gy,gz) in units of mGal.
c
      real km2m
      data gamma/6.67e-11/,si2mg/1.e5/,pi/3.14159265/,km2m/1.e3/
      ierror=0
      rx=xp-xq
      ry=yp-yq
      rz=zp-zq
      r=sqrt(rx**2+ry**2+rz**2)
      if(r.eq.0.)pause 'SPHERE:  Bad argument detected.'
      r3=r**3
      tmass=4.*pi*rho*(a**3)/3.
      gx=-gamma*tmass*rx/r3
      gy=-gamma*tmass*ry/r3
      gz=-gamma*tmass*rz/r3
      gx=gx*si2mg*km2m
      gy=gy*si2mg*km2m
      gz=gz*si2mg*km2m
      return
      end