Simulate image data representative of a real experiment

As part of my research, I need to validate a particular experimental approach of image capture and particle trajectory analysis by simulating "fake" data where I have control of the parameters.

More specifically, my task is to:

1. generate multiple particles' stochastic trajectories;
2. create a gaussian blur of the simulated point particles, where the gaussian blur is related to the radius parameter I specifiy;
3. discretize that gaussian-blurred particle into a grid corresponding to some pixel resolution related to the real experimental setup; and
4. output each time step in the trajectory as a pixelated image.

I can easily complete Step 1, but my issues reside with Step 2-4. The image I've attached will hopefully provide a useful diagram that may better descripe my goal. The portion in blue is what I desire to output.

The experimental setup involves a microscope viewing top-down the trajectory of particles on a flat plane, so my "fake data" images need to represent that setup.

I've spent quite a while looking through Matlab's imaging capabilitites; however, either due to my ignorance of imaging or to my ignorance of Matlab, I have been unable to come up with an approach that meets my needs.

I appreciate any input anyone can offer.

[J. Alex Lee]

are you able to create a movie/image stack of the unblurred simulated particle trajectories, and all you need is the blur? it's not clear where exactly you are stuck...you can look into imfilter() and conv2() for the blurring.

for these experiments, is it important for you that blur alone won't simulate things like the bright ring (i forget the name of it) when particles are in specific depth relative to focal plane?

[Ryan Muoio]

Yes, exactly. I'm having trouble blurring more than anything else. I want to be able to control the gaussian blur on each individual particle by inputting the radius of the particle. And then I want to create a specific pixelation after the blurring has occurred, which will bin the gaussian into the prescribed pixels.

I'm not sure how imfilter() or conv2() will do what I desire. But, again, I'm not experienced in this functionality of Matlab. I aim to perform this simulation to help our group with something specific, but I'm not at all specialized in this area.

In some sense, we want to create the "perfect" data, limited only by pixel resolution. So, we don't mind that other experimental factors, uncertainties, and imperfections are left unconsidered.

[Ryan Muoio]

Wow! Thanks for your feedback and examples. I will modify the code to match my needs, but I think you've given me a great foundation.

Since you techincally "commented" rather "answered" my question, I don't know how to accept your answer.

[Bjorn Gustavsson]

OK, since it was useful I'll move it to "answer"...

You can look for the keyword "splatting methods" in papers about calculating images of 3-D projections for tomographic imaging. We used a variant of that looks pretty much like this (with additional geometry-handling and such).

[Ryan Muoio]

Thanks! I accepted your answer.

I'll also look into "splatting methods." I've never heard of those before. Sounds interesting.

Answer 1

[Bjorn Gustavsson]

If it is enough for you to put particles at discrete pixel positions you could do something like this:

nP = 123;
dIm  = spalloc(1024,1024,nP); % Or bigger if you want finer "accuracy"
x = randi(1024,[nP,1]);
y = randi(1024,[nP,1]);
dx = 7; % Just some arbitrary Gaussian widths
dy = 6; % For you to adjust
[X,Y] = meshgrid(-15:15);
fK = exp(-X.^2/dx^2-Y.^2/dy^2);
fK = fK/sum(fK(:));
Im = conv(full(dIm),fK,'same');

If you want different particles blurred with Gaussians with different widths you might get a good enough result if you separate your particles into different size-groups:

nP = [12,23,34,45,56]; % Number of particles in each group
wG = [ 2, 4, 6, 8,12]; % Gaussian 1/e half-widths in pixels
sz = [1024,1024];
Im = zeros(sz);
dIm = Im;
[X,Y] = meshgrid(-25:25);
for i1 = 1:numel(nP)
  x = randi(1024,[nP(i1),1]);
  y = randi(1024,[nP(i1),1]);
  idx = sub2ind(sz,y,x);
  dIm(idx) = 1;
  fK = exp(-(X.^2+Y.^2)/wG(i1)^2);
  fK = fK/sum(fK(:));
  Im = Im + conv2(dIm,fK,'same');
  dIm(idx) = 0;
end

This way you get some separation into different blurrings of your particles, not prefectly continuos range of sizes, but you can at least start to get some variability in, likewise you will not get particles distributed uniformly over the image field - only uniformly over the pixel centres. But this should be a simple first step. If you need you can refine this by playing tricks with distributing particles between four neighbouring pixels - or do it with more analytical aproaches (but the utility of that comes down to how well you will know the image blurring in practice...)

HTH

Answer 2

[J. Alex Lee]

Bjorn's answer contains the conv2() route to blurring

fK = exp(-X.^2/dx^2-Y.^2/dy^2);
fK = fK/sum(fK(:));
Im = conv(full(dIm),fK,'same');

And if you have image processing, you could do it somewhat simpler as

Im = imgaussfilt(dIm,[dy,dx]); % or dx,dy, depending on which way is x

As for "pixelating", if you just mean binning the gray values, you could just use rounding like

s = 8; % round to nearest 6th gray value
ImP = round(Im/8)*8

[Ryan Muoio]

Thanks, J. Alex Lee! I will apply Bjorn Gustavsson's foundational code and add your pixelation idea. Hopefully, everything will work.