How to 'normalize' a grayscale image?

Question

My math is a bit rusty. I'm trying to equalize a histogram of a 2D array which represents grayscale values in the range 0-255 (values may not be whole numbers because of how they are computed).

I found this article on Wikipedia, but I don't quite understand the formulas they present.

ni, n and L I can compute, but I'm not quite sure how to implement this cdf function. Might this function be of use?

Here's what I've got so far:

static double[,] Normalize(double[,] mat)
{
    int width = mat.GetLength(0);
    int height = mat.GetLength(1);
    int nPixels = width*height;

    double sum = 0;
    double max = double.MinValue;
    double min = double.MaxValue;
    var grayLevels = new Dictionary<double, int>();

    foreach (var g in mat)
    {
        sum += g;
        if (g > max) max = g;
        if (g < min) min = g;
        if (!grayLevels.ContainsKey(g)) grayLevels[g] = 0;
        ++grayLevels[g];
    }
    double avg = sum/nPixels;
    double range = max - min;

    var I = new double[width,height];

    // how to normalize?

    return I;
}

Solution 1

You can use this function I just wrote:

public static Bitmap ContrastStretch(Bitmap srcImage, double blackPointPercent = 0.02, double whitePointPercent = 0.01)
{
    BitmapData srcData = srcImage.LockBits(new Rectangle(0, 0, srcImage.Width, srcImage.Height), ImageLockMode.ReadOnly,
        PixelFormat.Format32bppArgb);
    Bitmap destImage = new Bitmap(srcImage.Width, srcImage.Height);
    BitmapData destData = destImage.LockBits(new Rectangle(0, 0, destImage.Width, destImage.Height),
        ImageLockMode.WriteOnly, PixelFormat.Format32bppArgb);
    int stride = srcData.Stride;
    IntPtr srcScan0 = srcData.Scan0;
    IntPtr destScan0 = destData.Scan0;
    var freq = new int[256];

    unsafe
    {
        byte* src = (byte*) srcScan0;
        for (int y = 0; y < srcImage.Height; ++y)
        {
            for (int x = 0; x < srcImage.Width; ++x)
            {
                ++freq[src[y*stride + x*4]];
            }
        }

        int numPixels = srcImage.Width*srcImage.Height;
        int minI = 0;
        var blackPixels = numPixels*blackPointPercent;
        int accum = 0;

        while (minI < 255)
        {
            accum += freq[minI];
            if (accum > blackPixels) break;
            ++minI;
        }

        int maxI = 255;
        var whitePixels = numPixels*whitePointPercent;
        accum = 0;

        while (maxI > 0)
        {
            accum += freq[maxI];
            if (accum > whitePixels) break;
            --maxI;
        }
        double spread = 255d/(maxI - minI);
        byte* dst = (byte*) destScan0;
        for (int y = 0; y < srcImage.Height; ++y)
        {
            for (int x = 0; x < srcImage.Width; ++x)
            {
                int i = y*stride + x*4;

                byte val = (byte) Clamp(Math.Round((src[i] - minI)*spread), 0, 255);
                dst[i] = val;
                dst[i + 1] = val;
                dst[i + 2] = val;
                dst[i + 3] = 255;
            }
        }
    }

    srcImage.UnlockBits(srcData);
    destImage.UnlockBits(destData);

    return destImage;
}

static double Clamp(double val, double min, double max)
{
    return Math.Min(Math.Max(val, min), max);
}

The defaults mean that the darkest 2% of pixels will become black, the lightest 1% will become white, and everything in between will be stretched to fill the color space. This is the same as the default for ImageMagick.

This algorithm has the fun side effect that if you use values above 50% then it will invert the image! Set to .5, .5 to get a black & white image (2 shades) or 1, 1 to get a perfect inversion.

Assumes your image is already grayscale.

Solution 2

Found something that you might find useful

http://sonabstudios.blogspot.in/2011/01/histogram-equalization-algorithm.html

Hope that helps

Solution 3

Calculating the cumulative distribution function involves a couple of steps.

First you get the frequency distribution of your grayscale values.

So something like:

freqDist = new int[256];

for each (var g in mat)
{
    int grayscaleInt = (int)g;
    freqDist[grayscaleInt]++;
}

Then to get your CDF, something like:

cdf = new int[256];
int total = 0;

for (int i = 0; i < 256; i++)
{
    total += freqDist[i];
    cdf[i] = total;
}

Solution 4

I can help you to understand your link,

first, counting value which represent image, shows in that link,

Value   Count   Value   Count   Value   Count   Value   Count   Value   Count
   52       1      64       2      72       1      85       2     113       1
   55       3      65       3      73       2      87       1     122       1
   58       2      66       2      75       1      88       1     126       1
   59       3      67       1      76       1      90       1     144       1
   60       1      68       5      77       1      94       1     154       1
   61       4      69       3      78       1     104       2   
   62       1      70       4      79       2     106       1
   63       2      71       2      83       1     109       1

it is means, the image is created with values above, nothing else.

second, sums the value cumulatively from 52 to 154

Value   cdf Value   cdf Value   cdf Value   cdf Value   cdf
   52     1    64    19    72    40    85    51   113    60
   55     4    65    22    73    42    87    52   122    61
   58     6    66    24    75    43    88    53   126    62
   59     9    67    25    76    44    90    54   144    63
   60    10    68    30    77    45    94    55   154    64
   61    14    69    33    78    46   104    57 
   62    15    70    37    79    48   106    58
   63    17    71    39    83    49   109    59

it is means,

value 52 have 1 cdf cause it is initial value, 
value 55 have 4 cdf cause it has 3 count in image plus 1 cdf from 52, 
value 58 have 6 cdf cause it has 2 count in image plus 4 cdf from 55,
and so on.. till..
value 154 have 64 cdf cause it has 1 count in image plus 63 cdf from 144.

then, calculating histogram equalization formula for each image values based on the function

cdf(v) is represent current cdf value from current image value,

in this case, if h(v) = 61 so cdf(v) = 14

cdfmin is represent initial cdf value, in this case, 1 cdf from value 52

happy coding.. ^^

Solution 5

Here's my implementation:

private static byte[,] Normalize(byte[,] mat)
{
    int width = mat.GetLength(0);
    int height = mat.GetLength(1);
    int nPixels = width*height;

    var freqDist = new int[256];

    foreach (var g in mat)
    {
        ++freqDist[g];
    }

    var cdf = new int[256];
    int total = 0;

    for (int i = 0; i < 256; ++i)
    {
        total += freqDist[i];
        cdf[i] = total;
    }

    int cdfmin = 0;

    for (int i = 0; i < 256; ++i)
    {
        if (cdf[i] > 0)
        {
            cdfmin = cdf[i];
            break;
        }
    }

    var I = new byte[width,height];
    double div = (nPixels - cdfmin) / 255d;

    for (int y = 0; y < height; ++y)
    {
        for (int x = 0; x < width; ++x)
        {
            I[x, y] = (byte)Math.Round((cdf[mat[x, y]] - cdfmin) / div);
        }
    }

    return I;
}

I changed it from using doubles to bytes to work better with the histogram (freqDist).

View more solutions

Comments

[Blender]

This is a comment, not an answer.

[hanut]

I cant comment yet but I wanted to help.

[mpen]

Both of your answers were helpful, but yours was the trickier part. Got it implemented now, seems to be working well. Thank you!

[Admin]

Thanks for posting back your results. Any reason beyond convenience for the unsafe scope?

[mpen]

Speed. This was a long time ago, but I believe it's very slow to access pixel data in a safe context.