CODEWITHSUNDEEP
matlaboptimizationimage-processingnearest-neighbor
akhilc
akhilc

Reputation: 181

How to compare all patches in two images in matlab?

I have two images, let's say A and B, sizes of which need not be same. A be of dimension 255×255 and B be of dimension 100×100. And dimension of a patch is 5×5 for my problem. What I need to do is compare all the patches in A to all the patches in B.

Every patch will be overlapping with some neighboring patches. To clear this point, first patch in A will be A(1:5,1:5) (MATLAB notation). Second patch A(2:6,1:5) and so on, all the way till A(251:255,1:5) at the end of first row, and all the way to A(251:255,251:255) for the last patch in A.

I have to compare each of these patches against all the patches in B. As you can see, there are 251*251 patches in A and 96*96 patches in B. So there are a lot comparisons to be made. And comparison is just Euclidean distance, i.e, I'll just take the difference of two patches and sum the squares. For each patch comparison I'll get a scalar value as result.

I implemented this in MATLAB but it's taking several minutes to execute. So please suggest me the fastest way to implement this. This section is bottlenecking my entire project. The code which I wrote is given below. I'm not an expert clearly so please forgive the mistakes.

row = size(A,1);    
col = size(A,2);

row2 = size(B,1);    
col2 = size(B,2);

patch_long = zeros(5,5,(row2-4)*(col2-4));

idx = 1;

for i = 1:row2-4    
    for j = 1:col2-4

        patch_long(:,:,idx) = B(i:i+4,j:j+4);

        idx = idx+1;

    end

end

%// I rearranged 'B' matrix as 3d matrix with each patch arranged like 
%// slide behind one by one

parfor m = 1:row-4    
    for n = 1:col-4

        temp1 = bsxfun(@minus,(A(m:m+4,n:n+4)),patch_long);
        temp2 = sum(sum(temp1.^2));

        count = sum(temp2 <=threshold);

        if count > 1 
            % Do xyzend
        end

    end
end

%// count counts how many patches in 'B' are close to a particular patch in 'A'.

Upvotes: 3

Views: 2285

Answers (2)

Rody Oldenhuis
Rody Oldenhuis

Reputation: 38032

As mentioned by Mikhail, you actually did a terrific job. There is really not all that much you can improve.

There are only a few minor things you can do:

  1. swap the processing order of the matrices
  2. swap the loop order

Because A is smaller than B, swapping the processing order will cause bsxfun to do more of the heavy lifting. Since that is "closer to the silicon" than a JIT'ed MATLAB loop, this will be beneficial.

Doing a simple tic...toc test with A = rand(255) and B = rand(10) (I don't have all day! :) gives

Elapsed time is 3.460361 seconds.  %// original processing order
Elapsed time is 1.629848 seconds.  %// inverted processing order

For your case, the difference will be less of course, but this could be significant.

Swapping the loop order (loop over rows in the inner loop, columns on the outer loop) works because MATLAB stores data in column-major order. Swapping the order then means the pointer offset calculation can be re-used in the inner loop, thus requiring less operations. Normally, this can make up to ~30% difference, however, recent versions of MATLAB have greatly improved upon this (still not sure how), and you don't use 1 column, but 5...so this will probably not be as significant.

Doing the same test, with the different loop order, gives: 

Elapsed time is 3.462599 seconds. %// original loop order
Elapsed time is 3.272909 seconds. %// inverted loop order

Combining everything: 

Elapsed time is 1.571496 seconds. %// Both optimizations implemented

My code: 

%// random data for testing
A = rand(255);
B = rand(10);
threshold = rand;

[rowA, colA] = size(A);
[rowB, colB] = size(B);

patch_long = zeros(5,5, (rowA-4)*(colA-4));

tic

%// re-arrange A into 3D array of 5×5 patches 
idx = 1;
for jj = 1:colA-4     
    c = jj:jj+4;      %// re-use column indices
    for ii = 1:rowA-4 %// inner loop is over the rows

        patch_long(:,:,idx) = A(ii:ii+4, c);        
        idx = idx+1;

    end

end

%// The heavy lifting
for n = 1:colB-4     %// let bsxfun do more work; loop over smallest array here
    c = n:n+4;       %// re-use column indices
    for m = 1:rowB-4 %// inner loop is over the rows

        temp1 = bsxfun(@minus, B(m:m+4, c), patch_long);
        count = sum( sum(temp1.^2) <= threshold );

    end
end

toc

Note that I didn't use parfor, just to keep everything copy-pastable for everyone. Using parfor on the real problem gives me

Elapsed time is 387.802798 seconds. %// original version, but without parfor
Elapsed time is 362.991760 seconds. %// both optimizations, WITH parfor

but I think this is mostly because my work PC is pretty bad (AMD-A6 3650 APU).

Upvotes: 3

Shai
Shai

Reputation: 114866

You can use im2col to extract all patches

pA = im2col( A, [5 5], 'sliding' );
pB = im2bol( B, [5 5], 'sliding' );
% compute squared difference
D = sum( bsxfun( @minus, permute( pA, [2 3 1] ), permute( pB, [3 2 1] ) ).^2, 3 );

By the way, if you do not need ALL the distances and willing to compromise for an approximation for the k nearest neighbors, you might find PatchMatch very useful:
This is an approximate k-nearest neighbors algorithm tailored for image patches. It is very efficient in terms of memory (space) usage and very fast.

Upvotes: 4

Related Questions