function tracer_command ( command )
disp ( command );

global to;


current_point = get ( to.axis, 'CurrentPoint' );
step_points = [current_point(1,1) current_point(1,2)]
x0 = current_point(1,1); y0 = current_point(1,2);
hold on
plot ( x0, y0, 'r*' );


depth = str2num(get(to.depth_edit,'string')) * -1;


h1 = mcvgt ( to.cdf, 'h1', [0 0], [-1 -1] );
h1 = h1';
h2 = mcvgt ( to.cdf, 'h2', [0 0], [-1 -1] );
h2 = h2';

ang = mcvgt ( to.cdf, 'ang', [0 0], [-1 -1] );
ang = ang';

velocities = [];

for time_step = 1:(length(to.jd)-1)

	% 
	% flush the event queue so we can see what is going on.
	drawnow;

	set ( to.title_label, 'string', jd2str(to.jd(time_step) ) );
	u = zslice(to.cdf, 'u', time_step, depth );
	v = zslice(to.cdf, 'v', time_step, depth );

	%
	% Find the nearest xy point
	near_ind = nearxy(to.x(:),to.y(:),x0,y0);

	if ( isnan(to.depth(near_ind)) )
		;
	else
	
		%
		% compute the ij indices
		[r,c] = size(to.x);
		i = mod(near_ind-1,r)+1;
		j = floor((near_ind-1)/r)+1;
		
	
		%
		% calculate the corner points of the grid cell.  Need
		% to use the form factors h1 and h2 to do this.
	    corner_pts = [ (to.x(near_ind) - h1(near_ind)/1000/2) ...
					   (to.x(near_ind) - h1(near_ind)/1000/2) ...
					   (to.x(near_ind) + h1(near_ind)/1000/2) ...
					   (to.x(near_ind) + h1(near_ind)/1000/2); ... 
		               (to.y(near_ind) - h2(near_ind)/1000/2) ...
					   (to.y(near_ind) + h2(near_ind)/1000/2) ...
					   (to.y(near_ind) - h2(near_ind)/1000/2) ...
					   (to.y(near_ind) + h2(near_ind)/1000/2) ];
	
		rot_matr = [cos(ang(near_ind)*-1) sin(ang(near_ind)*-1); ...
				   -sin(ang(near_ind)*-1) cos(ang(near_ind)*-1)];
		corner_pts(1,:) = corner_pts(1,:) - to.x(near_ind);
		corner_pts(2,:) = corner_pts(2,:) - to.y(near_ind);
		corner_pts_rotated = rot_matr * corner_pts;
		corner_pts_rotated(1,:) = corner_pts_rotated(1,:) + to.x(near_ind);
		corner_pts_rotated(2,:) = corner_pts_rotated(2,:) + to.y(near_ind);
			
		x_corner = corner_pts_rotated(1,:);
		y_corner = corner_pts_rotated(2,:);
	
		%
		% Find the nearest corner neighbor.  This will tell us
		% which indices to use for the uv calculation.
		% The ordering is SW, NW, SE, NE.
		near_corner_ind = nearxy ( x_corner, y_corner, x0, y0 );

		max_ind = length(to.x);

		switch ( near_corner_ind )
	
			% SW
			case 1
				u_ind = [ ( near_ind - r ) ...
						  ( near_ind ) ...
						  ( near_ind - r + 1 ) ...
						  ( near_ind + 1 ) ];
				ifind = find((u_ind>=1) | u_ind<=max_ind);
				u_ind = u_ind(ifind);
				u_cell = u(u_ind);

				v_ind = [ ( near_ind - 1 ) ...
						  ( near_ind + r - 1 ) ...
						  ( near_ind ) ...
						  ( near_ind + r ) ];
				ifind = find((v_ind>=1) | v_ind<=max_ind);
				v_ind = v_ind(ifind);
				v_cell = v(v_ind);
	
				
			% NW
			case 2
				u_ind = [ ( near_ind ) ...
						  ( near_ind + r ) ...
						  ( near_ind + 1 ) ...
						  ( near_ind + r + 1 ) ];
				ifind = find((u_ind>=1) | u_ind<=max_ind);
				u_ind = u_ind(ifind);
				u_cell = u(u_ind);

				v_ind = [ ( near_ind - 1 ) ...
						  ( near_ind + r - 1 ) ...
						  ( near_ind ) ...
						  ( near_ind + r ) ];
				ifind = find((v_ind>=1) | v_ind<=max_ind);
				v_ind = v_ind(ifind);
				v_cell = v(v_ind);
	
	
				
			% SE
			case 3
				u_ind = [ ( near_ind - r ) ...
						  ( near_ind ) ...
						  ( near_ind - r + 1 ) ...
						  ( near_ind + 1 ) ];
				ifind = find((u_ind>=1) | u_ind<=max_ind);
				u_ind = u_ind(ifind);
				u_cell = u(u_ind);

				v_ind = [ ( near_ind ) ...
						  ( near_ind + r  ) ...
						  ( near_ind + 1 ) ...
						  ( near_ind + r + 1 ) ];
				ifind = find((v_ind>=1) | v_ind<=max_ind);
				v_ind = v_ind(ifind);
				v_cell = v(v_ind);
	
	
	
			% NE
			case 4
				u_ind = [ ( near_ind ) ...
						  ( near_ind + r ) ...
						  ( near_ind + 1 ) ...
						  ( near_ind + r + 1 ) ];
				ifind = find((u_ind>=1) | u_ind<=max_ind);
				u_ind = u_ind(ifind);
				u_cell = u(u_ind);
	
				v_ind = [ ( near_ind ) ...
						  ( near_ind + r  ) ...
						  ( near_ind + 1 ) ...
						  ( near_ind + r + 1 ) ];
				ifind = find((v_ind>=1) | v_ind<=max_ind);
				v_ind = v_ind(ifind);
				v_cell = v(v_ind);

		end
	
	
		%
		% Now compute the velocity at that point.
		% For now, just let this be the average.
		% And convert it to km/day
		u_point = mean(u_cell)/1000*86400;
		v_point = mean(v_cell)/1000*86400;
	
		psliceuv ( x0, y0, [u_point+sqrt(-1)*v_point],1,1,'black');
	
	
		%
		% Compute the next point.
		next_x = x0 + u_point * (to.jd(time_step+1) - to.jd(time_step)); 
		next_y = y0 + v_point * (to.jd(time_step+1) - to.jd(time_step)); 
	
		plot ( next_x, next_y, 'k*' );
	
		step_points = [step_points; [next_x next_y]];
		x0 = next_x;
		y0 = next_y;

	end
	
end

plot ( step_points(:,1), step_points(:,2), 'k*' );

