Calculating number of particles in each 64 x 64 blocks

Hi,

I have an image, consisting of particles of different shapes over a black background. I would like to:

1) Count the number of particles, the area and diameter of each particle in each 64 x 64 pixels blocks across the image; 2) Find the SNR in each 64 x 64 pixels blocks.

I would be happy if you could direct me.

Answer 1

[John BG]

I am answering your question in 2 parts:

Since both parts require a bit of insight I kindly ask you to thumbs-up vote each part separately.

In this first answer let's find a way to count the particles, in a second answer we will find a way to analyze SNR, CNR, CIR and other parameters that ARE POSSIBLE to calculate, ignore the comment above telling you SNR is not possible here, it is, but you have to be aware of clutter and what happens when particles too close each other, too close to background noise, or simply, the initial processing splits a particle into apparent 2 or more particles, and then you have to merge them back.

Let's start with the head counting:

if you have the Image Toolbox, just check the starter tutorial that shows how to enhance contrast (imadjust grarythresh bwareaopen) and count (bwconncomp, field: NumObjects) rice grains

http://uk.mathworks.com/help/images/image-enhancement-and-analysis.html

following, manual mode, i mean, without the image processing toolbox:

Squelch=149     % play with Squelch values, recommended range [120 240], to include or exclude particles not right on 'surface' of image
A=imread('counting particles.jpg')
A2=A(:,:,1)+A(:,:,2)+A(:,:,3)    % a bit more contrast

[A_sz1,A_sz2]=size(A2)   % sz2=571 is x, sz1=574 is y

figure(1);imshow(A)

figure(1) is your original image, no need to repeat it

figure(2);imshow(A2)

observe slightly increased contrast, and from a different point of view you can see the particles you want to count as peaks:

Now let's make the peaks sharper to be able to count them:

[X,Y]=meshgrid([1:1:A_sz2],[1:1:A_sz1])
A3=surface(X,Y,A2)
A4=A3.ZData
[A3_sz1,A3_sz2]=size(A3.ZData)
spots=zeros(A3_sz1,A3_sz2)
for k=1:1:A3_sz1   % sz1, amount of lines, sz1 is y axis
    [pks,locs,widths,proms]=findpeaks(double(A4(k,:)),'MinPeakHeight',Squelch)
%     figure(3);plot(A4(k,:));grid on;hold on;plot(locs,pks,'rx')  
    spots(k,locs)=255
end

Basically, the for loop catches peaks line by line. When plotting all the lines watching sideways we don't really see anything, like a thick bush, but for instance 20 lines sideways look like this:

Note that the function findpeaks without the option 'MinPeakHeight' set to Squelch=149 or higher, for this picture, it takes a lot of noise. Since the illumination conditions are set, or we assume you do not have control over how much light or the kind of light you use, you don't really want to count the blurred background spots that could be anything else.

If you really want to count all the particles, you have illuminate the volume with enough power to 'see' all the contours, even if overlapped, and then process the image to discern them all. If not enough light, we have to carry on eliminating small peaks too close to noise.

Play with Squelch to see what happens, I already did, Squelch<100 will catch a lot of noise

figure(6);surf(X,Y,spots,gradient(spots))

