C++, gamma correction

This does a brightness adjustment of the image using a simple gamma correction, with the gamma value determined separately for each component to match the target average.

The high level steps are:

1. Read image and extract histogram for each color component.
2. Perform a binary search of the gamma value for each component. A binary search is performed on the gamma values, until the resulting histogram has the desired average.
3. Read the image a second time, and apply the gamma correction.

All image input/output uses PPM files in ASCII. Images were converted from/to PNG using GIMP. The code was run on a Mac, image conversions were done on Windows.

Code:

#include <cmath>
#include <string>
#include <vector>
#include <sstream>
#include <fstream>
#include <iostream>

static inline int mapVal(int val, float gamma)
{
    float relVal = (val + 1.0f) / 257.0f;
    float newRelVal = powf(relVal, gamma);

    int newVal = static_cast<int>(newRelVal * 257.0f - 0.5f);
    if (newVal < 0)
    {
        newVal = 0;
    }
    else if (newVal > 255)
    {
        newVal = 255;
    }

    return newVal;
}

struct Histogram
{
    Histogram();

    bool read(const std::string fileName);
    int getAvg(int colIdx) const;
    void adjust(const Histogram& origHist, int colIdx, float gamma);

    int pixCount;
    std::vector<int> freqA[3];
};

Histogram::Histogram()
  : pixCount(0)
{
    for (int iCol = 0; iCol < 3; ++iCol)
    {
        freqA[iCol].resize(256, 0);
    }
}

bool Histogram::read(const std::string fileName)
{
    for (int iCol = 0; iCol < 3; ++iCol)
    {
        freqA[iCol].assign(256, 0);
    }

    std::ifstream inStrm(fileName);

    std::string format;
    inStrm >> format;
    if (format != "P3")
    {
        std::cerr << "invalid PPM header" << std::endl;
        return false;
    }

    int w = 0, h = 0;
    inStrm >> w >> h;
    if (w <= 0 || h <= 0)
    {
        std::cerr << "invalid size" << std::endl;
        return false;
    }

    int maxVal = 0;
    inStrm >> maxVal;
    if (maxVal != 255)
    {
        std::cerr << "invalid max value (255 expected)" << std::endl;
        return false;
    }

    pixCount = w * h;

    int sumR = 0, sumG = 0, sumB = 0;
    for (int iPix = 0; iPix < pixCount; ++iPix)
    {
        int r = 0, g = 0, b = 0;
        inStrm >> r >> g >> b;
        ++freqA[0][r];
        ++freqA[1][g];
        ++freqA[2][b];
    }

    return true;
}

int Histogram::getAvg(int colIdx) const
{
    int avg = 0;
    for (int val = 0; val < 256; ++val)
    {
        avg += freqA[colIdx][val] * val;
    }

    return avg / pixCount;
}

void Histogram::adjust(const Histogram& origHist, int colIdx, float gamma)
{
    freqA[colIdx].assign(256, 0);

    for (int val = 0; val < 256; ++val)
    {
        int newVal = mapVal(val, gamma);
        freqA[colIdx][newVal] += origHist.freqA[colIdx][val];
    }
}

void mapImage(const std::string fileName, float gammaA[])
{
    std::ifstream inStrm(fileName);

    std::string format;
    inStrm >> format;

    int w = 0, h = 0;
    inStrm >> w >> h;

    int maxVal = 0;
    inStrm >> maxVal;

    std::cout << "P3" << std::endl;
    std::cout << w << " " << h << std::endl;
    std::cout << "255" << std::endl;

    int nPix = w * h;

    for (int iPix = 0; iPix < nPix; ++iPix)
    {
        int inRgb[3] = {0};
        inStrm >> inRgb[0] >> inRgb[1] >> inRgb[2];

        int outRgb[3] = {0};
        for (int iCol = 0; iCol < 3; ++iCol)
        {
            outRgb[iCol] = mapVal(inRgb[iCol], gammaA[iCol]);
        }

        std::cout << outRgb[0] << " " << outRgb[1] << " "
                  << outRgb[2] << std::endl;
    }
}

int main(int argc, char* argv[])
{
    if (argc < 5)
    {
        std::cerr << "usage: " << argv[0]
                  << " ppmFileName targetR targetG targetB"
                  << std::endl;
        return 1;
    }

    std::string inFileName = argv[1];

    int targAvg[3] = {0};
    std::istringstream strmR(argv[2]);
    strmR >> targAvg[0];
    std::istringstream strmG(argv[3]);
    strmG >> targAvg[1];
    std::istringstream strmB(argv[4]);
    strmB >> targAvg[2];

    Histogram origHist;
    if (!origHist.read(inFileName))
    {
        return 1;
    }

    Histogram newHist(origHist);
    float gammaA[3] = {0.0f};

    for (int iCol = 0; iCol < 3; ++iCol)
    {
        float minGamma = 0.0f;
        float maxGamma = 1.0f;
        for (;;)
        {
            newHist.adjust(origHist, iCol, maxGamma);
            int avg = newHist.getAvg(iCol);
            if (avg <= targAvg[iCol])
            {
                break;
            }
            maxGamma *= 2.0f;
        }

        for (;;)
        {
            float midGamma = 0.5f * (minGamma + maxGamma);

            newHist.adjust(origHist, iCol, midGamma);
            int avg = newHist.getAvg(iCol);
            if (avg < targAvg[iCol])
            {
                maxGamma = midGamma;
            }
            else if (avg > targAvg[iCol])
            {
                minGamma = midGamma;
            }
            else
            {
                gammaA[iCol] = midGamma;
                break;
            }
        }
    }

    mapImage(inFileName, gammaA);

    return 0;
}

The code itself is fairly straightforward. One subtle but important detail is that, while the color values are in the range [0, 255], I map them to the gamma curve as if the range were [-1, 256]. This allows the average to be forced to 0 or 255. Otherwise, 0 would always remain 0, and 255 would always remain 255, which might never allow for an average of 0/255.

To use:

1. Save the code in a file with extension .cpp, e.g. force.cpp.
2. Compile with c++ -o force -O2 force.cpp.
3. Run with ./force input.ppm targetR targetG target >output.ppm.

Sample output for 150, 100, 100:

Other samples are above SE image size limit.

Sample output for 75, 91, 110:

Other samples are above SE image size limit.