彩色图像直方图均衡——HSI模型

前言

参考书籍:《数字图像处理(第三版)》冈萨雷斯
参考博客: 参考博客
限定范围: 处理 24 位色深 BMP 格式文件

RGB转HSI公式

RGB 值先除以 255 255 255 ,转到 [ 0 , 1 ] [0,1] [0,1] 区间。
θ = arccos ⁡ 1 2 [ ( R − G ) + ( R − B ) ] [ ( R − G ) 2 + ( R − B ) ( G − B ) ] 1 / 2 H = { θ , B ≤ G 2 π − θ , B > G S = 1 − 3 R + G + B [ min ⁡ ( R , G , B ) ] I = 1 3 ( R + G + B ) \begin{aligned} \theta&=\arccos\frac{\displaystyle\frac 12[(R-G)+(R-B)]}{[(R-G)^2+(R-B)(G-B)]^{1/2}}\\ H&=\begin{cases} \theta,&B\leq G\\ 2\pi-\theta,&B>G \end{cases}\\ S&=1-\frac{3}{R+G+B}[\min(R,G,B)]\\ I&=\frac 13(R+G+B) \end{aligned} θHSI=arccos[(RG)2+(RB)(GB)]1/221[(RG)+(RB)]={θ,2πθ,BGB>G=1R+G+B3[min(R,G,B)]=31(R+G+B)

HSI转RGB公式

0 ≤ H < 2 π 3 0\leq H< \displaystyle\frac{2\pi}3 0H<32π 时:
B = I ( 1 − S ) R = I [ 1 + S cos ⁡ H c o s ( π 3 − H ) ] G = 3 I − ( R + B ) \begin{aligned} B&=I(1-S)\\ R&=I[1+\frac{S\cos H}{cos(\displaystyle\frac{\pi}3-H)}]\\ G&=3I-(R+B) \end{aligned} BRG=I(1S)=I[1+cos(3πH)ScosH]=3I(R+B)
2 π 3 ≤ H < 4 π 3 \displaystyle\frac{2\pi}3 \leq H< \displaystyle\frac{4\pi}3 32πH<34π 时:
H = H − 2 π 3 R = I ( 1 − S ) G = I [ 1 + S cos ⁡ H c o s ( π 3 − H ) ] B = 3 I − ( R + G ) \begin{aligned} H&=H-\displaystyle\frac{2\pi}3\\ R&=I(1-S)\\ G&=I[1+\frac{S\cos H}{cos(\displaystyle\frac{\pi}3-H)}]\\ B&=3I-(R+G) \end{aligned} HRGB=H32π=I(1S)=I[1+cos(3πH)ScosH]=3I(R+G)
4 π 3 ≤ H < 2 π \displaystyle\frac{4\pi}3 \leq H< 2\pi 34πH<2π 时:
H = H − 4 π 3 G = I ( 1 − S ) B = I [ 1 + S cos ⁡ H c o s ( π 3 − H ) ] R = 3 I − ( G + B ) \begin{aligned} H&=H-\displaystyle\frac{4\pi}3\\ G&=I(1-S)\\ B&=I[1+\frac{S\cos H}{cos(\displaystyle\frac{\pi}3-H)}]\\ R&=3I-(G+B) \end{aligned} HGBR=H34π=I(1S)=I[1+cos(3πH)ScosH]=3I(G+B)
转换后的 RGB 值在 [ 0 , 1 ] [0,1] [0,1] 区间,乘以 255 255 255 转到原区间。

代码

#include<iostream>
#include<windows.h>
#include<cmath>
#include<vector>
#define pi acos(-1)
using namespace std;

vector<vector<double>> model;

int main()
{
    //fp为原图像文件指针,fp1为处理后图像文件指针,该文件可不存在
    FILE* fp = fopen("cry.bmp", "rb");
    FILE* fp1 = fopen("test.bmp", "wb");
    if (!fp)
    {
        cout << "原图像文件不存在,请保证原图像与代码在同一目录下" << endl;
        return -1;
    }
    BITMAPFILEHEADER* fileheader = new BITMAPFILEHEADER;
    BITMAPINFOHEADER* infoheader = new BITMAPINFOHEADER;
    fread(fileheader, 14, 1, fp);
    fread(infoheader, 40, 1, fp);
    int bits = infoheader->biBitCount;
    if (bits != 24)
    {
        cout << "该文件不符合要求" << endl;
        return -1;
    }
    int width = infoheader->biWidth, height = infoheader->biHeight;
    int lineByte = (width * 3 + 3) / 4 * 4;
    fileheader->bfSize = 14 + 40 + lineByte * height;
    infoheader->biSizeImage = lineByte * height;
    fwrite(fileheader, 14, 1, fp1);
    fwrite(infoheader, 40, 1, fp1);
    model.resize(height, vector<double>(3 * width));
    unsigned char* origin = new unsigned char[lineByte * height];
    memset(origin, 0, lineByte * height);
    double counter[256] = { 0 };
    double R, G, B;
    double H, S, I;
    int x = 0, y = 0;
    double eps = 1e-6;
    double Smin = 10, Smax = -1, Imin = 10, Imax = -1;
    double Rmin = 10, Rmax = -1, Gmin = 10, Gmax = -1, Bmin = 10, Bmax = -1;
    double low, top, bottom, theta;
    fread(origin, 1, lineByte * height, fp);
    fclose(fp);

    //rgb转hsi
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            model[i][3 * j] = *(origin + i * lineByte + 3 * j + 2);
            model[i][3 * j + 1] = *(origin + i * lineByte + 3 * j + 1);
            model[i][3 * j + 2] = *(origin + i * lineByte + 3 * j);
            R = model[i][3 * j] / 255;
            G = model[i][3 * j + 1] / 255;
            B = model[i][3 * j + 2] / 255;
            model[i][3 * j + 2] = (R + G + B) / 3;
            if (R == G && G == B && B == R) model[i][3 * j + 1] = 0;
            else
            {
                low = min(min(R, G), B);
                model[i][3 * j + 1] = 1 - 3 * low / (R + G + B);
            }
            if (model[i][3 * j + 1] == 0) model[i][3 * j] = 0;
            else {
                top = 0.5 * (2 * R - G - B);
                bottom = sqrt((R - G) * (R - G) + (R - B) * (G - B));
                theta = acos(top / (bottom +eps));
                if (B <= G) model[i][3 * j] = theta;
                else model[i][3 * j] = 2 * pi - theta;
            }
            Smin = min(Smin, model[i][3 * j + 1]);
            Imin = min(Imin, model[i][3 * j + 2]);
        }

    //对s、i进行归一化
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            model[i][3 * j + 1] -= Smin;
            model[i][3 * j + 2] -= Imin;
            Smax = max(Smax, model[i][3 * j + 1]);
            Imax = max(Imax, model[i][3 * j + 2]);
        }
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            if (Smax) model[i][3 * j + 1] /= Smax;
            if (Imax) model[i][3 * j + 2] /= Imax;
            counter[(int)(model[i][3 * j + 2] * 255 + 0.5)]++;
        }

    for (int i = 1;i < 256;i++)
        counter[i] += counter[i - 1];

    //hsi转rgb
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            model[i][3 * j + 2] = counter[(int)(model[i][3 * j + 2] * 255 + 0.5)] / width / height;
            H = model[i][3 * j];
            S = model[i][3 * j + 1];
            I = model[i][3 * j + 2];

            if (H < 2 * pi / 3)
            {
                B = I * (1 - S);
                R = I * (1 + S * cos(H) / cos(pi / 3 - H));
                G = 3 * I - R - B;
            }
            else if (H < 4 * pi / 3)
            {
                H -= 2 * pi / 3;
                R = I * (1 - S);
                G = I * (1 + S * cos(H) / cos(pi / 3 - H));
                B = 3 * I - R - G;
            }
            else
            {
                H -= 4 * pi / 3;
                G = I * (1 - S);
                B = I * (1 + S * cos(H) / cos(pi / 3 - H));
                R = 3 * I - G - B;
            }
            model[i][3 * j] = R;
            model[i][3 * j + 1] = G;
            model[i][3 * j + 2] = B;
            Rmin = min(Rmin, R);
            Gmin = min(Gmin, G);
            Bmin = min(Bmin, B);
        }

    //对r、g、b进行归一化
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            model[i][3 * j] -= Rmin;
            model[i][3 * j + 1] -= Gmin;
            model[i][3 * j + 2] -= Bmin;
            Rmax = max(Rmax, model[i][3 * j]);
            Gmax = max(Gmax, model[i][3 * j + 1]);
            Bmax = max(Bmax, model[i][3 * j + 2]);
        }
    for (int i = 0;i < height;i++)
        for (int j = 0;j < width;j++)
        {
            if (Rmax) model[i][3 * j] /= Rmax;
            if (Gmax) model[i][3 * j + 1] /= Gmax;
            if (Bmax) model[i][3 * j + 2] /= Bmax;
            *(origin + i * lineByte + 3 * j) = (unsigned char)(model[i][3 * j + 2] * 255 + 0.5);
            *(origin + i * lineByte + 3 * j + 1) = (unsigned char)(model[i][3 * j + 1] * 255 + 0.5);
            *(origin + i * lineByte + 3 * j + 2) = (unsigned char)(model[i][3 * j] * 255 + 0.5);
        }
    fwrite(origin, lineByte * height, 1, fp1);
    fclose(fp1);
    return 0;
}

