function [xi, yi, zi] = grid2d(X, Y, Z, V, xInc, yInc, zInc)
% GRID2D   - Grids FSI data onto an evenly spaced 3-D mesh
%
% Use as: [XI, YI, ZI, VI] = grid2d(X, Y, Z, V, xInc, yInc, zInc)
%
% Inputs: X,Y,Z,V = Vectors of 3-d data. X,Y, and Z specify
%         the spatial coordinates. V contains the data values.
%         xInc, yInc, zInc = Scalar specifying the distance 
%         between nodes on the nesh grid.
%
% Output: XI,YI,ZI,VI = 3d arrays of gridded data. XI,YI and ZI
%         specify the spatial coordinates while VI contains the data values
%
% Note: VI is determined by binning the data within a each z level while
%       looping through each z-level.

% Brian Schlining
% 25 Jul 2000

fprintf(1, 'Gridding data ...');

% Calculate the 3-d grid points
xMin = min(X);
xMax = max(X);
yMin = min(Y);
yMax = max(Y);
zMin = min(Z);
zMax = max(Z);
[xi, yi, zi] = meshgrid([xMin:xInc:xMax], [yMin:yInc:yMax], [zMin:zInc:zMax]);
vi = ones(size(xi))* NaN;

nMax = length([zMin:zInc:zMax]);
n    = 0;

for zn = min(Z):zInc:max(Z)                   % Loop through all depths by 2m increments
   fprintf(1, '.');
   n = n + 1;
   good = find((Z >= zn) & (Z < (zn + zInc)));  % Find all points in this 2m bin 
   vi(:,:,n) = griddata(X(good), Y(good), V(good), xi(:,:,n), yi(:,:,n));
end

zi = zi + zInc/2;


