How to find number of intensity levels in an image

I'm working with SEM images, for finding pore size distribution, I need to provide Number of intensity levels as an input, so how to determine it?

Answer 1

[DGM]

You can use unique() to find the number of unique gray levels if that's what you're after. Consider the following image with 8 gray levels:

inpict = imread('grad8-1.png');
nu = numel(unique(inpict)) 

nu = 8

But be aware that things like antialiasing and compression may mean that the number of nominal gray levels isn't the same as the actual number of unique values in the image.

inpict = imread('grad8-2.png'); % antialiased
nu = numel(unique(inpict)) 

nu = 114

inpict = imread('grad8-1.jpg'); % lossy compresssion
nu = numel(unique(inpict)) 

nu = 228

[William Rose]

@Deep, @DGM, very nice answer.

[Walter Roberson]

JPEG images that have distinct boundaries are almost always smoothed out near the boundaries, introducing unique intensity values that are not present in the original. If you need to worry about the number of distinct intensity levels, then you should actively avoid JPEG.

[DGM]

I was thinking about this while working, but I guess IA already beat me to it. As I said, if for some reason you want to find the number of nominal gray levels in a quantized image like those above, using unique() won't necessarily give you the right answer. You might try looking at the distribution of gray levels to see if you can find a pattern. Consider the similar 32-level gray ramp, as subject to JPG compression:

A = imread('grad32-1.jpg'); % an image (grayscale, uint8)

[counts,~] = imhist(A); % get the bin counts
counts = padarray(counts,[1 0],0,'both'); % end bins need neighbors
hs = stem(counts,'marker','none'); hold on % plot the histogram

