所有简易、非,       在对互相网格做统计以前

一、引言

自我初次接触HDR方面的学识,有描述不得法的地点烦请见谅。

为便宜小说讲述,引用部分百度中的文章对HDR图像举行简要的描述。

高动态范围图像(High-Dynamic
Range,简称HDR),比较平常的图像,可以提供更加多的动态范围和图像细节,依照分裂的揭露时间的LDR(Low-Dynamic
Range)图像,利用每个暴光时间相对应最佳细节的LDR图像来合成最后HDR图像,可以更好的反映人真实环境中的视觉效果。

具体确实存在的亮度差,即最亮的实体亮度,和微小的物体亮度之比为108
而人类的双眼所能看到的范围是105左右,可是一般的屏幕,照相机能表示的唯有256种分歧的亮度,计算一般的显示屏,照相机能表示的唯有256种不相同的亮度机在象征图象的时候是用8bit(256)级或16bit(65536)级来区分图象的亮度的,但那区区几百或几万不可能再次出现真实自然的黄石意况。HDR文件是一种越发图形文件格式,它的每一个像素除了常备的RGB音讯,还有该点的实际亮度音信。普通的图形文件每个象素唯有0
-255的灰度范围,那实则是不够的。想象一下太阳的发光强度和一个纯黑的物体之间的灰度范围或者说亮度范围的不一致,远远当先了256个级别。因而,一张普通的白昼风景图片,看上去白云和阳光可能都呈现是一致的灰度/亮度,都是纯白色,但骨子里白云和日光之间实际的亮度无法同样,他们中间的亮度差距是伟人的。因而,普通的图形文件格式是很不可依赖的,远远没有记录到具体世界的其实景况。而HDR格式则记录了很广范围内的亮度音讯。

不过最终,HDR图像要在屏幕中呈现,照旧须求对其数额进行处理的,如何处理即能丰裕利用这个数量,又能使得图像的显得尽量不丢失细节,是多年来不少图像工小编探讨的根本。

简短的说,就是当今有一堆离散的数量,数据的遍布范围或者很广,怎么样把那些离散的数码隐射到0到255里边。

二、相关算法的兑现

最简便的当然是线性隐射,先算出离散数据的最大值和纤维值,然后将数据线性的拉升至0到255里头,这种直白的操作往往不知所措取得满足的作用,会促成大量细节丢失,表现在视觉上就是一大块粉蓝色或者一大块白色的,如下图所示:

 图片 1  
   图片 2

     上边两幅图要么是暗部太暗,要么是亮部太亮,全体相比度太强,导致细节新闻多量丢失。

     针对这一难题,很四个人提议了无数卓殊不错的化解方案,比如依据全局操作符的,其中正文小编已毕其中的按照急忙双边滤波技术的HDR显示进度。

    
本文对应的参照随想地址: Fast
Bilateral Filtering for the Display of High-Dynamic-Range
Images

    
诗歌的细路很粗略,首先她将原本的HDR数据分解成四个层:base layer 和
detail layer,然后下落base layer的相比度,不转移detail
layer的数量,在将那两层合并。

    
其中:base layer的多寡用 HDR原始数据开展两岸滤波获取。

    
算法的简短流程入下所示:

    1、input
intensity= 1/61*(R*20+G*40+B)

    2、r=R/(input
intensity), g=G/input intensity, B=B/input intensity

    3、log(base)=Bilateral(log(input
intensity))

    4、log(detail)=log(input
intensity)-log(base)

    5、log
(output intensity)=log(base)*compressionfactor+log(detail) –
log_absolute_scale

    6、R
output = r*10^(log(output intensity)), etc.

    
上述进度中的变量compressionfactor,log_absolute_scale原文小编的提议取值为:

     
compressionfactor = targetContrast/(max(log(base)) – min(log(base)))
对于众多图像,targetContrast使用log(5)能获取较为理想的值。

     
而log_absolute_scale= max(log(base))*compressionfactor;

   
在进展双边的时候,SigmaS一般取值为0.02*马克斯(Width,Height)相比适当,而SigmaR取值0.4较为理想(这里是指多少量化到了0-1里面的)。

  所以都取优化的参数,则上述进度可以自动举行。

   
作者提到上述log操作都是以10为底开展的,我以为以e为底实际效果也没啥区其他。

  三、效果

    根据那一个思路编制程序后,确实能得到很不利的效应,比如上述两幅图像,依据前面讲的参数取值,解码后得到的图像如下:

   
图片 3 
图片 4

   
可知,图像的细节较为周密的反映出来了。

    当然,自动的参数不肯定能疏通最好的作用,比如依然那两幅图,手工选拔部分参数,可以调出如下效果:

   
图片 5 
图片 6

  更加是率先幅图,很有种蒙太奇的感觉。

   
在看望几张长出新在随笔中的图像的结果:

   
图片 7 
图片 8       

   
图片 9 
图片 10

   图片 11 图片 12

 
 图片 13 图片 14

   
          线性解码图                                             
双边滤波解码图

 
有些图线性解码啥都看不到,双边滤波解码后细节表现的就很显著了。

   
HDR格式的原本数据的解码可以凭借FreeImage来已毕,FreeImage就好像已经讲那几个多少量化到了0和1里头(不肯定不利)。一段不难的贯彻代码如下:

