Vectorization of array comparison

Question

I have surjective functions created by matching one element in an array MatchesX.trainIdx to one or more elements in a second array MatchesX.queryIdx.

To obtain only the bijective elements of said funciton I run the same function forward

Matches1=Matcher.match(Descriptors1,Descriptors2);

and then backwards

Matches2=Matcher.match(Descriptors2,Descriptors1);

and then look for the elements occuring in both function in following fashion:

k=1;
DoubleMatches=Matches1;

for i=1:length(Matches1)
    for j=1:length(Matches2)
        if((Matches1(i).queryIdx==Matches2(j).trainIdx)&&(Matches1(i).trainIdx==Matches2(j).queryIdx))
            DoubleMatches(k)=Matches1(i);
            k=k+1;
       end
   end
end

DoubleMatches(k:end)=[];

This of course does the work, but it is rather unelegant and seems to bother the JIT accelerator (calc time with accel on and accel off is the same).

Can you think of a way to vectorize this expresion? Is there any other way of avoiding the JIT from "striking"?

Thanks a lot and sorry about the strange structs, I'm working with MEX-functions. Let me know if rewriting the code in "normal" arrays would help

Answer 1

Access to data in multi-dimensional structures is notoriously slow in MATLAB, so transforming your data to an ordinary array will certainly help:

kk = 1;
DoubleMatches = Matches1;

%// transform to regular array
Matches1queryIdx = [Matches1.queryIdx];
Matches1trainIdx = [Matches1.trainIdx];

Matches2queryIdx = [Matches2.queryIdx];
Matches2trainIdx = [Matches2.trainIdx];

%// loop through transformed data instead of structures
for ii = 1:length(Matches1queryIdx)
    for jj = 1:length(Matches1queryIdx)
        if((Matches1queryIdx(ii)==Matches2trainIdx(jj)) && ...
                (Matches1trainIdx(ii)==Matches2queryIdx(jj)))
            DoubleMatches(kk) = Matches1(ii);
            kk = kk+1;
        end
    end
end

DoubleMatches(kk:end)=[];

There is also a solution that is almost entirely vectorized:

matches = sum(...
    bsxfun(@eq, [Matches1.queryIdx], [Matches2.trainIdx].') & ...
    bsxfun(@eq, [Matches1.trainIdx], [Matches2.queryIdx].'));

contents = arrayfun(@(x)..
    repmat(Matches1(x),1,matches(x)), 1:numel(matches), ...
    'Uniformoutput', false);

DoubleMatches2 = [contents{:}]';

Note that this can be a lot more memory intensive (it has O(N²) peak memory footprint, as opposed to O(N) for the others, although the data type at peak memory is logical and thus 8x smaller than double...). Better do some checks beforehand which one you should use.

A little test. I used the following dummy data:

Matches1 = struct(...
    'queryIdx', num2cell(randi(25,1000,1)),...
    'trainIdx', num2cell(randi(25,1000,1))...
);

Matches2 = struct(...
    'queryIdx', num2cell(randi(25,1000,1)),...
    'trainIdx', num2cell(randi(25,1000,1))...
);

and the following test:

%// Your original method
tic    
    kk = 1;
    DoubleMatches = Matches1;

    for ii = 1:length(Matches1)
        for jj = 1:length(Matches2)
            if((Matches1(ii).queryIdx==Matches2(jj).trainIdx) && ...
                    (Matches1(ii).trainIdx==Matches2(jj).queryIdx))
                DoubleMatches(kk) = Matches1(ii);
                kk = kk+1;
            end
        end
    end

    DoubleMatches(kk:end)=[];

toc

DoubleMatches1 = DoubleMatches;


%// Method with data transformed into regular array
tic

    kk = 1;
    DoubleMatches = Matches1;

    Matches1queryIdx = [Matches1.queryIdx];
    Matches1trainIdx = [Matches1.trainIdx];

    Matches2queryIdx = [Matches2.queryIdx];
    Matches2trainIdx = [Matches2.trainIdx];

    for ii = 1:length(Matches1queryIdx)
        for jj = 1:length(Matches1queryIdx)
            if((Matches1queryIdx(ii)==Matches2trainIdx(jj)) && ...
                    (Matches1trainIdx(ii)==Matches2queryIdx(jj)))
                DoubleMatches(kk) = Matches1(ii);
                kk = kk+1;
            end
        end
    end

    DoubleMatches(kk:end)=[];

toc

DoubleMatches2 = DoubleMatches;


% // Vectorized method
tic

    matches = sum(...
        bsxfun(@eq, [Matches1.queryIdx], [Matches2.trainIdx].') & ...
        bsxfun(@eq, [Matches1.trainIdx], [Matches2.queryIdx].'));

    contents = arrayfun(@(x)repmat(Matches1(x),1,matches(x)), 1:numel(matches), 'Uniformoutput', false);

    DoubleMatches3 = [contents{:}]';

toc

%// Check if all are equal
isequal(DoubleMatches1,DoubleMatches2, DoubleMatches3)

Results:

Elapsed time is 6.350679 seconds. %// (  1×) original method
Elapsed time is 0.636479 seconds. %// (~10×) method with regular array
Elapsed time is 0.165935 seconds. %// (~40×) vectorized
ans =
     1                            %// indeed, outcomes are equal