function x = sph_gradx ( z, rx, ry ) ;
% function x = sph_gradx ( z, rx, ry ) ;
%
% rx and ry must be specified in radian;
%
% SEE : sph_grad, sph_div, sph_curl

  global RADUS RADIAN DEGREE

%---- check and reshape input if necessary
  [zy,zx] = size(z); nx = length(rx); ny = length(ry);
  if ( zx ~= nx );
    z = reshape ( z, nx, ny )';
  end;

%---- calculations                   
  ix2=[2:nx 1]; ix1=[nx 1:(nx-1)];

  dx	= rx(ix2) - rx(ix1);
  dx(1) = dx(1) + 2*pi; dx(nx) = dx(nx) + 2*pi;

  acosdx = RADUS * cos( ry' ) * dx ;

  x = ( z(:,ix2) - z(:,ix1) ) ./ acosdx ;
  % check end points to ensure polar values are reasonable
  if cos(ry( 1))<1e-5; x( 1,:) = mean(x(   2,:)) * ones(1,nx); end
  if cos(ry(ny))<1e-5; x(ny,:) = mean(x(ny-1,:)) * ones(1,nx); end

%---- reshape output if necessary
  if ( zx ~= nx );
    x = reshape ( x', zy, zx );
  end;