public Bitmap LoadHdrFormFreeImage(string FileName)
{
    Bitmap Bmp = null;
    FREE_IMAGE_FORMAT fif = FREE_IMAGE_FORMAT.FIF_UNKNOWN; ;
    if (FreeImage.IsAvailable() == true)
    {
        fif = FreeImage.GetFileType(FileName, 0);
        if (fif != FREE_IMAGE_FORMAT.FIF_HDR)
        {
            MessageBox.Show("不是Hdr格式的图像.");
            return null;
        }
        fif = FreeImage.GetFIFFromFilename(FileName);
        FIBITMAP Dib = FreeImage.Load(fif, FileName, FREE_IMAGE_LOAD_FLAGS.DEFAULT);
        uint Bpp = FreeImage.GetBPP(Dib);

        if (Bpp != 96)
        {
            MessageBox.Show("无法支持的Hdr格式.");
            FreeImage.Unload(Dib);
            return null;
        }
        uint Width = FreeImage.GetWidth(Dib);                        //  图像宽度
        uint Height = FreeImage.GetHeight(Dib);                      //  图像高度
        uint Stride = FreeImage.GetPitch(Dib);                       //  图像扫描行的大小,必然是4的整数倍
        IntPtr Bits = FreeImage.GetBits(Dib);

        float* Data = (float*)Bits;
        int Speed, Index;
        byte* Pixel;
        float Value;

        if (RawData != null) Marshal.FreeHGlobal((IntPtr)RawData);
        RawData = (float*)Marshal.AllocHGlobal((int)Width * (int)Height * 3 * sizeof(float));
        CopyMemory(RawData, Data, (int)Width * (int)Height * 3 * sizeof(float));

        Bmp = new Bitmap((int)Width, (int)Height, PixelFormat.Format24bppRgb);
        BitmapData BmpData = Bmp.LockBits(new Rectangle(0, 0, Bmp.Width, Bmp.Height), ImageLockMode.ReadWrite, PixelFormat.Format24bppRgb);
        Pixel = (byte*)BmpData.Scan0;

        for (int Y = 0; Y < Height; Y++)
        {
            Speed = Y * BmpData.Stride;
            Index = Y * (int)Width * 3;
            for (int X = 0; X < Width; X++)
            {
                Value = (Data[Index + 2] * 255);
                if (Value > 255)
                    Value = 255;
                else if (Value < 0)
                    Value = 0;
                Pixel[Speed] = (byte)Value;
                Value = (Data[Index + 1] * 255);
                if (Value > 255)
                    Value = 255;
                else if (Value < 0)
                    Value = 0;
                Pixel[Speed + 1] = (byte)Value;
                Value = (Data[Index + 0] * 255);
                if (Value > 255)
                    Value = 255;
                else if (Value < 0)
                    Value = 0;
                Pixel[Speed + 2] = (byte)Value;
                Index += 3;
                Speed += 3;
            }
        }
        FreeImage.Unload(Dib);
        Bmp.UnlockBits(BmpData);
        Bmp.RotateFlip(RotateFlipType.RotateNoneFlipY);
        return Bmp;
    }
    else
        return null;
}

  以上选拔的是线性解码。

  附一个解码的调用程序:http://files.cnblogs.com/Imageshop/ReadHdrTest.rar 

     更加多的源码可参照:http://people.csail.mit.edu/sparis/code/src/tone_mapping

                            http://people.csail.mit.edu/fredo/PUBLI/Siggraph2002/

    
一些大面积的用于测试的HDR图像可以从那边下载:http://www.pauldebevec.com/Research/HDR/

     同样,那种 tone
mapping算法也足以用在平常的RGB图像上,效果如下所示:

 
 图片 15  图片 16

              原图                                    增添base
layer的比较度

  
图片 17 
图片 18

             原图                                    下降base
layer的相比度

  
对于常见的总体偏暗的图像,该办法也能获得卓殊不错的效应,一些测试结果如下:

  
图片 19  
图片 20

   图片 21  
图片 22

   图片 23  
图片 24

   
很多任何的算法也能起到接近上述的功效,但是他俩一般很简单暴发halo现象。  

   
相信对于偏亮的常见照片,也足以有平等的处理能力(未找到合适的测试图片)。

  

 图片 25

*********************************作者:
laviewpbt   时间: 2013.11.18   联系QQ:  33184777
 转发请保留本行音讯************************

   

2.4.伪代码

       
下图为两岸近似的伪代码。1.早先化所有的下采样齐次值,令所有下采样权值为0。2.找到全图的细小亮度点(这一步是为了持续将亮度域的限定平移到从0开头)。3.对此高分辨率图像的任一个像素点(X,Y)和亮度I(X,Y),(a)同样总结其齐次格局;(b)统计下采样后的坐标(那里运用了取整,因而,在原式高分辨率图像中的相邻空间,会被映射到低分辨率下的同一空间坐标,但在亮度域可能不在同一坐标上);(c)更新下采样空间的坐标值,即在低分辨率空间上添加高分辨率的齐次值(之所以用齐次值,是因为齐次坐标相加后,取出原始数据即成功了对相应数据的加权平均:(w1i1,w1)+(w2i2,w2)=(w1i1+w2i2,w1+w2),数据对应(w1i1+w2i2)/(w1+w2))。4.在采样域做高斯卷积,获得卷积后的齐次坐标。5.上采样回高分辨率图像,对于未知的数据点,使用三线性差值。

图片 26