代码说明

model 数组中当存的是 RGB 值时按照 R、G、B 的顺序存储,存的是 HSI 值时按照 H、S、I 的顺序存储,注意 BMP 文件的 RGB 值是按照 B、G、R 的顺序存储,最后存到 origin 数组中需调换顺序。
eps 是一个极小的量,不会影响最终结果,但是它是必要的,它能防止在计算 θ \theta θ 时出现 0 / 0 0/0 0/0 的状况,而且通常 0 / 0 0/0 0/0 不会报错,但输出图像时会出现很多无法察觉的噪声。
这是原图像:
在这里插入图片描述
这是不加 eps 时的图像(不同编译器不同,有的编译器输出结果变化特别大,此为 g++,MSVC比较逆天):
在这里插入图片描述
放大到阿尼亚的额头处:
在这里插入图片描述
加了 eps 后恢复正常:
在这里插入图片描述
注意归一化过程中做除法前判断分母是否为 0 0 0 ,为 0 0 0 的情况有原图所有像素 RGB 均相等(即实际上为灰度图时)、所有像素亮度值相等和其他的一些可能原因,总之需注意。

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
要在Python中进行彩色图像直方图均衡化,可以采用以下步骤: 1. 将彩色图像转换为HSI(Hue, Saturation, Intensity)模型。可以使用OpenCV库中的cvtColor()函数将RGB图像转换为HSI图像HSI模型中的亮度通道(Intensity)对应图像的灰度信息。 2. 对亮度通道进行直方图均衡化。可以使用OpenCV库中的equalizeHist()函数对亮度通道进行直方图均衡化操作。 3. 将处理后的亮度通道与原始的色调(Hue)和饱和度(Saturation)通道组合,得到均衡化后的HSI图像。 4. 将均衡化后的HSI图像转换回RGB模型。可以使用OpenCV库中的cvtColor()函数将HSI图像转换为RGB图像。 下面是一个简单的示例代码: ```python import cv2 # 读取彩色图像 image = cv2.imread('image.jpg') # 将图像转换为HSI模型 hsi_image = cv2.cvtColor(image, cv2.COLOR_BGR2HSV) # 提取亮度通道 intensity_channel = hsi_image[:,:,2] # 对亮度通道进行直方图均衡化 equalized_intensity = cv2.equalizeHist(intensity_channel) # 将均衡化后的亮度通道与原始的色调和饱和度通道组合 hsi_image[:,:,2] = equalized_intensity # 将均衡化后的HSI图像转换回RGB模型 equalized_image = cv2.cvtColor(hsi_image, cv2.COLOR_HSV2BGR) # 显示结果 cv2.imshow('Original Image', image) cv2.imshow('Equalized Image', equalized_image) cv2.waitKey(0) cv2.destroyAllWindows() ``` 在上面的示例代码中,我们首先读取了彩色图像,然后将其转换为HSI模型。接下来,我们提取了亮度通道,对其进行直方图均衡化,然后将处理后的亮度通道与原始的色调和饱和度通道组合,最后将均衡化后的HSI图像转换回RGB模型。最后,我们显示了原始图像均衡化后的图像。 请注意,上述代码仅为示例,你可以根据自己的需求进行适当的修改和优化。<span class="em">1</span><span class="em">2</span><span class="em">3</span> #### 引用[.reference_title] - *1* [(Python数字图像处理)彩色图像处理---色调和彩色校正以及直方图均衡化](https://blog.csdn.net/qq_44926189/article/details/121178739)[target="_blank" data-report-click={"spm":"1018.2226.3001.9630","extra":{"utm_source":"vip_chatgpt_common_search_pc_result","utm_medium":"distribute.pc_search_result.none-task-cask-2~all~insert_cask~default-1-null.142^v93^chatsearchT3_2"}}] [.reference_item style="max-width: 50%"] - *2* *3* [灰度图像和彩色图像直方图均衡化(python实现)](https://blog.csdn.net/zhao9428686/article/details/105394328)[target="_blank" data-report-click={"spm":"1018.2226.3001.9630","extra":{"utm_source":"vip_chatgpt_common_search_pc_result","utm_medium":"distribute.pc_search_result.none-task-cask-2~all~insert_cask~default-1-null.142^v93^chatsearchT3_2"}}] [.reference_item style="max-width: 50%"] [ .reference_list ]

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值