Matrix computation

Question

I have a 384*512*3 RGB matrix. There are only 512 unique colors which repeat themselves with different weights. From them, I have to select half, and the other half must be replaced with the closest elements from the first one.

I thought of looping through the image, and search the closest color for the current one. After finding it, I replace one with the other.

But I have 3 loops 1:384, 1:512, 1:256. With the first two I loop through the RGB matrix, the third is used to loop to the matrix containing the final colors. This takes some time to compute.

What can be done to speed it up?

The loops look like this:

dim=size(RGB);
for i=1:dim(1)
    for j=1:dim(2)
        aux=[RGB(i,j,1) RGB(i,j,2) RGB(i,j,3)];
        minim=RGB_dist(uint8(V_colors_K(1,1:3)),aux);
        index=1;
            for k=1:K
            %index=1;
            if (minim>RGB_dist(uint8(V_colors_K(k,1:3)),aux))
                minim=RGB_dist(uint8(V_colors_K(k,1:3)),aux);
                index=k

            end
            RGB(i,j,1)=V_colors_K(index,1);
            RGB(i,j,2)=V_colors_K(index,2);
            RGB(i,j,3)=V_colors_K(index,3);
        end
    end
end

V_colors_K represent the half colors selected to be the final ones.

There is some minor improvement I could think about. The minimum distance is not needed if the color is in the good half.

Here's the algorithm to be more precise:

Definition 1. Function D(c1,c2) is the distance between two color vectors c1 and c2, such as Euclidean distance.

Definition 2. Function P(c) is the pixel number of color c. Definition 3. Base-color cbase of an initial color set Q is the color that satisfies the equation

Definition 4. A weighted product of color c, V(c), is defined as

where wp is the weight of the number of pixels, and wd is the weight of a color distance.

Given the first color cbase, our method calculates the weighted products of other colors and selects the first K-1 largest product. The corresponding K-1 colors with the base-color are used to form an initial palette. The left N–K colors are merged with the closest colors in the initial palette to produce the final palette.

RGB_dist function:

function[distance]=RGB_dist(x,y)

    distance=sqrt(sum((double(x-y)).^2*[3;4;2],2));

end

I do have a function that works on a whole matrix, it calculates the distance between all pairs.

function[D]=RGB_dist_full_2(x)

I = nchoosek(1:size(x,1),2);

D = squareform(RGB_dist(x(I(:,1),:), x(I(:,2),:)))

end

Then, I would need to get the minimum on each column.

Answer 1

If I'm reading this right, you are applying the operation RGB_dist() pairwise to every color in V_colors_k and every pixel in RBG. If RGB_dist() is a linear function, like a dot product then you can apply it to the whole matrix all at once. For example if it was a dot product then you could replace the entire inner loop with:

DISTS = V_colors_K * RGB(i,j,:)';
k = find( DISTS == min(DISTS(:)) );
RGB(i,j,:) = V_colors_K(k,:);

Without knowing what's in RBG_dist() I can't give you a better answer. The general answer that I can give is: Matlab loops are slow, if you want it to run fast you need to remove all loops and use only matrix operations (which are lightening fast). This process of replacing loops with matrix ops is called vectorizing your code, and it can be tricky. Without knowing what you're doing inside RBG_dist() it's impossible to come up with a full vectorization for you.

My usual workflow in matlab is to write something the intuitive way with loops, like you did, then once it gives me the correct results I go back and figure out how to vectorize it (ie replace all loops with matrix operations) to make it fast. Vectorizing is tricky, it's like a linear-algebra puzzle and the speeding-up always take me WAY longer than writing the code in the first place.

UPDATE:

The best I got uses one loop over your basis colors. You were super close with your RGB_dist function since that line will work perfectly well on whole matrices:

[RGBwidth RGBHeight RGBdepth] = size(RGB);
minDists=inf( [RGBwidth RGBheight] );
bestKs=zeros( [RGBwidth RGBheight] );
for k=1:K
    % make matrix out of the color k, the same shape as RGB
    color_K_mat = premute(repmat(V_colors_K(k,:), [RGBwidth 1 RGBheight]), [3 1 2]);

    % compute the distance from each pixel's color to color k
    dists = sqrt(sum((RGB-color_K_mat).^2, 3));

    % create a binary mask showing which pixels are closer to this color than to any previous one
    mask = (dists < minDists);

    % update your results matrices
    bestKs = not(mask)*bestKs + mask*k
    minDists = min(bestKs, dists);
end

% now bestKs(i,j) gives you the index k of the closest color for pixel (i,j)
% and minDists(i,j) gives you the distance to that color

In theory it should be possible to vectorize out even this loop, but that's a bigger puzzle and I have my own work to do :P

Answer 2

to reduce the inner loop operations:

1. define a vector the size of V_colors_K, with all the values assigned as aux (lets call it v_aux).
2. then calculate RGB_dist on v_aux vs V_colors_K.
3. take the minimum result (it represents the minimum distance color).

Answer 3

Use kmeans:

img = im2double( img );
[IDX,C] = kmeans( reshape( img, [], 3 ), 256 ); %// cluster into 256 clusters
cimg = ind2rgb( reshape( IDX, size(img(:,:,1)) ), C );