Cross correlation of images, adding a loop to go through multiple images

Hello,

I have a code that cross correlates an image 1 and matches it to an image 2. It then outputs the cross-correlation matrix (X).

How can I put a loop around this code for a stack of 3 images (like my_mat_file.mat that I have attached), and save all the cross correlation matrices that result from the loop.

This is my code for cross-correlation of images: (I have also attached my .mat data file to be used with this code)

%%
clc;
clear;
close all;
load('my_mat_file.mat');
%%
template=ims(:,:,1); 
full=ims(:,:,2);
S_full = size(full);
S_temp = size(template);
%%
X=normxcorr2(template, full);
r=max(X(:));
[i,j]=find(X==r);
%%
figure;
subplot(2,2,1), imagesc(full), title('image2');
subplot(2,2,2), imagesc(template), title('image1');
subplot(2,2,3), imagesc(X), rectangle('Position', [j-20 i-20 40 40]);
R = zeros(S_temp);
shift_a = [0 0];
shift_b = [i j] - S_temp;
R((1:S_full(1))+shift_a(1), (1:S_full(2))+shift_a(2)) = full;
R((1:S_temp(1))+shift_b(1), (1:S_temp(2))+shift_b(2)) = template;
subplot(2,2,4), imagesc(R);
%%
figure;
mesh(X);

Answer 1

[Image Analyst]

I don't see a stack of 3 images and a template image in there. ims is a uint16 RGB image. It doesn't make sense to have the template be the full size red channel, and then the other gray scale image be the green channel image. Since your images are so similar, it would just give the trivial result that the closest match is where there is no shift at all. Your template should be smaller.

[Gucci]

yes, the first 2 images are very similar as the particle only moves a little bit, however this is what I am meant to run for upto 10000 images in the ims struct

eg ims(:,:,1) ... all the way to ims(;,:,10000), while incrementing each time,

so for example, how does image 1 relate to image 2,

then how does image 2 relate to image 3, ...... and so on.

How can I just add a for loop or something around this, and save each output of X?

By the end of it I should have 9999 matrices of 'X'... but in my_mat_file, I should be able to save 2 'X' matrices.

[Image Analyst]

OK, so it's just like a movie. You can find the centroid of each image independently. Why are you trying to do a correlation? I don't see any need for that.

[Gucci]

@Image Analyst this is because am comparing traditional centroid tracking with correlation, I have updated the code below for correlation:

I have tried to put a loop, and I am supposed to get 2 values of C and R, as there are 3 images in the video: (The C and R values for image 1 to image 2..... and then The C and R values for image 2 to image 3).

However my loop is wrong, because I am getting 4 values for C and R... How can I fix the loop, and then store the 2 values of C and R in a struct?

%%
clc;
clear;
close all;
load('my_mat_file');

for tempImg=1:2
for fullImg=2:3

template=ims(:,:,tempImg);
full=ims(:,:,fullImg);
S_temp = size(template);
S_full = size(full);

%%
X=normxcorr2(template, full);
r=max(X(:));
[i,j]=find(X==r);

%%
B=im2gray(X);
 temp=graythresh(B);
 bin_im=im2bw(B,temp);

%% centroid of cross-correlation matrix
tot_mass = sum(bin_im(:));
[ii,jj] = ndgrid(1:size(bin_im,1),1:size(bin_im,2));
R = sum(ii(:).*bin_im(:))/tot_mass;
C = sum(jj(:).*bin_im(:))/tot_mass;
  %%
    end
end

[Image Analyst]

You need to index your arrays to save them for separate time points.

Try this:

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;
s = load('my_mat_file.mat')
ims = s.ims;
[rows, columns, numberOfImages] = size(ims);
hFig = figure;
for t = 1 : numberOfImages
    % Get this time point image.
    thisImage = ims(:, :, t);
    subplot(2, 1, 1);
    imshow(thisImage);
    colormap(gray(256));
    impixelinfo;
    caption = sprintf('Time point %d', t);
    title(caption, 'FontSize', fontSize);
    drawnow;
    hFig.WindowState = 'maximized';
    % Find centroid.
    mask = bwareafilt(thisImage > 12000, 1);
    subplot(2, 1, 2);
    imshow(mask);
    props = regionprops(mask, thisImage, 'WeightedCentroid');
    rowCentroid(t) = props.WeightedCentroid(2);
    colCentroid(t) = props.WeightedCentroid(1);
    % Plot it in red crosshairs
    subplot(2, 1, 1);
    hold on;
    plot(colCentroid(t), colCentroid(t), 'r+', 'LineWidth', 2, 'MarkerSize', 100)
    hold off;