Now the peaks are sharper, but there are particles that have different isolated white pixels counting as separate particles. Let's merge close enough sharp peaks into single peaks with the following:

    F=griddedInterpolant(X',Y',spots','cubic')
    [K,H]=surfature(X,Y,spots)

support function surfature is available from MATLAB CENTRAL here:

<http://uk.mathworks.com/matlabcentral/fileexchange/11168-surface-curvature?s_tid=srchtitle?

There are other ways to 'drop' an envelope on really close each other peaks so they are not counted as separate peaks. Human eye can tell such dots belong to same particle. I found extrema2 really useful to such purpose:

    [spots_max,imax,zmin,imin] = extrema2(spots)   % put support function extrema.m in same folder
    [x_spots,y_spots]=ind2sub(size(spots),imax)        % getting coordinates of particles

extrema2.m is available from MATLAB CENTRAL, it requires support function extrema.m in same folder

http://uk.mathworks.com/matlabcentral/fileexchange/12275-extrema-m--extrema2-m?s_tid=srchtitle

Now let's count peaks

part_count=size(spots_max)                               % counting particles
part_count=part_count(1)

With Squelch=149 I count 943 particles.

Depending upon your CPU, operative system, .. the for loop may take a minute or so.

I understood the point that you said 64x64 pixel basic block, but I saved the shifting window coding that is easy to write, but takes time writing it, you can implement the 64x64 shifting window if you want, but bear in mind that it will take more than you think solving apparently small issues like avoiding repetition of points.

If you shift the window, instead of moving it 64 to 64, there may be additional delay. I repeat, If you want you can implement a smaller processing window, once you find a way to count the spots in a window regardless of window size, does it make sense?

If you find this answer of any help to solve your question, regarding the particles count, please click on the thumbs-up vote link above, thanks in advance.

I will go through the SNR request next week, approx, ok?

John [e-mail]

Answer 2

[Walter Roberson]

Threshold. Subdivide. Regionprops.

The SNR cannot be calculated for those sub-images, as you have no idea what is signal and what is noise. Calculating SNR requires a version of the signal with no noise, whether that version is measured physically or is generated by a software model.

For example when we were doing NMR or MRS, we could put in a pure water sample: with the known composition, any deviation from a clean signal could only be due to noise.

[John BG]

Walter have a look to the sideways screenshots I have attached below.

It is possible to calculate SNR because there are sharp peaks and the question doesn't say 'some are friends some are foes, have a guess'. What happens is that illumination conditions do not get to the bottom of the volume, and there are many particles that are really close to noise level, and there may be other things that 'look like' particles 'down there', so the only safe thing to do is to draw a threshold, which is what the question is about: counting peaks that are high enough, don't you agree?

[John BG]

ok, perhaps it's air, not water, and tethys is attempting to count floating particles in air sample.

1st, on Earth vacuum is expensive to create.

Do you suggest to have a sealed void transparent tank as a reference?

That reference sample would be? black, pitch black, there you have your reference.

In radio instruments calibration you calibrate short open just before measuring, because the temperature of the instrument and DUT may affect the measurements, yet modern instruments have it all in tables, and most of times you just retrieve reference calibration data made available in the instrument firmware, no need to every measurement plug short/open, or so to speak, hang around with a glass of clean water any more or, a transparent bottle of void.

Regards

John

[e-mail]

[Walter Roberson]

I see no reason in the question to justify an assumption that there is any fluid or any illumination condition issue involved.

Starting from the information that we have been given (nothing more than "image" and "particles"), you have multiplied hypotheses unnecessarily.

• Are the particles 2D or 3D? We do not know.
• Is there a fluid involved? We do not know, and we do not know how thick it might be, and we do not know if the image is being taken from inside the fluid or from outside with a refraction boundary to take into account
• Is the potential fluid diffusing the light in any significant manner? We do not know -- just as we do not know the index of refraction of the potential fluid
• Are the particles generating their own light, or are they being distantly illuminated (e.g., starlight on a dust field in space), or are they being illuminated by room lights, or are they being illuminated by an instrument light specifically for that purpose? We do not know -- and we do not know the spectral properties of the potentially illumination either, and those make a difference if we are contenting with any significant refraction.
• Are the wavelengths involved optical? Is, for example, this an x-ray image? We do not know.
• Are the particles firmly enough edged that we can study distributions of shapes and sizes to find "expected" size and shape and use that to assume boundaries between touching particles in order to count them individually? Or are the particles soft boundary (e.g., droplets of liquid) so each blob should be counted as one? We do not know.
• What, in all of this, is leading to noise anyhow? Refraction and diffraction do not cause noise, they cause blurring. Considering the word "image" was used, we can probably assume that shot noise exists in the system, but do we need to worry about dark current (thermal) noise in this system, or... ? http://qsimaging.com/ccd_noise.html Was it a fast exposure, slow exposure? Do we need to be taking into account Brownian movement of the particles?

You simply cannot estimate noise from an image (or sub-image) in isolation. You might be able to estimate noise if given a sufficiently detailed model of the system -- but "sufficiently detailed" in the case of dark current / thermal noise includes a pretty recent reading taken with no data source (because it either is absent or access to it has been cut off) in order to measure sufficiently recently coefficients that have quantum sources and so might have long term statistical properties but whose values now are uncertain enough to make significant differences in the measurements.

All that stuff you posted... it might be relevant to some set-ups, but we have no idea if it is relevant to this set-up.

Answer 3

[Brattv]

Hi, The Mathworks tutorial called Correcting Nonuniform Illumination might have some good answers to several of your questions. Best of luck!

[Image Analyst]

Theoretically they should not have used background subtraction to correct for the non-uniform background. It's not x-rays or fluorescence so they should have used background division not subtraction to flatten the background. Think about it - if you have 90% as many photon coming from the corner of your image, you need to divide the darker/bad value by 0.90 to bring it up to what it should be.