I'm making an algorithm which extracts information from images - most of these informations are energy-like data. This is basically done by running through the image with a kernel (size is given as parameter), and get the squared sums of values in that kernel.
This is done in three scales, where the three kernel (patch) sizes are (at this point): smallestPatchSize, smallestPatchSize * 3, smallestPatchSize * 9 (with overlapping kernels in the second and third case). This is done for several color channels, gradient filters, etc (17 of these altogether).
My question is whether it is possible to vectorize the code below; obviously this part of the code requires much more time to run than any other. I'm quite a beginner in Matlab and still try to get the hang of vectorization but this one beats me :)
for dim = 1:17
for i = 1:smallestPatchSize:(size(image,1) - currentPatchSize)
for j = 1:smallestPatchSize:(size(image,2) - currentPatchSize)
% calculate the position in the energy matrix
% which has different dimensions than the input pictures
iPosition = (i - 1 + smallestPatchSize) / smallestPatchSize;
jPosition = (j - 1 + smallestPatchSize) / smallestPatchSize;
% calculate the energy values and save them into the energy matrix
energy(iPosition, jPosition, dim) = sum(sum(abs(...
filters(i:i+currentPatchSize, j:j+currentPatchSize,dim)))) ^ 2;
end
end
end
Thanks in advance - this is my first question @ StackOverflow :)
You're taking local sums-of-values of blocks. A sum-filter is separable, i.e. if you want to do a sum over m-by-n, you can break it down into taking first the m-by-1 sum, and then the 1-by-n sum on the result. Note that with a full convolution you'll take the sum a little bit more often, but it's still faster than using blkproc
or the more recent blockproc
.
Thus, for even smallestPatchSize
, you can write your inner loop as:
tmp = conv2(...
conv2(abs(filter),ones(currentPatchSize,1),'same'),...
ones(1,currentPatchSize,'same').^2;
%# tmp contains the energy for a sliding window filter.
energy = tmp((currentPatchSize/2+1):(currentPatchSize+1):size(image,1)-(currentPatchSize/2+1),...
(currentPatchSize/2+1):(currentPatchSize+1):size(image,2)-(currentPatchSize/2+1));
If you have access to the image processing toolkit you can use blkproc:
B = blkproc(A,[m n],fun)
In your case:
[m,n] = size(filters) ./ patchSize;
fun = @(x) sum( sum(x^2) ); % or whatever
B should end up being the right size. You can also use overlapping regions if you want.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With