互相近似伪代码

 3 双边滤波源码完结

 1 #include <time.h>
 2 #include <stdio.h>
 3 #include <stdlib.h>
 4 #include <math.h>
 5 #include "bf.h"
 6 
 7 
 8 //--------------------------------------------------
 9 /**
10 双边滤波器图像去噪方法
11 \param   pDest       去噪处理后图像buffer指针
12 \param   pSrc        原始图像buffer指针
13 \paran   nImgWid     图像宽
14 \param   nImgHei     图像高
15 \param   nSearchR    搜索区域半径
16 \param   nSigma      噪声方差 
17 \param   nCoff       DCT变换系数个数
18 
19 return   void  
20 */
21 //--------------------------------------------------
22 void BilateralFilterRGB(float *pDest, BYTE *pSrc,int nImgWid, int nImgHei, int nSearchR, int nSigma, int nCoff)
23 {
24 
25     BYTE *pSrcR  = new BYTE[nImgHei * nImgWid];                  // 新建图像缓冲区R
26     BYTE *pSrcG  = new BYTE[nImgHei * nImgWid];                  // 新建图像缓冲区G
27     BYTE *pSrcB  = new BYTE[nImgHei * nImgWid];                  // 新建图像缓冲区B
28 
29 
30     for(int i = 0;i < nImgHei; i++)
31     {
32         for (int j = 0; j < nImgWid; j++)
33         {
34             pSrcR[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 0];         
35             pSrcG[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 1];         
36             pSrcB[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 2];         
37         }
38     }    
39 
40 
41     float* pDestR = new float[nImgHei * nImgWid];         // 去噪后的图像数据R通道
42     float* pDestG = new float[nImgHei * nImgWid];         // 去噪后的图像数据G通道
43     float* pDestB = new float[nImgHei * nImgWid];         // 去噪后的图像数据B通道
44     memset(pDestR,255,nImgHei * nImgWid * sizeof(float)); // 传递变量初始化
45     memset(pDestG,255,nImgHei * nImgWid * sizeof(float));
46     memset(pDestB,255,nImgHei * nImgWid * sizeof(float));
47 
48 
49     float *pII = new float[nImgHei * nImgWid ];                       // 积分图
50     float *pWeights  = new float[nImgHei * nImgWid ];                       // 每一个窗的权重之和,用于归一化  
51 
52 
53     float mDctc[MAX_RANGE];                                     // DCT变换系数
54     float mGker[MAX_RANGE];                                     // 高斯核函数的值
55 
56     for (int i = 0; i< MAX_RANGE; ++i)
57     {
58         mGker[i] = exp(-0.5 * pow(i / nSigma, 2.0));                       // 首先计算0-255的灰度值在高斯变换下的取值,然后存在mGker数组中
59     }
60 
61     CvMat G = cvMat(MIN_RANGE, MAX_RANGE, CV_32F, mGker);           // 转换成opencv中的向量
62     CvMat D = cvMat(MIN_RANGE, MAX_RANGE, CV_32F, mDctc);   // 把mDctc系数数组转换成opencv中的向量形式
63 
64     cvDCT(&G,&D,CV_DXT_ROWS);                                           // 然后将高斯核进行离散的dct变换
65     mDctc[0] /= sqrt(2.0);                                              // 离散余弦变换的第一个系数是1/1.414;
66 
67 
68     bf(pSrcR, pDestR, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR);
69     bf(pSrcG, pDestG, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR);
70     bf(pSrcB, pDestB, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR);
71 
72 
73     for(int i=0; i < nImgHei;i++)
74     {
75         for (int j=0; j < nImgWid; j++)
76         {
77             pDest[i * nImgWid * 3 + j * 3 + 0] = pDestR[i * nImgWid + j];           
78             pDest[i * nImgWid * 3 + j * 3 + 1] = pDestG[i * nImgWid + j];                   
79             pDest[i * nImgWid * 3 + j * 3 + 2] = pDestB[i * nImgWid + j];           
80         }
81     }
82 
83 
84     free(pSrcR);      // 释放临时缓冲区
85     free(pSrcG);
86     free(pSrcB);
87     free(pDestR);
88     free(pDestG);
89     free(pDestB);
90 
91 }

  1 //--------------------------------------------------
  2 /**
  3 双边滤波器图像去噪方法
  4 \param   pDest       去噪处理后图像buffer指针
  5 \param   pSrc        原始图像buffer指针
  6 \paran   nImgWid     图像宽
  7 \param   nImgHei     图像高
  8 \param   pII         积分图缓冲区
  9 \param   pWeights    归一化权重系数 
 10 \param   nCoff       DCT变换系数
 11 \param   nPatchWid   图像块宽度
 12 
 13 return   void  
 14 */
 15 //--------------------------------------------------
 16 void bf (BYTE *pSrc, float *pDest, float *pDctc, float *pII, float *pWeights, int nImgHei, int nImgWid, int nCoff, int nPatchWid)  
 17 {  
 18     if (pDest == NULL || pSrc ==NULL)
 19     {
 20         return;
 21     }
 22 
 23     float *cR  = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float));
 24     float *sR  = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float));
 25     float *dcR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float));
 26     float *dsR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float));
 27 
 28     
 29     //  余弦变换查找表
 30     int fx     = 0;
 31     int ck     = 0;
 32     int r2     = pow(2.0*nPatchWid + 1.0, 2.0);         // 这个是图像块的像素点的总个数
 33 
 34     float c0   = pDctc[0];                              // 这是DCT的第一个系数
 35     float c0r2 = pDctc[0] * r2;                         // 这个用于初始化归一化权重之和
 36 
 37 
 38     for (ck = 1; ck < nCoff; ck++)                            // 注意到这里面的系数并不包含C0,所以后面才有把C0额外的加上
 39     {                                                                     // 一个矩阵,用于存放各种系数和数值,便于查找
 40         int ckr = (ck - 1) * MAX_RANGE;                                 // 数组是从0开始的,所以这里减1
 41 
 42         for (fx = 0; fx<MAX_RANGE; fx++)                              // fx就是图像的灰度值
 43         {
 44             float tmpfx = PI * float(fx) * float(ck) / MAX_RANGE;   // ck其实就相当于那个余弦函数里面的t,fx相当于u,PI/MAX相当于前面的那个系数,这都是余弦变换相关的东西
 45             
 46             cR[ckr + fx]  = cos(tmpfx);                                      // 存储余弦变换,这个用在空间域上      
 47             sR[ckr + fx]  = sin(tmpfx);                                      // 存储正弦变换
 48             dcR[ckr + fx] = pDctc[ck]  * cos(tmpfx);                         // 存储余弦变换和系数的乘积,这个则用在强度范围上
 49             dsR[ckr + fx] = pDctc[ck]  * sin(tmpfx);                         // 存储正弦变换和系数的乘积       
 50         }        
 51     }
 52 
 53       
 54     float *pw  = pWeights;                          // 新建一个归一化权重的中间变量进行数据的初始化
 55     float *pwe = pWeights + nImgWid * nImgHei;      // 限定最大位置
 56 
 57     while (pw < pwe) 
 58     {
 59         *pw++ = c0r2;                   // 赋初值,让它们都等于第一个DCT系数乘以图像块中总的像素点的个数
 60     }
 61     
 62     for (int ck = 1; ck < nCoff; ck++) 
 63     {        
 64         int ckr = (ck-1)*MAX_RANGE; // 数组是从0开始的,所以这里减1
 65 
 66         add_lut(pWeights, pSrc, pII,cR + ckr, dcR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid);  // add cos term to pWeights
 67         add_lut(pWeights, pSrc, pII,sR + ckr, dsR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid);  // add cos term to pWeights  
 68         
 69     } 
 70 
 71     ImgRectSumO(pDest, pSrc, pII, c0, nImgHei, nImgWid, nPatchWid, nPatchWid);  //加上C0的变换之后,再初始化滤波后的函数     
 72     
 73     for (int ck = 1; ck < nCoff; ck++) 
 74     {        
 75         int ckr = (ck-1)*MAX_RANGE;
 76         
 77         add_f_lut(pDest, pSrc, pII, cR + ckr, dcR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid);  // add cosine term to dataf
 78         add_f_lut(pDest, pSrc, pII, sR + ckr, dsR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid);  // add sine term to dataf
 79   
 80     }
 81 
 82 
 83     float *pd = pDest + nPatchWid * nImgWid;
 84     pw        = pWeights + nPatchWid * nImgWid;
 85 
 86     float *pdie  = pd + nImgWid - nPatchWid;
 87     float *pdend = pDest + (nImgHei - nPatchWid) * nImgWid;
 88     
 89     while (pd < pdend) 
 90     {
 91         pd += nPatchWid; //把边都给去掉了
 92         pw += nPatchWid; //这个也是把边去掉了
 93 
 94         while (pd < pdie) 
 95         {
 96             *pd++ /= ( (*pw++)*255);     // 之所以要除以255,是为了把它化到[0,1]之间,便于显示
 97         }
 98 
 99         pd += nPatchWid;                 // 把边都给去掉了
100         pw += nPatchWid;                 // 这个也是把边去掉了
101         pdie += nImgWid;                 // 过渡到下一行       
102     }
103    
104     free(cR); 
105     free(sR); 
106     free(dcR); 
107     free(dsR);   // 释放缓冲区
108 }

