Determine the intensity value, "T", in a 2D image for which 99.9% of all intensity values are less than "T"

% I can read in an image and calculate a threshold for which the pixels are above some value, i.e., 6000.

Threshold=Q>6000;

% I can count the number of those pixels above that value

PixelCount = sum(Threshold(:))

% ratio the PixelCount over the total # of pixels in an image with dimensions X1, Y1 to determine a fraction, or percent of pixels above that value

Fraction=PixelCount/(X1*Y1)

But if I wish to calculate the intensity value that represents a desired fraction, i.e. 0.001, - instead of choosing an arbitrary number, like 6000........... how do I program to extract the intensity value that represents cutoff for which 99.9% of pixels have a lower intensity ?

I should mention these are uint16 images.

[Ed Principe]

I think I have something that works for me..... Q is an image from an image stack extracted in a for loop wih index 'p':

Q=l.images1(:,:,p);

SortedQ=sort(Q(:));

figure;plot(SortedQ);

title(['Sorted Pixel Intensities:',F,': image #: ', num2str(p),'/',num2str(l.noimages)],'Interpreter', 'none');

SortedQ(round(0.999*X1*Y1))

Answer 1

[Akira Agata]

I believe prctile function will be helpful to this task, like:

% Read sample gray scale image
I = imread('cameraman.tif');

% Calculate threshold value (99.9% cut-off)
th = prctile(I(:),99.9);

Check the result:

>> nnz(I<th)/numel(I)

ans =

0.9990

[Ed Principe]

Thanks! But It appears that function requires both the Statistics and Machine Learning toolboxes. :-(

[Akira Agata]

OK, then how about the following workaround?

Isort = sort(I(:));
pt = round(numel(x)*0.999);
th = Isort(pt);

If the numel(I) (= total number of pixel) is sufficiently large, the result will be the same.

[Ed Principe]

Yes, that works, that is essentially what I came up with in my comment on the original post above. Thank you very much!!