一、简介

这里使用一种简单的方式来计算图像中的像素值直方图分布。计算过程如下所述：

第一种方式：

1、首先将图像变为一维数组（reshape），并将数组中的数值进行排序。
2、基于排序的连续性，统计相同的像素值个数，直到最后一个像素被统计，则计算结束。、

第二种方式：

这种方式可能更简单一些，借用std中的unorder_map数据结构，来完成图像像素值的直方图统计。

二、实现代码

第一种方式：

//标准文件
#include <errno.h>
#include <iostream>
#include <fstream>
#include <queue>
#include <map>

//GDAL
#include "gdal_priv.h"

//Eigen
#include <Eigen/Dense>

typedef Eigen::MatrixXd Matrix;

//图像直方图统计
void HistogramStatistics(const Matrix& images, std::map<double,int>& hists)
{
    if (!images.size()) return;

    Eigen::VectorXd imgVec = images.reshaped();
    size_t imgCount = imgVec.size();
    Eigen::VectorXi indices = Eigen::VectorXi::LinSpaced(Eigen::Sequential, imgCount, 0, imgCount);
    std::sort(indices.data(), indices.data() + indices.size(), [&imgVec](int i1, int i2)
        { return imgVec[i1] < imgVec[i2]; });

    int label = 0;
    float curVal = imgVec(indices(0));
    hists[curVal] = 1;
    for (int i = 1; i < indices.size(); ++i)
    {
        if (imgVec(indices(i)) == curVal)
        {
            hists[curVal]++;
            continue;
        }
        else
        {
            curVal = imgVec(indices(i));
            hists[curVal] = 1;
        }
    }
}

int main(int argc, const char* argv[])
{
    //--------------------------------获取数据集-----------------------------
    std::cout << "--------------------------------获取数据集-----------------------------" << std::endl;
    if (argc != 2) {
        return EINVAL;
    }
    //它可能是一个URL，一个文件名等存储事物，它依赖于驱动程序
    const char* pszFilename = argv[1];

    GDALDatasetUniquePtr poDataset;     //该指针已经包含GDALClose函数
    GDALAllRegister();
    const GDALAccess eAccess = GA_ReadOnly;
    //如果GDALOpen()返回NULL，则意味着打开失败，并且已经通过CPLError()发出了错误消息。
    poDataset = GDALDatasetUniquePtr(GDALDataset::FromHandle(GDALOpen(pszFilename, eAccess)));
    if (!poDataset)
    {
        std::cerr << "read failure！！！" << std::endl;
        //GDALClose((GDALDatasetH)poDataset.get());
        return -1;
    }
    std::cout << "read successfully！！！" << std::endl;

    //--------------------------------获取栅格数据-----------------------------
    std::cout << "--------------------------------获取栅格数据-----------------------------" << std::endl;

    //第一种方式
    double* pafScanline;
    GDALRasterBand* poBand;
    poBand = poDataset->GetRasterBand(1);
    int   nXSize = poBand->GetXSize();
    int   nYSize = poBand->GetYSize();
    pafScanline = (double*)CPLMalloc(sizeof(double) * nXSize * nYSize);
    poBand->RasterIO(GF_Read, 0, 0, nXSize, nYSize,
        pafScanline, nXSize, nYSize, GDT_Float64,
        0, 0);

    Matrix img(nXSize, nYSize);
    for (int i = 0; i < nXSize; ++i)        //行
    {
        for (int j = 0; j < nYSize; ++j)         //列
        {
            int index = i * nXSize + j;
            img(i, j) = pafScanline[index];
        }
    }
    CPLFree(pafScanline);       //释放缓冲区

    std::map<double, int> hists;
    HistogramStatistics(img,hists);

    //--------------------------------输出结果-----------------------------
    std::cout << "--------------------------------输出结果-----------------------------" << std::endl;
    for (auto& it : hists)
    {
        std::cout << it.first << ":" << it.second << std::endl;
    }


    std::cout << "--------------------------------输出结果-----------------------------" << std::endl;
    std::cout << "compute successfully！！！" << std::endl;

    return 0;
}

第二种方式：

//标准文件
#include <errno.h>
#include <iostream>
#include <fstream>
#include <queue>
#include <map>
#include <unordered_map>

//GDAL
#include "gdal_priv.h"

//Eigen
#include <Eigen/Dense>

typedef Eigen::MatrixXd Matrix;

//图像直方图统计
void HistogramStatistics(const Matrix& images, std::map<double,int>& hists)
{
    if (!images.size()) return;

    Eigen::VectorXd imgVec = images.reshaped();
    size_t imgCount = imgVec.size();
    Eigen::VectorXi indices = Eigen::VectorXi::LinSpaced(Eigen::Sequential, imgCount, 0, imgCount);
    std::sort(indices.data(), indices.data() + indices.size(), [&imgVec](int i1, int i2)
        { return imgVec[i1] < imgVec[i2]; });

    int label = 0;
    float curVal = imgVec(indices(0));
    hists[curVal] = 1;
    for (int i = 1; i < indices.size(); ++i)
    {
        if (imgVec(indices(i)) == curVal)
        {
            hists[curVal]++;
            continue;
        }
        else
        {
            curVal = imgVec(indices(i));
            hists[curVal] = 1;
        }
    }
}

int main(int argc, const char* argv[])
{
    //--------------------------------获取数据集-----------------------------
    std::cout << "--------------------------------获取数据集-----------------------------" << std::endl;
    if (argc != 2) {
        return EINVAL;
    }
    //它可能是一个URL，一个文件名等存储事物，它依赖于驱动程序
    const char* pszFilename = argv[1];

    GDALDatasetUniquePtr poDataset;     //该指针已经包含GDALClose函数
    GDALAllRegister();
    const GDALAccess eAccess = GA_ReadOnly;
    //如果GDALOpen()返回NULL，则意味着打开失败，并且已经通过CPLError()发出了错误消息。
    poDataset = GDALDatasetUniquePtr(GDALDataset::FromHandle(GDALOpen(pszFilename, eAccess)));
    if (!poDataset)
    {
        std::cerr << "read failure！！！" << std::endl;
        //GDALClose((GDALDatasetH)poDataset.get());
        return -1;
    }
    std::cout << "read successfully！！！" << std::endl;

    //--------------------------------获取栅格数据-----------------------------
    std::cout << "--------------------------------获取栅格数据-----------------------------" << std::endl;

    //第一种方式
    double* pafScanline;
    GDALRasterBand* poBand;
    poBand = poDataset->GetRasterBand(1);
    int   nXSize = poBand->GetXSize();
    int   nYSize = poBand->GetYSize();
    pafScanline = (double*)CPLMalloc(sizeof(double) * nXSize * nYSize);
    poBand->RasterIO(GF_Read, 0, 0, nXSize, nYSize,
        pafScanline, nXSize, nYSize, GDT_Float64,
        0, 0);

    std::unordered_map<double, int> hists;
    for (int i = 0; i < nXSize; ++i)        //行
    {
        for (int j = 0; j < nYSize; ++j)         //列
        {
            int index = i * nXSize + j;
            hists[pafScanline[index]]++;
        }
    }
    CPLFree(pafScanline);       //释放缓冲区

    //--------------------------------输出结果-----------------------------
    std::cout << "--------------------------------输出结果-----------------------------" << std::endl;
    for (auto& it : hists)
    {
        std::cout << it.first << ":" << it.second << std::endl;
    }


    std::cout << "--------------------------------输出结果-----------------------------" << std::endl;
    std::cout << "compute successfully！！！" << std::endl;

    return 0;
}