//--------------------------------------------------
/**
余弦函数首系数积分图
\param   pDest       去噪处理后图像buffer指针
\param   pSrc        原始图像buffer指针
\paran   nImgWid     图像宽
\param   nImgHei     图像高
\param   pII         积分图缓冲区
\param   c0          DCT变换第一个系数
\param   nPatchWid   图像块宽度
\param   nPatchHei   图像块高度

return   void  
*/
//--------------------------------------------------
void ImgRectSumO (float* pDest, const BYTE* pSrc, float* pII, float c0,const int nImgHei, const int nImgWid, const int nPatchWid, const int nPatchHei) 
{
    if (pDest == NULL || pSrc ==NULL)
    {
        return;
    }

    int ri1 = nPatchWid + 1;                                // 中心像素点
    int rj1 = nPatchHei + 1;
    int ri21 = 2 * nPatchWid + 1;                           // 图像块的宽度
    int rj21 = 2 * nPatchHei + 1;

    const BYTE *pSrcImg    = pSrc;                          // 原始图像数据
    const BYTE *pSrcImgEnd = pSrc + nImgHei * nImgWid;      // 最大的像素点位置    
    const BYTE *pSrcImgWid = pSrcImg + nImgWid;             // 第一行的最大的像素点位置,下面循环用的的变量
    float       *pIITmp     = pII;                           // 积分图数据
    float       *ppII       = pIITmp;                        // 积分图中间变量


    *pIITmp++ = c0 * float( *pSrcImg++ );                    // 第一个系数乘以图像数据

    while (pSrcImg < pSrcImgWid)                             // 遍历第一行进行求积分图,其实这个是图像数据的积分图
    {
        *pIITmp++ = (*ppII++) + c0 * float(*pSrcImg++);      
    }  

    while (pSrcImg < pSrcImgEnd)                             // 遍历所有的像素点的变量
    {    

        pSrcImgWid   += nImgWid;                             // 进行第二行的循环
        ppII          = pIITmp;                              // 积分图中间变量
        *pIITmp++     = c0 * float(*pSrcImg++);              // 第二行第一个像素点的处理

        while (pSrcImg < pSrcImgWid) 
        {
            (*pIITmp++) = (*ppII++) + c0 * float(*pSrcImg++); // 求第二行的积分图
        }        

        float *pIIWid = pIITmp;                        // 当前位置像素点
        pIITmp = pIITmp - nImgWid;                     // 上一行的起始点位置
        float *pIIup  = pIITmp - nImgWid;              // 上上一行的起始点位置

        while (pIITmp < pIIWid)                        // 遍历整个行的每一个数据
        {
            (*pIITmp++) = (*pIITmp) + (*pIIup++);      // 行与行之间的积分图
        }

    }


    float *pres = pDest + ri1 * nImgWid;                    // 最小的行位置
    float *pend = pDest + (nImgHei - nPatchWid) * nImgWid;  // 最大的行位置

    float *pii1 = NULL;
    float *pii2 = NULL;
    float *pii3 = NULL;
    float *pii4 = NULL;                               // 求积分图的四个点

    pii1 = pII + ri21 * nImgWid + rj21;               // 右下角
    pii2 = pII + ri21 * nImgWid;                      // 左下角
    pii3 = pII + rj21;                                // 右上角
    pii4 = pII;                                       // 左上角

    while (pres < pend) 
    {
        float *pe = pres + nImgWid - nPatchHei;
        pres += rj1;

        while (pres < pe)                             // 这个只是求图像数据的积分图
        {
            (*pres++) = (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++);
        }

        pres += nPatchHei;
        pii1 += rj21;
        pii2 += rj21;
        pii3 += rj21;
        pii4 += rj21;
    }

}

