findarray(a,b,flag)

function ind = findarray(a,b,flag)
%
% FINDARRAY    Find one array within another
%   I = FINDARRAY(A,B,'first')  for the array B returns an index array of 
%   the same size as B containing the smallest absolute index in A for each
%   element of B which is a member of A and 0 if there is no such index. 
%   FINDARRAY(A,B) is the same as FINDARRAY(A,B,'first').
%   
%   I = FINDARRAY(A,B,'last')  for the array B returns an index array of 
%   the same size as B containing the highest absolute index in A for each
%   element of B which is a member of A and 0 if there is no such index. 
%
%   I = FINDARRAY(A,B,'all') returns a NDIMS(B)+1 dimensional array such 
%   that I(:,...,k) contains the k-th absolute index in A for each element
%   of B, 0 otherwise.
% 
% See also find, ismember
% 
%  Example:
%
%  >> findarray(pascal(3),magic(2),'all')
%
%  ans(:,:,1) =
%       1     6
%       0     5
%  ans(:,:,2) =
%       2     8
%       0     0
%  ans(:,:,3) =
%      3     0
%      0     0
%  ans(:,:,4) =
%      4     0
%      0     0
%  ans(:,:,5) =
%      7     0
%      0     0

% Mukhtar Ullah
% [e-mail]
% First Created: December 12, 2004
% Laste modified: November 16, 2010

if nargin < 3
    flag = 'first';
else
    flag = validatestring(flag, {'first','last','all'});
end

[aa,bb] = ndgrid(a,b);
iseq = aa ~= bb;

if iseq
    ind = zeros(size(b)); 
    return; 
end

ia = (1:numel(a))';
ia = ia(:, ones(1,numel(b)));

switch flag
    case 'first'
        ia(iseq) = inf;
        ind = reshape(min(ia), size(b));
        ind(isinf(ind)) = 0;
    case  'last'
        ia(iseq) = 0;
        ind = reshape(max(ia), size(b));
    case 'all'
        ia(iseq) = inf;
        ia = unique(sort(ia),'rows');
        ia(isinf(ia)) = 0;
        if ~ia(end,:), ia(end,:) = []; end
        [siz{1:ndims(b)}] = size(b);
        ind = reshape(ia.',siz{:},[]);
end