% hist3D
% usage: H = hist3D(X,Y,Z,m,n,r);
% find a three dimensional correlation histrogram of three matrix quantities 
% the returned nxmxr matrix array H satisfies 
%
%         H(i,j,k) = # of X,Y and Zvalues satifying  
%         x0 +(i-1)*dx <= X < x0 + i*dx 
%         y0 +(j-1)*dy <= Y < y0 + j*dy
%         z0 +(k-1)*dy <= Z < z0 + k*dy
%
% where
%          x0 = minX, dx = m/(maxX-minX) 
%          y0 = minY, dy = n/(maxY-minY)
%          z0 = minY, dz = r/(maxZ-minZ)

function  H = hist2D(X,Y,Z,m,n,r);

X = X(:);
Y = Y(:);
Z = Z(:);
mX = min(X);
MX = max(X);
mY = min(Y);
MY = max(Y);
mZ = min(Z);
MZ = max(Z);

H = zeros(m,n,r);
X = ceil((m/(MX-mX))*(X-mX));
J = find(X<1);
X(J) = 1;
J = find(X>m);
X(J) = m;

Y = ceil((n/(MY-mY))*(Y-mY));
J = find(Y<1);
Y(J) = 1;
J = find(Y>n);
Y(J) = n;

Z = ceil((n/(MZ-mZ))*(Z-mZ));
J = find(Z<1);
Z(J) = 1;
J = find(Z>n);
Z(J) = n;


for i=1:m
   for j=1:n
      for k=1:r
         H(i,j,k)= sum((X ==i)&(Y==j)&(Z==k));
         end;
      end;
   end;