//--------------------------------------------------
/**
余弦积分图加函数
\param   pNomalWts   像素点的归一化权重矩阵
\param   pSrc        原始图像buffer指针
\param   pII         积分图
\paran   nImgWid     图像宽
\param   nImgHei     图像高
\param   pII         积分图缓冲区
\param   pCosMtx     余弦函数矩阵
\param   pCosDctcMtx 余弦函数与DCT变换系数乘积矩阵
\param   nPatchWid   图像块宽度
\param   nPatchHei   图像块高度

return   void  
*/
//--------------------------------------------------

void add_lut (float* pNomalWts, const BYTE* pSrc, float* pII, float* pCosMtx, float* pCosDctcMtx, const int nImgHei, const int nImgWid, const int nPatchWid, const int nPatchHei) 
{ 

    int ri1  = nPatchWid + 1;      // kernel相关
    int rj1  = nPatchHei + 1;
    int ri21 = 2 * nPatchWid + 1;  // 这是kernel的宽度和长度
    int rj21 = 2 * nPatchHei + 1;

    const BYTE *pi  = pSrc;           // 这个是图像数据
    float *pii       = pII;            // 这个是中间的缓冲变量,积分图的缓冲区
    const BYTE *piw = pi + nImgWid;   // 也是赋初值的时候使用的中间变量
    float *pii_p     = pii;            // 直线积分图指针的指针


    *pii++ = pCosMtx[*pi++];  //图像数据值所对应的查找表中的余弦或者正弦变换


    while (pi < piw)
    {            
        *pii++ = (*pii_p++) + pCosMtx[*pi++]; // 第一行的积分图
    }    

    const BYTE *piend = pSrc + nImgHei * nImgWid;       // 限定循环位置

    while (pi < piend)                                     // 遍历整个图像数据
    {  

        piw    += nImgWid;                            // 第二行
        pii_p   = pii;                             // 指向积分图指针的指针
        *pii++  = pCosMtx[*pi++];              // 获取第一个图像数据值所对应的查找表中的余弦或者正弦变换

        while (pi < piw) 
        {
            (*pii++) = (*pii_p++) + pCosMtx[*pi++]; // 求这一行的积分图
        }

        float *piiw = pii;
        pii = pii - nImgWid;
        float *pii_p1 = pii - nImgWid;

        while (pii < piiw) 
        {
            (*pii++) = (*pii) + (*pii_p1++);                  // 行与行之间求积分图
        }

    }  


    float *pres = pNomalWts + ri1 * nImgWid;                   // 定位要处理点的像素的那一行
    float *pend = pNomalWts + (nImgHei - nPatchWid) * nImgWid; // 限定了边界

    float *pii1 = NULL;
    float *pii2 = NULL;
    float *pii3 = NULL;
    float *pii4 = NULL;

    pii1 = pII + ri21 * nImgWid + rj21;  // 获得积分图中那个最大的数据右下角
    pii2 = pII + ri21 * nImgWid;         // 积分图左下角
    pii3 = pII + rj21;                   // 积分图右上角
    pii4 = pII;                          // 积分图左上角

    pi = pSrc + ri1 * nImgWid;                  // 定位要处理的像素的位置


    while (pres < pend)                         // 设定高度做循环
    {
        float *pe = pres + nImgWid - nPatchHei; // 限定了宽度的范围
        pres = pres + rj1;                      // 定位了要处理的像素点的归一化系数矩阵的位置
        pi = pi + rj1;                          // 定位了要处理的像素点

        while (pres < pe)                       // 设定宽度做循环
        {
            (*pres++) = (*pres) + pCosDctcMtx[*pi++] * ( (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++) );    // 得到的就是归一化系数矩阵之和
        }

        pres += nPatchHei;                     
        pi   += nPatchHei;                     
        pii1 += rj21;                          
        pii2 += rj21;                          
        pii3 += rj21;                          
        pii4 += rj21;                          // 略过不处理的边界区域
    }    

}

