Thilina's Blog

I might be wrong, but…

Process spatial operations without ‘for’ loops in MATLAB

While working on some image enhancement applications which the intensity of a given pixel is also a function of its location (spatial coordinates), I used to use ‘for’ loops to iterate each pixel at once and calculate its intensity based on (x,y) value obtained from those two ‘for’ loops. This is very annoying since it takes quite some time to process two ‘for’ loops as well as it is a belief that the MATLAB ‘for’ loops are used to slow. In this article I share my experience how I avoided the requirement of two ‘for’ loops by using the properties of complex numbers for a photo effect which I described in a previous article of mine.

Normally in an image there is an mxnx3 matrix which represents red, green and blue intensities respectively. In my case I append another 4th matrix of coordinates which makes the image matrix as mxnx4. (Even-though there is a memory drawback in this design, it is very easy to handle rather than two nested ‘for’ loops).

This 4th matrix contains complex numbers which their real value is equal to row index and its imaginary value is equal to its column index. And now we do have another advantage. Using the properties of complex numbers we can easily calculate the distance and angle between two pixels in interest which is very useful in image effects such as ‘focal BW’.

First of all let’s look how to create this coordinate matrix. Say the matrix which contains your image data is ‘img’ with the datatype of ‘uint8’.

First cast matrix into double.

img = cast(img, ‘double’);

Next step is to get the dimensions of the image;

[r c d] = size(img);

Then create a zero matrix of size rxc.

cord = zeros(r,c,1);

Now fill the coordinate matrix using repmat command as below.

u = 1:r;
v = 1:c;
p = repmat(v,[r 1]);
q = repmat(u’, [1 c]);
cord = q + (1i*p);

Now your coordinate matrix is ready. Append this to your image matrix;

img(:,:,4) = cord;

Now wipe unwanted ‘cord’ matrix from memory.

clear cord;

Now you can perform any coordinate based example as in matrix handling without using for loops. Following example demonstrate how to use this method for focal BW effect. This is the traditional code I used in my previous article.

clear all;
close all;
clc;
%% Parameters
scale = 0.5;
grade = 5;
radius = 1.7;
%% File select
[FileName,PathName] = uigetfile({'*.jpg';'*.jpeg';'*.JPG'},'SelectImage');
file = [PathName,FileName];
img = cast(imresize(imread(file),scale),'double');
[r c d] = size(img);
IMG = zeros(r,c,d);
%% Pick point in interest
h = figure; imshow(cast(img,'uint8'));
[X Y] = getpts(h);
c0 = floor(X(1));
r0 = floor(Y(1));
c_ = max(c0,(c-c0));
r_ = max(r0,(r-r0));
L = sqrt(r_^2 + c_^2);
close(h);
%% Operation
frame = img;
bwframe = cast(rgb2gray(cast(frame,'uint8')),'double');
for i = 1:r
    for j = 1:c
        Ln = sqrt((r0-i)^2 + (c0-j)^2);
        fn = (L-Ln)./L;
        f = (radius*fn^grade);
        if f>1
            f = 1;
        end
        IMG(i,j,1) = (f*frame(i,j,1))+((1-f)*bwframe(i,j));
        IMG(i,j,2) = (f*frame(i,j,2))+((1-f)*bwframe(i,j));
        IMG(i,j,3) = (f*frame(i,j,3))+((1-f)*bwframe(i,j));
    end
end
%% Display
imG = cast(IMG,'uint8');
imshow(imG);

And the modified code will be now as follows.

clear all;
close all;
clc;
%% Parameters
scale = 0.5;
grade = 5;
radius = 1.7;
%% File select
[FileName,PathName] = uigetfile({'*.jpg';'*.jpeg';'*.JPG'},'SelectImage');
file = [PathName,FileName];
img = cast(imresize(imread(file),scale),'double');
[r c d] = size(img);
IMG = zeros(r,c,d);
cord = zeros(r,c,1);
u = 1:r;
v = 1:c;
p = repmat(v,[r 1]);
q = repmat(u', [1 c]);
cord = q + (1i*p);
img(:,:,4) = cord;
clear cord;
%% Pick point in interest
h = figure; imshow(cast(img(:,:,1:3),'uint8'));
[X Y] = getpts(h);
c0 = floor(X(1));
r0 = floor(Y(1));
c_ = max(c0,(c-c0));
r_ = max(r0,(r-r0));
L = sqrt(r_^2 + c_^2);
close(h);
%% Operation
frame = img(:,:,1:3);
bwframe = cast(rgb2gray(cast(frame,'uint8')),'double');
fn = (L - abs(img(r0,c0,4)-img(:,:,4)))/L;
f = ((radius*fn).^grade);
f(f>1) = 1;
IMG(:,:,1) = (f.*frame(:,:,1))+((1-f).*bwframe);
IMG(:,:,2) = (f.*frame(:,:,2))+((1-f).*bwframe);
IMG(:,:,3) = (f.*frame(:,:,3))+((1-f).*bwframe);
%% Display
imG = cast(IMG,'uint8');
imshow(imG);

Results are same as previous code are below.

bwFocus_000

Result

Hope you got an idea of how to avoid traditional for loop in spatial coordinates based operation. In this case you need to choose your option between memory and speed. Thank you very much for reading.

2013 October 18 - Posted by | Image Processing, MATLAB | , , , ,

1 Comment »

  1. Good job but it would have been nice if you would have mentioned how much was the speed-up.

    Comment by Hell | 2014 November 13 | Reply


Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: