how to track beam light motion in X-ray microscope?

In X-ray microscopy we acquire several frames at each X-ray energy point and merge them together with an aim to improving signal to noise ratio of the final image. The evil comes from the fact that beam is unstable and in the state of flux. Thus, due to the beam motion there is a limit to the number of frames we can acquire and merge to increase SNR. Is there any way to track beam motion using MATLAB and eventually take the motion into account during data acquisition? Any help is much appreciated. Thank you in advance.

Two flat field images taken at different period of time are here: http://imgur.com/UAP768E

I am not sure that beam motion in the two flat field images is noticeable visually though, because the beam moves around the image center following an unknown trajectory with a shift of a few pixels. Also, I am not sure it is possible to outline beam boundaries precisely because of light diffusion.

[Walter Roberson]

I wonder if we could get a few images so we could see the effect of the beam motion ?

[Alex]

I would be grateful for your advice, I posted a link to two flat field images. Thank you.

[Alex]

Please, let me know if you would need more details (experimental setup, images, etc.)

Answer 1

[Image Analyst]

What do you get if you simply find the weighted centroid of each image?

[Alex]

We are working with 16 bit images. I took two flat fields (min intensity ca. 8000, max ca. 30000 ) with the acquisition time of 0.6 s and deadtime of around 1 hr and did the following:

Centroid

   bw = flat > 10000;  % thresholding
  imagesc(flat); colormap gray
  axis image
  L = bwlabel(bw,8);
  s = regionprops(L, 'Centroid');
  hold on
  for k = 1:numel(s)
      plot(s(k).Centroid(1), s(k).Centroid(2), 'r*')
  end
  hold off

Weighted Centroid

bw = flat > 10000;  % thresholding
imagesc(flat); colormap gray
axis image
s = regionprops(L, 'PixelIdxList', 'PixelList');
idx = s(1).PixelIdxList;
sum_region1 = sum(flat(idx));
x = s(1).PixelList(:, 1);
y = s(1).PixelList(:, 2);

xbar = sum(x .* double(flat(idx))) / sum_region1;
ybar = sum(y .* double(flat(idx))) / sum_region1;

hold on
for k = 1:numel(s)
    idx = s(k).PixelIdxList;
    pixel_values = double(flat(idx));
    sum_pixel_values = sum(pixel_values);
    x = s(k).PixelList(:, 1);
    y = s(k).PixelList(:, 2);
    xbar = sum(x .* pixel_values) / sum_pixel_values;
    ybar = sum(y .* pixel_values) / sum_pixel_values;

      plot(xbar, ybar, '*')
  end
  hold off

Here you can see both results: http://imgur.com/eGOlCPs

I am not sure why there are so may connected objects, I expected only one - the beam centred in the centre of the image...Nethertheless, we can see that * in the center of the image shifts as we move from flat 1 to flat 2 . Do you think it is a reliable way to track beam motion?

[Image Analyst]

You can call bwareaopen() to get rid of all the small noise specks at the top and bottom of the thresholded area. Then you can do

s = regionprops(L, flat, 'WeightedCentroid');

to get the weighted centroid directly without calculating it manually.

[Alex]

Thank you for your time. It looks cool now http://imgur.com/ufE3anh Is there any other method to test out for comparison purpose?

[Image Analyst]

Why do you have no decimal places of accuracy? You should get subpixel resolution on the centroid. For another method, you could try to fit the image to a model, like a 2D quadratic, and look for the peak. Try this: http://www.mathworks.com/matlabcentral/fileexchange/34765-polyfitn

[Alex]

Thank you, it was a bad idea to use round to coordinates. From your expertise, do you think it is possible eventually to describe the motion using some function (polynomial, etc.)? Also, do you think it is possible to predict the motion? Thank you very much!

[Image Analyst]

It may or may not be possible to describe the motion, with a polynomial or anything else. I don't know what causes the shift in the beam. Is it a periodic/cyclic motion, or is it chaotic? Try plotting the centroid and see what its path looks like.