//--------------------------------------------------
/**
积分图
\param   pDest   像素点的归一化权重矩阵
\param   pSrc        原始图像buffer指针
\param   pII         积分图
\paran   nImgWid     图像宽
\param   nImgHei     图像高
\param   pII         积分图缓冲区
\param   pCosMtx     余弦函数矩阵
\param   pCosDctcMtx 余弦函数与DCT变换系数乘积矩阵
\param   nPatchWid   图像块宽度
\param   nPatchHei   图像块高度

return   void  
*/
//--------------------------------------------------

void add_f_lut (float* pDest,const BYTE* pSrc,float* pII,float* pCosMtx,float* pCosDctcMtx,const int nImgHei, const int nImgWid,const int nPatchWid, const int nPatchHei) 
{  

    int ri1 = nPatchWid + 1;      // kernel相关,目标是指向中心像素点
    int rj1 = nPatchHei + 1;      // 目标是指向中心像素点
    int ri21 = 2 * nPatchWid + 1; // 这是kernel的宽度和长度
    int rj21 = 2 * nPatchHei + 1; // 其实就是指向搜索区域的右下角,也就是最大位置的那个像素点,就是边界了

    const BYTE *pi = pSrc;       // 这个是图像数据,图像的灰度值
    float *pii = pII;             // 这是积分图

    const BYTE *piw = pi + nImgWid;  // 指向第一行的末尾位置,用于循环控制
    float *pii_p = pii;               // 指向积分图指针的指针

    *pii++ = float(*pi) * pCosMtx[*pi]; // 灰度值乘以余弦系数
    ++pi;      

    while (pi < piw)                    // 第一行遍历
    {
        *pii++ = (*pii_p++) + float(*pi) * pCosMtx[*pi];
        ++pi;
    }

    const BYTE *piend = pSrc + nImgHei * nImgWid;

    while (pi < piend) 
    {        
        piw   += nImgWid;
        pii_p  = pii;
        *pii++ = float(*pi) * pCosMtx[*pi];
        ++pi;

        while (pi < piw)
        {
            (*pii++) = (*pii_p++) + float(*pi) * pCosMtx[*pi];
            ++pi;
        }

        float *piiw = pii;
        pii = pii - nImgWid;   // 上一行起点
        float *pii_p1 = pii - nImgWid; // 上上一行的起始点

        while (pii < piiw)             // 循环一行
        {
            (*pii++) += (*pii_p1++);  //其实就是在列的位置上加上上一行
        }

    }


    float *pres = pDest + ri1*nImgWid;                     
    float *pend = pDest + (nImgHei - nPatchWid) * nImgWid; 
    float *pii1 = NULL;
    float *pii2 = NULL;
    float *pii3 = NULL;
    float *pii4 = NULL;

    pii1 = pII + ri21 * nImgWid + rj21;   // 积分图搜索区域最大位置的元素,堪称右下角的元素
    pii2 = pII + ri21 * nImgWid;          // 可以看成搜索区域中左下角的像素位置
    pii3 = pII + rj21;                    // 搜索区域右上角像素位置
    pii4 = pII;                           // 搜索区域左上角像素位置
    pi   = pSrc + ri1 * nImgWid;          // 定位要处理的像素的位置的那一行

    while (pres < pend) 
    {
        float *pe = pres + nImgWid - nPatchHei; //限定了宽度的范围

        pres = pres + rj1; //定位了要处理的像素点的归一化系数矩阵的位置
        pi   = pi   + rj1;     //定位了要处理的像素点

        while (pres < pe)  //遍历整个行
        {
            (*pres++) = (*pres) + pCosDctcMtx[*pi++] * ( (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++) ); //这个其实是计算整个像素块
        }

        pres += nPatchHei;   
        pi   += nPatchHei;   

        pii1 += rj21; 
        pii2 += rj21;
        pii3 += rj21;
        pii4 += rj21;
    }
}

 

2.互相滤波器的急忙近似

       在对相互网格做总括从前,先介绍一下多头滤波器的全速近似方法(A
Fast Approximation of the Bilateral Filter using a Signal Processing
Approach
),它是两者网格的雏形,是一个将双方滤波器扩大到高维空间举行线性卷积的主意。而两端网格的撰稿人在那篇文章的根底上,将高维空间数据映射成3D
array,进而提出了双方网格。

