Improving performance of mapping algorithm

Almost. I want to know the pairings, which for my data match more or less uniquely. It does not have to be sorted according to distance.

Some data sets have points matching perfectly, others have some distance between matching points. Otherwise I could use the distance matrix and just look for all entries below a tolerance, e.g., 1E-10.

Here is the MATLAB version:

Lk = 4;
Lp = 1500;
pload = rand(Lk,Lp,3);
irs = cell(Lk,1);
ics = cell(Lk,1);
irs{1} = 1:Lp;
ics{1} = 1:Lp;
for k = 2:Lk
    disp('In for loop, k = '); disp(k)
    [irs{k}, ics{k}] = sortpermByPairwiseDistance(pload(1,:,:),pload(k,:,:));
end

function [irf,icf] = sortpermByPairwiseDistance(pl1,pl2) 
  Li = size(pl1,2);
  Lj = size(pl2,2);
  % get distance matrix A[i,j]; i relative points, j points to be mapped
  A = zeros(Li,Lj);
  for kk = 1:3
    A = A + (pl1(:,kk) - pl2(:,kk)').^2;
  end
  [irf, icf] = evalMat(A);
end

function [irf,icf] = evalMat(A) 
    Ar = size(A);
    %% find pairs of smallest distance
    [~,I] = sort(A(:));
    [ir, ic] = ind2sub(Ar,I);
    irf = zeros(Ar(1),1); icf = zeros(Ar(1),1);
    j = 1;
    for i=1:Ar(1)*Ar(2)
        %if ~(ismember(ir(i),irf) || ismember(ic(i),icf))
        if ~(any(irf(:)==ir(i)) || any(icf(:)==ic(i)))
           irf(j) = ir(i);
           icf(j) = ic(i);
           j = j + 1;
        end
    end
end

I assume these iterations with “any(…)” are not efficient with MATLAB, while matrix operations would be better.