% try to find the peaks (this doesn't always work)
[pk pkidx] = findpeaks(counts,'MinPeakProminence',20);
plot(pkidx,pk,'bo')

% compare the two methods of counting gray levels
nuniqlevels = numel(unique(A))

nuniqlevels = 256

ngraylevels = numel(pk)

ngraylevels = 32

So in this example, we were able to correctly estimate that the image had been quantized to 32 levels, despite all of the noise. That said, I doubt this would reliably work in practice. The test image has relatively few gray levels, minimal edge-affected area, and the quantization is uniform. Perhaps most importantly, each gray level is represented by a relatively large portion of the image area (despite the fact that they are not equal areas).

Consider the following two versions of cameraman.tif, both quantized to 32 levels. Unlike the example above, the number of pixels in each bin varies significantly, and there is a much larger portion of the image which is now adjacent to an edge feature between two gray levels.

If the image is clean, we can still get the number of gray levels as before:

figure
A = imread('cam32.png'); % an image (grayscale, uint8)

[counts,~] = imhist(A); % get the bin counts
counts = padarray(counts,[1 0],0,'both'); % end bins need neighbors
hs = stem(counts,'marker','none'); hold on % plot the histogram

% try to find the peaks (this doesn't always work)
[pk pkidx] = findpeaks(counts,'MinPeakProminence',20);
plot(pkidx,pk,'bo')

% correct estimates are easier on clean images
ngraylevels = numel(pk)

ngraylevels = 32

Note how there is relatively little content in the bins in the mid-dark and bright regions. If anything were to spread the peaks, they wouldn't be discernable anymore. What happens if this image is smoothed or compressed?

figure
A = imread('cam32.jpg'); % an image (grayscale, uint8)

[counts,~] = imhist(A); % get the bin counts
counts = padarray(counts,[1 0],0,'both'); % end bins need neighbors
hs = stem(counts,'marker','none'); hold on % plot the histogram

% try to find the peaks (this doesn't always work)
[pk pkidx] = findpeaks(counts,'MinPeakProminence',20);
plot(pkidx,pk,'bo')

% not correct
ngraylevels = numel(pk)

ngraylevels = 23

I'm still not sure if this is really part of the needed solution. I'm inclined to doubt that your images are 32-level uint8 images or something similar, so maybe it's worth asking how the images are being presented (their class, data range, etc), what "gray levels" means to you in the context of your work, and how you intend to use that information.

[Deep]

Thank you for the detailed explanation, it also helped me clear my confusion about graylevel intensity

[Image Analyst]

OK, can you give us some context as to why you think you need this? What will that enable you to do with your SEM images that you can't do unless you have the number of unique gray levels?

[Deep]

Actually, I came across an open source MATLAB code for finding SEM Image Porosity and Pore Size, In that code, Number of intensity level was asked as an input, Here's a small snippet of it

[Image Analyst]

You can see that the number N is not used in what you've shown so far. And I highly doubt it will be needed to find porosity and pore size.

[Deep]

That image was only a snapshot of a small part of code, where it was asked as input

By changing the value of N, porosity changes significantly, that's why it was neccesary to find it

Note: below attached is also a small snapshot of the code

[Image Analyst]

@Deep you totally misunderstood what the N in the multithresh means.

help multithresh

 MULTITHRESH Multi-level image thresholds using Otsu's method.
    THRESH = MULTITHRESH(A) computes a single threshold for image A using
    Otsu's method and returns it in THRESH. THRESH can be used to convert A
    into a two-level image using IMQUANTIZE.
 
    THRESH = MULTITHRESH(A, N) computes N thresholds for image A using the
    Otsu's method and returns them in THRESH. THRESH is a 1xN vector which
    can be used to convert A into an image with (N+1) discrete levels using
    IMQUANTIZE. The maximum value allowed for N is 20.
 
    [THRESH, METRIC] = MULTITHRESH(A,...) returns the effectiveness metric
    as the second output argument. METRIC is in the range [0 1] and a
    higher value indicates greater effectiveness of the thresholds in
    separating the input image into N+1 regions based on Otsu's objective
    criterion.
 
    Class Support
    -------------
    The input image A is an array of one of the following classes: uint8,
    uint16, int16, single, or double. It must be nonsparse. N is a positive
    integer between 1 and 20, and can be any numeric class. THRESH is a 1xN
    numeric vector of the same class as A. METRIC is a scalar of class
    double.
 
    Notes
    -----
    1. A can be an array of any dimension. MULTITHRESH finds the thresholds
       based on the aggregate histogram of the entire array. An RGB image
       is considered a 3D numeric array and the thresholds are computed
       for the combined data from all the three color planes.
 
    2. MULTITHRESH uses the range of the input A [min(A(:)) max(A(:))] as
       the limits for computing image histogram which is used in subsequent
       computations. Any NaNs in A are ignored in computation. Any Infs and
       -Infs in A are counted in the first and last bin of the histogram,
       respectively.
 
    3. For N > 2, MULTITHRESH uses search-based optimization of Otsu's
       criterion to find the thresholds. The search-based optimization
       guarantees only locally optimal results. Since the chance of
       converging to local optimum increases with N, it is preferable to
       use smaller values of N, typically N < 10.
 
    4. For degenerate inputs where the number of unique values in A is less
       than or equal to N, there is no viable solution using Otsu's method.
       For such inputs, MULTITHRESH returns THRESH such that it includes
       all the unique values from A and possibly some extra values that are
       chosen arbitrarily.
 
    5. The thresholds (THRESH) returned by MULTITHRESH are in the same
       range as the input image A. This is unlike GRAYTHRESH which returns
       a normalized threshold in the range [0, 1].
 
    Reference
    ---------
    N. Otsu, "A Threshold Selection Method from Gray-Level Histograms," IEEE
    Transactions on Systems, Man, and Cybernetics, Vol. 9, No. 1, pp. 62-66,
    1979.
 
    Example
    ---------
      % This example computes multiple thresholds for an image using
      % MULTITHRESH and applies those thresholds to the image using IMQUANTIZE
      % to get segment labels.
 
      I = imread('circlesBrightDark.png');
      imshow(I, [])
      title('Original Image');
 
      % Compute the thresholds
      thresh = multithresh(I,2);
 
      % Apply the thresholds to obtain segmented image
      seg_I = imquantize(I,thresh);
 
      % Show the various segments in the segmented image in color
      RGB = label2rgb(seg_I);
      figure, imshow(RGB)
      title('Segmented Image');
 
    See also GRAYTHRESH, IMBINARIZE, IMQUANTIZE, RGB2IND.

    Documentation for multithresh
       doc multithresh

As you can see, N+1 is the number of ranges you want to split your historgam up into. It returns the thresholds (dividing gray levels) between the ranges. So it is NOT the number of unique gray levels in your image. So I'm standing firm on my answer that the number of unique gray levels in your image is not something you need to compute according to the code you posted.

FWIW, this thread is a classic demonstration of the XY Problem.

https://en.wikipedia.org/wiki/XY_problem

In other words you thought you knew what you needed and asked for it, and people gave you that but once the real situation was explained, it turned out you didn't want what you asked for at all.

So just forget what everyone told you about computing the number of gray levels and just loop with N = 1 to 20 (the max multithresh allows) and see how your porosity changes.

[Deep]

Thank you for clearing the concept, I was also bit confused about the N value, because that too depends on the image.

I really appreciate that you came across the concept and cleared the concept

[William Rose]

@Image Analyst,

Thanks for explaining the XY problem clearly and for the link. That is useful concept. I will keep it in mind when I am seeking or trying to give help.

[Image Analyst]

I'm not sure why you need multiple thresholds (because I haven't see the image). It could be that you just need one set of thresholds - a lower one and an upper one - and then your mask is everything in between. SEM images are notoriously hard to threshold and I'm not sure any automatic method would be good. You may want to use an interactive one like my utility:

https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

It could be, if you're lucky, that you ask for 5 thresholds and pick the intensity range between the second and third threshold to mask out the stuff you want (the pores).

You might want to learn about granulometry:

https://en.wikipedia.org/wiki/Granulometry_(morphology)

Answer 2

[Image Analyst]

Just mentioning imhist in case you wanted the number of pixels at each gray level (distribution) instead the number of unique gray levels (of which I don't think I've ever used and can't imagine why that would even be necessary).

[Walter Roberson]

The number of distinct intensity levels can be used as part of huffman encoding.