3 双边滤波落成步骤

  有了上述辩解之后完成Bilateral Filter就相比较不难了,其实它也与普通的Gaussian Blur没有太大的不相同。那里根本概括3有的的操作: 基于空间距离的权重因子变化;基于相似度的权重因子的成形;最后filter颜色的臆想。

2.1.相互滤波器

       
首先介绍一下怎么着是双方滤波器,那里引用一些客人的下结论博客:三头滤波算法原理

       
简单地说,双边滤波器是一种高斯滤波器的壮大。传统的高斯滤波器有一个高斯核,通过对空中相邻的像素点以高斯函数位权值取均值,落成对图像的平缓。由于观念高斯滤波器只考虑了空荡荡的新闻,即使它能兑现对图像的得力平滑,但还要也搅乱了边缘消息。双边滤波器的规律即在价值观高斯滤波器的根基上添加了一个表征亮度差距的高斯核,即既考虑了上空的信息,由考虑了值域的音信。在灰度差异不大的限量,表征亮度高斯核的权值较大、接近于1,因而双方滤波退化为观念的高斯滤波器,完结对图像的平整;在边缘部分,即便相邻像素点的半空中接近,空间距离小,但由于边缘部分,灰度差距较大,由此在值域上相差较大,所以出席的新高斯核使得在该部分不实施平滑处理,保留图像的边缘。所以说,双边滤波是一种非线性(四个高斯核的乘积)的滤波方法,是整合图像的上空邻近度和像素值相似度的一种折中拍卖,同时考虑空域音信和灰度相似性,达到保边去噪的目的。

2 双边滤波原理

  滤波算法中,指标点上的像素值经常是由其所在位置上的四周的一个小一些邻居像素的值所决定。在2D高斯滤波中的具体完结就是对周围的早晚范围内的像素值分别赋以分歧的高斯权重值,并在加权平均后获取当前点的末段结出。而那里的高斯权重因子是利用多少个像素之间的半空中中远距离(在图像中为2D)关系来变化。通过高斯分布的曲线可以窥见,离目的像素越近的点对终极结果的孝敬越大,反之则越小。其公式化的叙说相似如下所述:

                     图片 27

   其中的c即为基于空间距离的高斯权重,而用来对结果举办单位化。

  高斯滤波在低通滤波算法中有不错的显现,可是其却有其它一个难题,那就是只考虑了像素间的长空地点上的关联,因而滤波的结果会丢掉边缘的音讯。那里的边缘紧若是指图像中根本的两样颜色区域(比如藏蓝色的天空,蓝色的毛发等),而Bilateral就是在Gaussian blur中进入了别的的一个权重分部来化解这一标题。Bilateral滤波中对此边缘的维持通过下述表明式来兑现:

                     图片 28

                    图片 29

  其中的s为基于像素间相似程度的高斯权重,同样用来对结果开展单位化。对两岸进行重组即可以赢得基于空间距离、相似程度综合考量的Bilateral滤波:

                  图片 30

  上式中的单位化分部综合了三种高斯权重于联合而得到,其中的c与s计算可以详细描述如下:

                  图片 31

   且有图片 32

                   图片 33

   且有图片 34

  上述给出的表明式均是在空间上的最好积分,而在像素化的图像中本来不可以这么做,而且也没须求这么做,因此在采取前必要对其举行离散化。而且也不须要对此每个局地像素从整张图像上开展加权操作,距离当先一定程度的像素实际上对现阶段的对象像素影响很小,可以忽略的。限定局地子区域后的离散化公就可以简化为如下情势:

                   图片 35

  上述申辩公式就组成了Bilateral滤波完成的底子。为了直观地明白高斯滤波与互相滤波的分别,我们得以从下列图示中看看按照。假如目标源图像为下述左右区域泾渭鲜明的带有噪声的图像(由程序自动生成),紫色框的主干即为目标像素所在的任务,那么当前像素处所对应的高斯权重与双方权重因子3D可视化后的造型如前边两图所示:

                   图片 36

  左图为本来的噪音图像;中间为高斯采样的权重;右图为Bilateral采样的权重。从图中可以看出Bilateral插足了一般程度分部将来可以将源图像左边这一个跟当前像素差值过大的点给滤去,那样就很好地涵养了边缘。为了越发形象地洞察两者间的界别,使用Matlab将该图在三种不一样方式下的惊人图3D绘制出来,如下:

                 图片 37

  上述三图从左到右依次为:双边滤波,原始图像,高斯滤波。从中度图中得以显明看出Bilateral和Gaussian三种艺术的分别,前者较好地保全了边缘处的梯度,而在高斯滤波中,由于其在边缘处的扭转是线性的,由此就选取连累的梯度突显出渐变的状态,而那彰显在图像中的话就是境界的散失(图像的示范可知于后述)。         

4.相互网格的上采样

        提出双方网格的撰稿人接着又拓展了上下一心的两端网格,公布了Bilateral
Guided
Upsampling
,介绍了何等运用双边网格达成部分图像操作算子。算法的要旨情想就是将一副高分辨率的图像通过下采样转换成一个两岸网格,而在两岸网格中,每个cell里提供一个图像变换算子,它的原理是在上空与值域相近的区域内,相似输入图像的亮度经算子变换后也应有是形似的,由此在每个cell里的操作算子可以视作是输入/输出的近乎曲线,也即一个仿射模型,而对此不一样的cell,小编通过给定输入和梦想输出去陶冶这一个两岸网格落成其仿射模型的全局分段平滑。再经过上采样,得到高分辨率的处理后的图像。