end
subplot(2, 1, 2);
plot(colCentroid, rowCentroid, 'r.-', 'LineWidth', 2, 'MarkerSize', 50);
grid on;
xlabel('X Centroid', 'FontSize', fontSize);
ylabel('Y Centroid', 'FontSize', fontSize);

[Gucci]

I understand but the code you have sent me I think is for centroid calculation, however I am supposed to find the cross-correlation matrix between each 2 images and then threshold and find the centroid of the cross-correlation matrix 'X'

From my code below: this means that I must:

1) first compare Image 1 and Image 2, then find the cross correlation matrix which is 'X', (image1 would be template and image 2 would be full)

2) perform otsu thresholding on X and then find the centroid of 'bin_im'

3) repeat steps 1 and 2 for Image 2 and image 3 (so this time, image 2 would be template and image 3 would be full)

so my result will need to have a [C, R] value for the 'X' of Image 1 and Image 2

and then another [C,R] value for the 'X' of Image 2 and Image 3.

I will then need to store each value of [C] in one row of a struct, and then each value of [R] in the next row, underneath the corresponding value of [C]

however in the code below, my for loops are wrong, instead of 2 [C,R] values, I get more... I hope this makes sense

%%
clc;
clear;
close all;
load('my_mat_file');
for tempImg=1:2
for fullImg=2:3
template=ims(:,:,tempImg);
full=ims(:,:,fullImg);
S_temp = size(template);
S_full = size(full);
%%
X=normxcorr2(template, full);
r=max(X(:));
[i,j]=find(X==r);
%% Otsu threshold
B=im2gray(X);
 temp=graythresh(B);
 bin_im=im2bw(B,temp);
%% centroid of cross-correlation matrix
tot_mass = sum(bin_im(:));
[ii,jj] = ndgrid(1:size(bin_im,1),1:size(bin_im,2));
R = sum(ii(:).*bin_im(:))/tot_mass;
C = sum(jj(:).*bin_im(:))/tot_mass;
  %%
    end
end

[Image Analyst]

R(fullImg) = sum(ii(:).*bin_im(:))/tot_mass;
C(fullImg) = sum(jj(:).*bin_im(:))/tot_mass;

[Gucci]

Thanks, this worked to store the values in structs, however there is still a problem with my loop...

The current loop I put will compare image 1 to image 2, and then image 1 to image 3 XXX

(so 'tempImg' is not incrementing for some reason)

instead of what it should do: image 1 to image 2, and then image 2 to image 3 ✓✓✓

How can I fix this? I am not sure what I did wrong in the loop

[Image Analyst]

You don't need two nested loops. You need only 1:

for k = 1 : size(ims, 3)
    template = ims(:, :, k-1); % Template is prior image.
    full = ims(:, :, k); % This is the "current" image.

[Gucci]

Thanks @Image Analyst, this works,

Is there anything I could do to the code to make it run faster (for a set of thousands of images) and still obtain the same results? Maybe I could preallocate the R and C variables, but i don't know how

Here is the updated code,

%%
clc;
clear;
close all;
load('my_mat_file.mat');
%%
for k = 2 : size(ims, 3)
    template = ims(:, :, k-1); % Template is prior image.
    full = ims(:, :, k); % This is the "current" image.
S_temp = size(template);
S_full = size(full);
%%
X=normxcorr2(template, full);
r=max(X(:));
[i,j]=find(X==r);
%%
B=im2gray(X);
 temp=graythresh(B);
 bin_im=im2bw(B,temp);
%% centroid of cross-correlation matrix
tot_mass = sum(bin_im(:));
[ii,jj] = ndgrid(1:size(bin_im,1),1:size(bin_im,2));
R(k) = sum(ii(:).*bin_im(:))/tot_mass;
C(k) = sum(jj(:).*bin_im(:))/tot_mass;

end

[Image Analyst]

You could get rid of this

r=max(X(:));
[i,j]=find(X==r);

since you don't seem to do anything with (the badly-named) i and j.

I still don't know why you're doing this. The cross correlation will just give you how much the weighted centroid is shifted from the prior image. You can get that a lot faster with the regionprops() code I first gave you. Your cross correlation matrix is going to be 4 times as big as your input matrices. Plus it's computationally intensive. It's going to take a long time if you use full size matrices. I've said it before, there's just no need to do it.