图片 38

Bilateral Guided Upsampling

3.4 Color Filtering

   其中的相似度分部的权重s主要基于五个Pixel之间的水彩差值总括面来。对于灰度图而言,这么些差值的界定是足以预言的,即[-255, 255],由此为了狠抓统计的频率大家可以将该有的权重因子预总计生成并存表,在行使时神速查询即可。使用上述完结的算法对几张带有噪声的图像进行滤波后的结果如下所示:

2.2.连忙近似

图片 39

三头网络接近例子

       
作者首先将双边滤波器的计算公式写成了齐次坐标的花样(关于齐次坐标的明亮,可以参照至于齐次坐标的驾驭(经典)):

图片 40

双面滤波齐次表示

       
接着,作者引入δ函数,将公式从二维空间扩大到三维空间(即一个2D的空间域和一个1D的值域/亮度域):

图片 41

引入δ函数

        小编引入一些函数来代表上式中的一些操作:

图片 42

高斯核

图片 43

δ函数

        最后,将双方滤波器的公式改写为:

图片 44

相互滤波的线性格局

图片 45

改正后的两岸滤波进程

       
小编提出,由于采样定理,能够透过对图像下采样做卷积处理,再上采样复苏原式分辨率,达到急速的目标。

图片 46

一对标注的求证

1 双边滤波简介

  双边滤波(Bilateral filter)是一种非线性的滤波方法,是整合图像的空中邻近度和像素值相似度的一种折衷处理,同时考虑空域音讯和灰度相似性,达到保边去噪的目标。具有简易、非迭代、局地的性状。

  双边滤波器的好处是足以做边缘保存(edge preserving),一般过去用的维纳滤波或者高斯滤波去降噪,都会较肯定地混淆边缘,对于频仍细节的掩护作用并不明确。双边滤波器顾名思义比高斯滤波多了一个高斯方差sigma-d,它是基于空间分布的高斯滤波函数,所以在边缘附近,离的较远的像素不会太多影响到边缘上的像素值,那样就保证了边缘附近像素值的保留。可是出于保存了过多的频仍音信,对于彩色图像里的反复噪声,双边滤波器无法彻底的滤掉,只好够对于低频新闻进行较好的滤波。

3.双方网格

       
有专家总计了二者滤波器的接近方法,使用一种3D数组的款型来表示那种在三维空间的双边滤波,提出了双方网格的辩解(Real-time
edge-aware image processing with the bilateral
grid
)。

图片 47

互相网格的双面滤波器已毕原理

       
有了前边双边滤波近似的论争,再来领会两者网格就显示简单得多了。考虑上图所示的1D输入图像。双边网格是在空间域和亮度域举办采样,划分成网格。双边的边也就是从那里来的,空间(space)和亮度(range)。离散后,每个点的坐标和亮度新闻取整到相应的网格内。每个值域内的亮度值可由加权平均获得。通过在网格内进行滤波等拍卖,再由上采样方法插值出处理后的原来图像。

       
由此,打造一个两端网格的主意是:首先,根据空域和值域,划分出网格结构,并开首化所有网格的节点为0:

图片 48

初叶化节点

       
对于高分辨率图像的坐标,分别除以空间采样间隔和亮度采样间隔并取整,再对应的节点处填入其原始的灰度值。考虑到取整函数会促成原图像的某一邻域内的像素点会下采样到同一个网格中,因而将灰度值写成齐次值(I(x,y),1)。并对在该网格内的具备齐次值叠加。那样,这几个网格对应的末段灰度最可以写成这个灰度值的平均值:(I(x1,y1)+…+I(xn,yn))/n。

图片 49

网格填充

       
那样,大家就将原先的图像数据转载为了一个3D的双边网格。任何图像操作函数都足以作用在这么些双方网格上,即一对一于将该双边网格左乘一个3D图像操作矩阵(例如,对于相互滤波器而言,这几个函数是一个高斯卷积核,包涵了半空中的方差和亮度的方差),得到一个在低分辨率情况下拍卖的结果。最终经过上采样,并开展插值,得到最终于高分辨率图像结果。上采样的原理是选项一个参考图,对其自由一个空中的像素进行空域和值域的采样,那里不进行取整操作,然后找到其在网格中的位置,由于尚未了取整,要求动用三线性插值的法子,来完结未知范围的亮度值的乘除,那些进度称作“slicing”。由插值上采样后,我们得到了最终的滤波结果。

3.1Spatial Weight

  那就是平时的Gaussian Blur中利用的测算高斯权重的方法,其首要通过四个pixel之间的相距并行使如下公式统计而来:

                 图片 50

1.前言

        目前在看Deep Bilateral Learning for Real-Time Image
Enhancement
,是一篇用CNN完毕实时图像增强的不二法门,可以在三哥大上飞快落成HDR。小说中涉嫌到相互网格(Bilateral
Grid)的盘算,查阅了大气的文献,对两端网格的驳斥做了弹指间私家的明亮。

3.1Similarity Weight

  与基于距离的高斯权重总结类似,只然而此处不再根据多个pixel之间的上空距离,而是基于其相似程度(或者多个pixel的值时期的偏离)。

            图片 51

   其中的代表多个像素值之间的相距,可以直接采纳其灰度值之间的差值或者RGB向量之间的欧氏距离。

相关文章