顺着列方向一行一行的更行整行的列直方图,中值滤波是一种经典的图像操作

     首要参考论文:Median Filter
in Constant
Time.pdf

   
 在图像处理中,局部算法一般的话,在很大程度上会得到比全局算法更为好的作用,因为她考虑到了图像领域像素的新闻,而不少片段算法可以依靠直方图获得加快。同时,一些常规的算法,比如中值滤波、最大值滤波、最小值滤波、表面模糊等等都得以由此一些直方图进行加快。而传统的获取局地直方图总结量很大,越发是半径增加时,耗时会成平方关系增添。一些有些算法唯有在半径较大时才会取得很好的作用,因而,必须找到一种适于的增速总结局地直方图的章程。

     
SSE图像算法优化种类五:超高速指数模糊算法的完毕和优化(10000*10000在100ms左右落到实处) 一文中,我早已说过优化后的ExpBlur比BoxBlur还要快,那多少个时候自己相比的BoxBlur算法是通过积分图+SSE达成的,我在09年此外一个博客账号上已经提供过一篇那一个文章彩色图像高速模糊之懒惰算法,里面也介绍了一种高效的图像模糊算法,那一个算法的施行时间基本也是和半径无关的。在当年的SSE优化学习之路上自家早就也设想过将该算法使用SSE达成,但立时认为这么些算法逐像素同时逐行都是上下看重的(单纯的逐像素信赖算法我在指数模糊里有关系怎么着用SSE优化),是无能为力用SSE处理的,一向没考虑,直到日前有心上人提出某个基于局地局方差的算法希望能提速,我又重新触发灵感,终于将那种算法也促成的下令集完结,并且测试速度比积分图快二倍,比解析opencv中Box
Filter的落成并提议更进一步加速的方案(源码共享)
那篇作品的速度快3倍,比opencv的cvSmooth函数快5倍,在一台老旧的I3台式机上处理3000*2000的灰度图达到了6ms的进程,本文分享该优化进程并提供灰度版本的优化代码供我们学习和议论。

   
参考代码:http://files.cnblogs.com/Imageshop/CTMF.rar

      在参考Median Filter in Constant
Time.pdf
一文附带的C的代码的基础上,本文提议了依据SSE加快的恒长任意半径局地直方图获得技能,可以大大加快算法的预计时间,更加是大半径时的涨潮更为明朗。

  在彩色图像高速模糊之懒惰算法一文附带的代码中(VB6的代码)是对准24位的图像,为了商讨方便,大家先贴出8位灰度的C++的代码:

   
中值滤波是一种经典的图像操作,越发适用于椒盐噪音的删除。同样,他也是USM锐化(表示疑虑,我记念是高斯滤波)、顺序处理、形态学操作(比如去孤点)等算法的基础。更高级其他利用包涵目的划分、语音和文字识别以及历史学图像处理等。

     
首要的优化思路是,沿着列方向一行一行的更行整行的列直方图,新的一条龙对应的列直方图更新时只要求裁减已经不再限制内的老大像素同时出席新进入的像素的直方图音信。之后,对于一行中的第三个像素点,累加半径辐射范围内的列直方图,获得改点的片段直方图,对于行中的其余的像素,则类似于革新行直方图,先减去不在范围内那列的列直方图,然后加上移入范围内的列直方图。由于选择了依照SSE函数的加快进度,直方图想加和相减的进度较一般的加减法有了10倍以上的涨潮,因而大大的升高了总体的实用性。

  1 /// <summary>
  2 /// 按照Tile方式进行数据的扩展,得到扩展后在原尺寸中的位置,以0为下标。
  3 /// </summary>
  4 int IM_GetMirrorPos(int Length, int Pos)
  5 {
  6     if (Pos < 0)
  7         return -Pos;
  8     else if (Pos >= Length)
  9         return Length + Length - Pos - 2;
 10     else    
 11         return Pos;
 12 }
 13 
 14 void FillLeftAndRight_Mirror_C(int *Array, int Length, int Radius)
 15 {
 16     for (int X = 0; X < Radius; X++)
 17     {
 18         Array[X] = Array[Radius + Radius - X];
 19         Array[Radius + Length + X] = Array[Radius + Length - X - 2];
 20     }
 21 }
 22 
 23 int SumofArray_C(int *Array, int Length)
 24 {
 25     int Sum = 0;
 26     for (int X = 0; X < Length; X++)
 27     {
 28         Sum += Array[X];
 29     }
 30     return Sum;
 31 }
 32 
 33 /// <summary>
 34 /// 将整数限幅到字节数据类型。
 35 /// </summary>
 36 inline unsigned char IM_ClampToByte(int Value)            //    现代PC还是这样直接写快些
 37 {
 38     if (Value < 0)
 39         return 0;
 40     else if (Value > 255)
 41         return 255;
 42     else
 43         return (unsigned char)Value;
 44     //return ((Value | ((signed int)(255 - Value) >> 31)) & ~((signed int)Value >> 31));
 45 }
 46 
 47 int IM_BoxBlur_C(unsigned char *Src, unsigned char *Dest, int Width, int Height, int Stride, int Radius)
 48 {
 49     int Channel = Stride / Width;
 50     if ((Src == NULL) || (Dest == NULL))                    return IM_STATUS_NULLREFRENCE;
 51     if ((Width <= 0) || (Height <= 0) || (Radius <= 0))        return IM_STATUS_INVALIDPARAMETER;
 52     if (Radius < 1)                                            return IM_STATUS_INVALIDPARAMETER;
 53     if ((Channel != 1) && (Channel != 3) && (Channel != 4))    return IM_STATUS_NOTSUPPORTED;
 54 
 55     Radius = IM_Min(IM_Min(Radius, Width - 1), Height - 1);        //    由于镜像的需求,要求半径不能大于宽度或高度-1的数据
 56 
 57     int SampleAmount = (2 * Radius + 1) * (2 * Radius + 1);
 58     float Inv = 1.0 / SampleAmount;
 59 
 60     int *ColValue = (int *)malloc((Width + Radius + Radius) * (Channel == 1 ? Channel : 4) * sizeof(int));
 61     int *ColOffset = (int *)malloc((Height + Radius + Radius) * sizeof(int));
 62     if ((ColValue == NULL) || (ColOffset == NULL))
 63     {
 64         if (ColValue != NULL)    free(ColValue);
 65         if (ColOffset != NULL)    free(ColOffset);
 66         return IM_STATUS_OUTOFMEMORY;
 67     }
 68     for (int Y = 0; Y < Height + Radius + Radius; Y++)
 69         ColOffset[Y] = IM_GetMirrorPos(Height, Y - Radius);
 70 
 71     if (Channel == 1)
 72     {
 73         for (int Y = 0; Y < Height; Y++)
 74         {
 75             unsigned char *LinePD = Dest + Y * Stride;
 76             if (Y == 0)
 77             {
 78                 memset(ColValue + Radius, 0, Width * sizeof(int));
 79                 for (int Z = -Radius; Z <= Radius; Z++)
 80                 {
 81                     unsigned char *LinePS = Src + ColOffset[Z + Radius] * Stride;
 82                     for (int X = 0; X < Width; X++)
 83                     {
 84                         ColValue[X + Radius] += LinePS[X];                                            //    更新列数据
 85                     }
 86                 }
 87             }
 88             else
 89             {
 90                 unsigned char *RowMoveOut = Src + ColOffset[Y - 1] * Stride;                //    即将减去的那一行的首地址    
 91                 unsigned char *RowMoveIn = Src + ColOffset[Y + Radius + Radius] * Stride;    //    即将加上的那一行的首地址
 92                 for (int X = 0; X < Width; X++)
 93                 {
 94                     ColValue[X + Radius] -= RowMoveOut[X] - RowMoveIn[X];                                            //    更新列数据
 95                 }
 96             }
 97             FillLeftAndRight_Mirror_C(ColValue, Width, Radius);                //    镜像填充左右数据
 98             int LastSum = SumofArray_C(ColValue, Radius * 2 + 1);                //    处理每行第一个数据                                
 99             LinePD[0] = IM_ClampToByte(LastSum * Inv);
100             for (int X = 1; X < Width; X++)
101             {
102                 int NewSum = LastSum - ColValue[X - 1] + ColValue[X + Radius + Radius];
103                 LinePD[X] = IM_ClampToByte(NewSum * Inv);
104                 LastSum = NewSum;
105             }
106         }
107     }
108     else if (Channel == 3)
109     {
110 
111     }
112     else if (Channel == 4)
113     {
114 
115     }
116     free(ColValue);
117     free(ColOffset);
118     return IM_STATUS_OK;
119 }

   
不过,过多的拍卖时间严重的限制住了中值滤波器的应用。由于其算法的非线性和不可分离性普通的优化技术并不适于。最原始的步骤就是获得图像像素的一个列表,然后开展排序,接着取中值。在相似的情状下,该算法的复杂度是O(r2logr),其中r为核的半径。当像素的可能取值是个常数时,比如对于8位图像,可以行使桶排序(巴克(Buck)er
sort),该排序使得算法的复杂度下跌为O(r2)。不过除此之外小半径的图景外,那样的改正任然是不行接受的。

     
 具体的进程本身用代码加以证实:

  之前没察觉到,就这么不难的代码用C写后能爆发速度也是很诱人的,3000*2000的图能不辱任务39ms,若是在编译选项里勾选编译器的“启用增强指令集:流式处理
SIMD 增加 2 (/arch:SSE2)”,
则系统会对上述所有浮点总计的有些应用相关的SIMD指令优化,如下图所示:

   
那里插一句,从自家个人的咀嚼上说,任何依照排序的中值滤波,都是不可以对大半径举办实时有效的拍卖的。

    1、一些公用的内存分配进度

                       
  图片 1

    在<A Fast
Two-Dimensional Median Filtering
Algorithm
>一文中,提出了根据直方图计算的迅猛中值滤波,其时间复杂度为O(r),这几个想法我在自身的Imageshop软件里也早就进行过(那些时候可不曾看过这么些随笔,可知这么些优化是个东风标致都能想到的一步)。后来意识,在Paint.net里的中值用的也是其一。具体来讲,那个算法通过总括在核内部的像素的直方图直接得到中值。当从一个像素移动到其余一个与之紧邻(水平或垂直方向)的像素时,只须要立异一部分的信息,如图1所示。为了立异直方图,2r+1次加法以及2r+1次减法要求实践,而从直方图中统计中值所急需的时间是肯定的,如代码段1所示。对于8位图像,直方图由256个元素构成,在平均上说,总计中值需求128次相比和127次加法。实际上,通过转移终止寻找的尺度大家可以总括任何其余百分比效果(见代码段1中的Percentile参数)。

    TMatrix *Row = NULL, *Col = NULL;
    unsigned char *LinePS, *LinePD;
    int  X, Y, K, Width = Src->Width, Height = Src->Height;
    int *RowOffset, *ColOffSet;

    unsigned short *ColHist    = (unsigned short *)IS_AllocMemory(256 * (Width + 2 * Radius) * sizeof(unsigned short), true);    
    if (ColHist == NULL) {Ret = IS_RET_ERR_OUTOFMEMORY; goto Done8;}
    unsigned short *Hist    = (unsigned short *)IS_AllocMemory(256 * sizeof(unsigned short), true);    
    if (Hist == NULL) {Ret = IS_RET_ERR_OUTOFMEMORY; goto Done8;}
    Ret = GetValidCoordinate(Width, Height, Radius, Radius, Radius, Radius, Edge, &Row, &Col);        //    获取坐标偏移量
    if (Ret != IS_RET_OK) goto Done8;

  那些时候3000*2000的图能不负众望25ms,,基本上接近我立异的OPENCV的代码的速度了。

         图片 2         
            图片 3

  其中的ColHist用于保存一行像素对应的列直方图
,注意这里的行是用的扩张后的行的轻重即:Width + 2 *
Radius。IS_AllocMemory是个里头使用了_mm_malloc定义的内存分配函数,重即使考虑SSE函数的16字节对齐问题。

  简单的叙述下各函数的职能先。

                                      图1
 黄氏算法,本图半径r=2

     
Hist变量用于保存每个像素点的局地直方图数据,任何根据局部直方图技术的函数最后都演化为对此该函数进行多种多样的盘算。
   

  IM_GetMirrorPos函数紧如果获取某一个岗位Pos根据镜像的方法处理时在Length倾向的坐标,那里Pos可以为负值,这一个首如果用来获得前期的坐标偏移的。      

int StopAmount = (2 * Radius + 1) * (2*Radius +1) * Percentile;
int Sum = 0;
for (int I = 0; I <= 255; I++)
{
    Sum += HistBlue(I);
    if (Sum >= StopAmount)      // 满足停止条件
    {
        Value = I;              // 得到中值
        break; 
    }
}

     
GetValidCoordinate是一个用于支援边界处像素点处理的函数,具体可详细附件中提交的代码。

  FillLeftAndRight_Mirror_C重如若用来获取两边镜像数据的(直接得到,不调用IM_GetMirrorPos函数),比如比如Array原始数据为
***abcdefgh*** (Length = 8, Radius =
3),则结果为dcbabcdefghgfe。

                                              代码1:  
从直方图中计算中值(Percentile=0.5)或擅自百分比率

      2、更新一行像素的列直方图

  SumofArray_C首假诺计量一个数组的富有的因素的和。

  为了使中值滤波的时日复杂性降低至线性以下,人们做出了诸多矢志不渝。Gel使用了基于树的算法将复杂度下跌为O(log2r),在同样篇杂文中,他们差不离的说了一种复杂度为O(log
r)的二维图像中值算法。大家那边指出的算法也是用直方图的乘除来代表如龟速的排序。如今,Weiss使用层次直方图技术取得了O(log
r)复杂度的职能,在他的章程中,纵然速度是快了,然而代码复杂难懂。本文提议一种简易而使得的算法,易于用CPU或其他硬件完结。

for (Y = 0; Y < Height; Y++)
{
    if (Y == 0)                                            //    第一行的列直方图,要重头计算
    {
        for (K = -Radius; K <= Radius; K++)                    
        {
            LinePS = Src->Data + ColOffSet[K] * Src->WidthStep;
            for (X = -Radius; X < Width + Radius; X++)
            {
                ColHist[X * 256 + LinePS[RowOffset[X]]]++;
            }
        }
    }
    else                                                //    其他行的列直方图,更新就可以了
    {
        LinePS = Src->Data + ColOffSet[Y - Radius - 1] * Src->WidthStep;        
        for (X = -Radius; X < Width + Radius; X++)        // 删除移出范围内的那一行的直方图数据
        {
            ColHist[X * 256 + LinePS[RowOffset[X]]]--;
        }

        LinePS = Src->Data + ColOffSet[Y + Radius] * Src->WidthStep;
        for (X = -Radius; X < Width + Radius; X++)        // 增加进入范围内的那一行的直方图数据
        {
            ColHist[X * 256 + LinePS[RowOffset[X]]]++;
        }
    }
  //  依次获取一行每个像素的局部直方图
   //  根据局部直方图获的结果
}

  IM_BoxBlur_C函数内部即为模糊的函数体,选择的优化思路完全和自由半径中值滤波(扩充至百分比滤波器)O(1)时间复杂度算法的原理、完毕及效能是一样的。当半径为R时,方框模糊的值等于以某点为主干,左右内外各扩大R像素的的限量内装有像素的综合,像素总个数为(2*R+1)*(2*R+1)个,当然大家也足以把他分成(2*R+1)列,每列有(2*R+1)个元素本例的优化措施大家就是把累加数据分为一列一列的,丰裕利用重复音信来落成速度提高。

 
  为更好的了然小说的算法,我们先来探视黄氏算法的阙如。更加注意到该算法行与行之间从未其它信息得到保留,而各种像素的处理至少有2r+1次加法和减法的直方图总结,那就是其复杂度为O(r)的原由。凭直觉,大家思疑应该有某种格局使得对各类像素,直方图只需加上一个定点的次数,从而获取O(1)的复杂度。正如我辈所看到的,通过保留行与行之间的新闻,这变得实惠。首先让大家来介绍下直方图的一部分特性上。对于不相邻的区域A和B,有下式制造:

  可知,这一部分和寻常的片段优化措施接近,没有何样出格的地点。

  我们定义一个Radius + Width + Radius的内存数据ColValue,用来保存列方向的累加数据,对于第一行数据,需要做特殊处理,也就是用原始的方式计算一行像素所有元素在列方向上的和,
详见78至于86行代码,当然这里只计算了中间Width范围内的数据,当不是第一行时,如下图左边所示,新的累加值只需减去移出的哪一行像素值同时加上移入的一行像素值,详见90到96
行代码。

  上面只计算了中间Width范围内的累加值,为了方便后续代码的编写以及使用SSE优化,下面的FillLeftAndRight_Mirror_C主要作用就是填充左边和右边分别填充数据,而且是按照镜像的方式进行填充。

                     
  H(A ∪ B) = H(A) + H(B)

  3、依次获得一行每个像素的一些直方图

       
在更新了上述累加值后,我们开首拍卖总括均值了,对于每行的率先个点,SumofArray_C总计了前2*R
+
1个列的累加值,那些累加值就是该点周边(2*R+1)*(2*R+1)个要素的累积和,对于一行的其余像素,其实就类似于行方向列累加值的换代,减去移出的加盟新进的列,如下图右边图所示,102到104行即落实了该进度。

  注意到八个直方图的增进是一个O(1)操作。他和直方图的因素个数有关,而直方图元素个数是由图像位深决定的。有了那或多或少,我们就足以支付一个新的算法。

    for (Y = 0; Y < Height; Y++)
    {
     //  更新一行像素的列直方图

        memset(Hist, 0, 256 * sizeof(unsigned short));        //    每一行直方图数据清零先
        LinePS = Src->Data + Y * Src->WidthStep;
        LinePD = Dest->Data + Y * Dest->WidthStep;
        for (X = 0; X < Width; X++)
        {
            if (X == 0)
            {
                for (K = -Radius; K <= Radius; K++)            //    行第一个像素,需要重新计算    
                    HistgramAddShort(ColHist + K * 256, Hist);
            }
            else
            {
              /*  HistgramAddShort(ColHist + RowOffset[X + Radius] * 256, Hist);    
                  HistgramSubShort(ColHist + RowOffset[X - Radius - 1] * 256, Hist);
           */
                HistgramSubAddShort(ColHist + RowOffset[X - Radius - 1] * 256, ColHist + RowOffset[X + Radius] * 256, Hist);  //    行内其他像素,依次删除和增加就可以了
            }

        //  根据局部直方图获的结果

            LinePS++;
            LinePD++;
        }
    }

           图片 4   
 图片 5

   
首先,对于每一列图像,大家都为其有限支撑一个直方图(对于8位图像,该直方图有256个因素),在全体的处理进程中,那些直方图数据都不可以不得到维护。每列直方图累积了2r+1个垂直方向上相邻像素的音信,开始的时候,那2r+1个像素是分别以率先行的各类像素为主导的。核的直方图通过积攒2r+1个相邻的列直方图数据得到。其实,大家所做的就是将核直方图分解成他对应的列直方图的集纳,在漫天滤波的经过中,这几个直方图数据在多个步骤内用恒定的时间保持最新。

  上面处理的进度实际上和2的进程的优化道理是近似的,只可是一个是行方向,一个是列方向,聪明者自然能清楚,稍微愚昧者请自己多么琢磨,自然有出现转机的每天。

  原理基本上就是这样,这个算法占用的额外内存很明显很少,但是不支持In-Place操作。
  为了追求速度,我们把整个过程能用SSE优化的地方都用SSE优化。
  首先是第79至86行的数据,这个很容易优化:

for (int Z = -Radius; Z <= Radius; Z++)
{
    unsigned char *LinePS = Src + ColOffset[Z + Radius] * Stride;
    int BlockSize = 8, Block = Width / BlockSize;
    for (int X = 0; X < Block * BlockSize; X += BlockSize)
    {
        int *DestP = ColValue + X + Radius;
        __m128i Sample = _mm_cvtepu8_epi16(_mm_loadl_epi64((__m128i *)(LinePS + X)));
        _mm_storeu_si128((__m128i *)DestP, _mm_add_epi32(_mm_loadu_si128((__m128i *)DestP), _mm_cvtepi16_epi32(Sample)));
        _mm_storeu_si128((__m128i *)(DestP + 4), _mm_add_epi32(_mm_loadu_si128((__m128i *)(DestP + 4)), _mm_unpackhi_epi16(Sample, _mm_setzero_si128())));
    }
    //  处理剩余不能被SSE优化的数据
}

   
考虑从某个像素向右移动一个像素的事态。对于眼前行,核最左侧的列直方图首先要求更新,而此时该列的列直方图中的数据照旧以上一行对应地方很是像素为主导测算的。由此要求减去最上一个像素对应的直方图然后增加其下部一像素的直方图音讯。那样做的成效就是将列直方图数据回落一行。这一步很强烈是个0(1)操作,唯有三回加法和五次减法,而于半径r无关。

  4、 按照局地直方图获的结果

  用_mm_loadl_epi64加载8个字节数据到内存,并用_mm_cvtepu8_epi16将其转移为16位的__m128i变量,前边再把低4位和高4位的数额分别转换成32位的,然后用_mm_add_epi32增进,注意后边我改换函数用了两种不相同的方法,因为那边的数目绝对都是正数,因而也是能够利用_mm_cvtepi16_epi32和_mm_unpackhi_epi16组合Zero实现。

   
第二步更新核直方图,其是2r+1个列直方图之和。那是经过减去最左侧的列直方图数据,然后再添加第一步所处理的那一列的列直方图数据获得的。这一步也是个O(1)操作,如图2所示。如前所述,加法、减法以及统计直方图的中值的耗时都是部分凭借于图像位深的盘算,而于滤波半径无关。

  根据分裂的算法需要,结合局地直方图新闻来取得结果,比如最大值算法可以用如下格局获取:

  再来看看92到95行代码,这些也很简单。

                                   
  图片 6

    for (K = 255; K >= 0; K--)
    {
        if (Hist[K] != 0)
        {
            LinePD[X] = K;
            break;
        }
    }
int BlockSize = 8, Block = Width / BlockSize;
__m128i Zero = _mm_setzero_si128();
for (int X = 0; X < Block * BlockSize; X += BlockSize)
{
    int *DestP = ColValue + X + Radius;
    __m128i MoveOut = _mm_unpacklo_epi8(_mm_loadl_epi64((__m128i *)(RowMoveOut + X)), Zero);
    __m128i MoveIn = _mm_unpacklo_epi8(_mm_loadl_epi64((__m128i *)(RowMoveIn + X)), Zero);
    __m128i Diff = _mm_sub_epi16(MoveIn, MoveOut);                        //    注意这个有负数也有正数的,有负数时转换为32位是不能用_mm_unpackxx_epi16体系的函数
    _mm_storeu_si128((__m128i *)DestP, _mm_add_epi32(_mm_loadu_si128((__m128i *)DestP), _mm_cvtepi16_epi32(Diff)));
    _mm_storeu_si128((__m128i *)(DestP + 4), _mm_add_epi32(_mm_loadu_si128((__m128i *)(DestP + 4)), _mm_cvtepi16_epi32(_mm_srli_si128(Diff, 8))));
}
//  处理剩余不能被SSE优化的数据

                图 2  算法的两步执行
 (a)左侧的列直方图的立异通过减最上部和加最下边像素消息获取

     关于直方图累加的代码如下:

  那里也是两遍性处理8个像素,那里自己利用了其它一种转移技术来把8位转换为16位,但是前边的Diff因为有正有负,要更换为32位就非得采纳_mm_cvtepi16_epi32来转换,不能够用unpack体系组合函数来贯彻,因为unpack不会移动符号位,我在此间吃了好四次亏了。

                     
                                                    (b)
核直方图的翻新通过减最右侧和加最右侧列直方图音信得到

/// <summary>
/// 无符号短整形直方图数据相加,Y = X + Y, 整理时间2014.12.28; 
/// </summary>
/// <param name="X">加数。</param>
/// <param name="Y">被加数,结果保存于该数中。</param>
/// <remarks>使用了SSE优化。</remarks>
void HistgramAddShort(unsigned short *X, unsigned short *Y)
{
    *(__m128i*)(Y + 0)        =    _mm_add_epi16(*(__m128i*)&Y[0],        *(__m128i*)&X[0]);        //    不要想着用自己写的汇编超过他的速度了,已经试过了
    *(__m128i*)(Y + 8)        =    _mm_add_epi16(*(__m128i*)&Y[8],        *(__m128i*)&X[8]);
    *(__m128i*)(Y + 16)        =    _mm_add_epi16(*(__m128i*)&Y[16],    *(__m128i*)&X[16]);
    *(__m128i*)(Y + 24)        =    _mm_add_epi16(*(__m128i*)&Y[24],    *(__m128i*)&X[24]);
    *(__m128i*)(Y + 32)        =    _mm_add_epi16(*(__m128i*)&Y[32],    *(__m128i*)&X[32]);
    *(__m128i*)(Y + 40)        =    _mm_add_epi16(*(__m128i*)&Y[40],    *(__m128i*)&X[40]);
    *(__m128i*)(Y + 48)        =    _mm_add_epi16(*(__m128i*)&Y[48],    *(__m128i*)&X[48]);
    *(__m128i*)(Y + 56)        =    _mm_add_epi16(*(__m128i*)&Y[56],    *(__m128i*)&X[56]);
    *(__m128i*)(Y + 64)        =    _mm_add_epi16(*(__m128i*)&Y[64],    *(__m128i*)&X[64]);
    *(__m128i*)(Y + 72)        =    _mm_add_epi16(*(__m128i*)&Y[72],    *(__m128i*)&X[72]);
    *(__m128i*)(Y + 80)        =    _mm_add_epi16(*(__m128i*)&Y[80],    *(__m128i*)&X[80]);
    *(__m128i*)(Y + 88)        =    _mm_add_epi16(*(__m128i*)&Y[88],    *(__m128i*)&X[88]);
    *(__m128i*)(Y + 96)        =    _mm_add_epi16(*(__m128i*)&Y[96],    *(__m128i*)&X[96]);    
    *(__m128i*)(Y + 104)    =    _mm_add_epi16(*(__m128i*)&Y[104],    *(__m128i*)&X[104]);
    *(__m128i*)(Y + 112)    =    _mm_add_epi16(*(__m128i*)&Y[112],    *(__m128i*)&X[112]);
    *(__m128i*)(Y + 120)    =    _mm_add_epi16(*(__m128i*)&Y[120],    *(__m128i*)&X[120]);
    *(__m128i*)(Y + 128)    =    _mm_add_epi16(*(__m128i*)&Y[128],    *(__m128i*)&X[128]);
    *(__m128i*)(Y + 136)    =    _mm_add_epi16(*(__m128i*)&Y[136],    *(__m128i*)&X[136]);
    *(__m128i*)(Y + 144)    =    _mm_add_epi16(*(__m128i*)&Y[144],    *(__m128i*)&X[144]);
    *(__m128i*)(Y + 152)    =    _mm_add_epi16(*(__m128i*)&Y[152],    *(__m128i*)&X[152]);
    *(__m128i*)(Y + 160)    =    _mm_add_epi16(*(__m128i*)&Y[160],    *(__m128i*)&X[160]);
    *(__m128i*)(Y + 168)    =    _mm_add_epi16(*(__m128i*)&Y[168],    *(__m128i*)&X[168]);
    *(__m128i*)(Y + 176)    =    _mm_add_epi16(*(__m128i*)&Y[176],    *(__m128i*)&X[176]);
    *(__m128i*)(Y + 184)    =    _mm_add_epi16(*(__m128i*)&Y[184],    *(__m128i*)&X[184]);
    *(__m128i*)(Y + 192)    =    _mm_add_epi16(*(__m128i*)&Y[192],    *(__m128i*)&X[192]);
    *(__m128i*)(Y + 200)    =    _mm_add_epi16(*(__m128i*)&Y[200],    *(__m128i*)&X[200]);
    *(__m128i*)(Y + 208)    =    _mm_add_epi16(*(__m128i*)&Y[208],    *(__m128i*)&X[208]);
    *(__m128i*)(Y + 216)    =    _mm_add_epi16(*(__m128i*)&Y[216],    *(__m128i*)&X[216]);
    *(__m128i*)(Y + 224)    =    _mm_add_epi16(*(__m128i*)&Y[224],    *(__m128i*)&X[224]);    
    *(__m128i*)(Y + 232)    =    _mm_add_epi16(*(__m128i*)&Y[232],    *(__m128i*)&X[232]);
    *(__m128i*)(Y + 240)    =    _mm_add_epi16(*(__m128i*)&Y[240],    *(__m128i*)&X[240]);
    *(__m128i*)(Y + 248)    =    _mm_add_epi16(*(__m128i*)&Y[248],    *(__m128i*)&X[248]);
}

  接着是FillLeftAndRight_Mirror_C函数的优化,改写如下:

   上述的实际效果就是核直方图向右移动,而列直方图向下移动。在盘算中,每个像素只需访问一次,并且被添加到一个直方图数据中。这样,最终一步就是计量中值了,如代码段1所示,这也是一个O(1)操作。

  _mm_add_epi16足以一遍性达成16个short类型的数目标加法,比传统的add指令快了很多倍。

void FillLeftAndRight_Mirror_SSE(int *Array, int Length, int Radius)
{
    int BlockSize = 4, Block = Radius / BlockSize;
    for (int X = 0; X < Block * BlockSize; X += BlockSize)
    {
        __m128i SrcV1 = _mm_loadu_si128((__m128i *)(Array + Radius + Radius - X - 3));
        __m128i SrcV2 = _mm_loadu_si128((__m128i *)(Array + Radius + Length - X - 5));
        _mm_storeu_si128((__m128i *)(Array + X), _mm_shuffle_epi32(SrcV1, _MM_SHUFFLE(0, 1, 2, 3)));
        _mm_storeu_si128((__m128i *)(Array + Radius + Length + X), _mm_shuffle_epi32(SrcV2, _MM_SHUFFLE(0, 1, 2, 3)));
    }
    //  处理剩余不能被SSE优化的数据
}

   
综上所述,所有的单像素操作(包含更新列以及核直方图、计算中值)都是
O(1)操作。现在,大家爱慕来说说开端化操作,即经过积累前r行的多寡来测算列直方图以及过去r列直方图数据测算第四个像素点的核直方图。那些进度是个O(R)操作。其余,从一行移动到下一行也侵夺了别的一个O(R)操作(表示不知底)。但是,既然那么些发轫化操作对每行只举办三遍,绝对于任何计量的话,那么些可以忽略不计。

   
 由于_mm_add_epi16是这对短整形数据举行的拍卖,由此,一般意况下改指令所能处理的半径不可能超过127,如果要求大于127,则须求修改进程序中的short类型为int,同时须求使用_mm_add_epi32指令,这样程序的速度会有所下降。

  镜像就是倒序的进程,直接利用SSE的shuffle函数很有益完结。

                  图片 7

  经过测试,在本人的I5的台式机中,1024*768图像在直方图更新上所须要的平分之间约为30ms,比较有的算法的基本就是部分岁月(比如上述的求最大值),可能大部分耗时并不在这里。

  统计数组的拉长和优化也造福。

    针对8位灰度图像,我们对上述算法举行一下总计。

   
 附件的代码中有个完全的测试工程,并有本人眼前持有的TMatrix结构的完全代码,我之后的小说都将以改结构为依托举行处理。

int SumofArray_SSE(int *Array, int Length)
{
    int BlockSize = 8, Block = Length / BlockSize;
    __m128i Sum1 = _mm_setzero_si128();
    __m128i Sum2 = _mm_setzero_si128();
    for (int X = 0; X < Block * BlockSize; X += BlockSize)
    {
        Sum1 = _mm_add_epi32(Sum1, _mm_loadu_si128((__m128i *)(Array + X + 0)));
        Sum2 = _mm_add_epi32(Sum2, _mm_loadu_si128((__m128i *)(Array + X + 4)));
    }
    int Sum = SumI32(_mm_add_epi32(Sum1, Sum2));
    //  处理剩余不能被SSE优化的数据
    return Sum;
}

   
(1)、对核最左边的列直方图执行一遍加法。

   
 代码还共享了家常便饭处理的函数,我很自信必然值得朋友去上学的。

  使用八个__m128i
变量首假使为着丰富利用XMM寄存器,其中SumI32函数如下,紧假设为了总计__m128i
内八个整数的累加值。

   
(2)、对同一列直方图实施三次减法,去除多余的像素新闻。

   
 那种光景看重的算法都有一个很沉重的老毛病,就是不可以相互,把图像分段处理,也会造成过多开首化耗时。

//    Convert 4 32-bit integer values to 4 unsigned char values.
void _mm_storesi128_4char(__m128i Src, unsigned char *Dest)
{
    __m128i T = _mm_packs_epi32(Src, Src);
    T = _mm_packus_epi16(T, T);
    *((int*)Dest) = _mm_cvtsi128_si32(T);
}

   
(3)、将革新后的列直方图数据加到核直方图中,那举行了256次加法。

   
  代码下载地址:http://files.cnblogs.com/files/Imageshop/BaseFile.rar

  代码不表明。

   
(4)、将无济于事的列直方图数据从核直方图中减去,那亟需256次减法。

图片 8

  最后就是100到104行的代码的转换。

   
(5)、为找到核直方图的中值,平均要求128次比较和127次加法。

 

int BlockSize = 4, Block = (Width - 1) / BlockSize;
__m128i OldSum = _mm_set1_epi32(LastSum);
__m128 Inv128 = _mm_set1_ps(Inv);
for (int X = 1; X < Block * BlockSize + 1; X += BlockSize)
{
    __m128i ColValueOut = _mm_loadu_si128((__m128i *)(ColValue + X - 1));
    __m128i ColValueIn = _mm_loadu_si128((__m128i *)(ColValue + X + Radius + Radius));
    __m128i ColValueDiff = _mm_sub_epi32(ColValueIn, ColValueOut);                            //    P3 P2 P1 P0                                                
    __m128i Value_Temp = _mm_add_epi32(ColValueDiff, _mm_slli_si128(ColValueDiff, 4));        //    P3+P2 P2+P1 P1+P0 P0
    __m128i Value = _mm_add_epi32(Value_Temp, _mm_slli_si128(Value_Temp, 8));                 //    P3+P2+P1+P0 P2+P1+P0 P1+P0 P0
    __m128i NewSum = _mm_add_epi32(OldSum, Value);
    OldSum = _mm_shuffle_epi32(NewSum, _MM_SHUFFLE(3, 3, 3, 3));                              //    重新赋值为最新值
    __m128 Mean = _mm_mul_ps(_mm_cvtepi32_ps(NewSum), Inv128);
    _mm_storesi128_4char(_mm_cvtps_epi32(Mean), LinePD + X);
}

    
上述计算量看起来相比较多。但是,大部分操作都可并行处理,从而显然的骤降了处理时间。更主要的,还有不少优化措施能下跌统计量,如下所示。

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

  从前认为这一个算法难以使用SSE优化,首要难题就在此地,可是在上学了Opencv的积分图技术时,这一个题目就一蹴而就了,进一步分析还发现Opencv的代码写的并不到家,还有立异的空间,见上述代码,使用三遍对相同数据移位就足以收获丰裕,由P3
P2 P1 P0变为P3+P2+P1+P0 P2+P1+P0 P1+P0 P0。我们有点分析一下。          
  

    1、并行化

 

  __m128i ColValueDiff =
_mm_sub_epi32(ColValueIn, ColValueOut);
那句代码获得了移入和移出种类的差值,我们计为;  P3 P2 P1 P0
(高位到没有)由于算法的增进特性,假如说OldSum的原始值为A3 A3 A3
A3,那么新的NewSum就应该是:

  现代处理器提供了SIMD指令可以用来加速我们的算法。下面描述的操作半数以上都是对直方图数据进行加和减操作。通过MMX,SSE2或Altivec指令可以并行处理多少个直方图操作。为了在一条指令中做越多的直方图处理,直方图的只可以用16位的数据类型。因而,大旨的大大小小最大为216(相应的核的轻重缓急是(256/2-1)),对于一般的使用充分了。这一个范围并不只是对大家的算法:这是一种优化的手段而已。

    A3+P3+P2+P1+P0  A3+P2+P1+P0  A3+P1+P0  A3+P0;

 
 此外一个足以运作并行的地方就是从图像中读取数据以及将其增加到相应的直方图中。同上述交替更新列和核直方图差别的是,我们可以率先更新整行的列直方图。然后采纳SIMD指令,大家可以相互的翻新多列直方图数据(分裂列直方图的立异之间从未其余关系,所以好相互)。然后再像平时一样更新核直方图。

__m128i Value_Temp = _mm_add_epi32(ColValueDiff, _mm_slli_si128(ColValueDiff, 4)); 这句的_mm_slli_si128得到中间结果 P2 P1 P0 0, 再和P3 P2 P1 P0相加得到

    Value_Temp =  P3+P2   P2+P1  P1+P0   P0

__m128i Value = _mm_add_epi32(Value_Temp, _mm_slli_si128(Value_Temp, 8));这句的_mm_slli_si128得到中间结果P1+P0 P0 0 0,再和P3+P2 P2+P1 P1+P0  P0相加得到:

    Value = P3+P2+P1+P0   P2+P1+P0   P1+P0   P0
好简单的过程。

  OldSum = _mm_shuffle_epi32(NewSum, _MM_SHUFFLE(3, 3, 3, 3)); 这一句为什么要这样写,恐怕还是读者自己体会比较好,似乎不好用文字表达。

   
那条优化手段对于有些高级语言是力不从心兑现的。像VC6,VC.NET那类可以一贯内嵌汇编的语言,就算可以兑现,也亟需作者具有很好的汇编语言基础,由此,实施的难度相比较大。有趣味的读者可以参考附件中的中的SSE代码。

   最终一个_mm_storesi128_4char是本身自己定义的一个将1个__m128i其中的4个32位整数转换为字节数并保存到内存中的函数,详见附件代码。

   
 2、缓存优化

  至于24位图像的优化,是个相比较为难的境地,因为SSE一次性处理4个32位数,而24位各类像素只有3个轻重,那种景色相似如故把她恢弘到4个轻重,比如说ColValue就能够改成4坦途的,然后累积和也亟需处理成4大路的,那本来须求自然的精细淫技,那里自己不想多谈,留点东西给自己。当然也足以设想先把24位分解到3个灰度内存,然后利用灰度的算法进行统计,前边在合成到24位,那些解释和合成的进程也是能够行使SSE加快的,若是您结合多线程,还能把3个灰度进度的处理并行化,这样除了多占用内存外,速度比至二级处理24位要快(因为直接处理算法不可以并行的,前后器重的原因)。其余就是最后在盘算列累积求平均值的长河也变得更为自然了,不会油不过生灰度这种要在__mm128i里面举办添加的长河,而是一贯得五个SSE变量累加。

  恒常时间的中值滤波算法须要在内存中为每列保持一个直方图,对于图像,那很不难就多达数百KB的大小,经常那高于前几天的处理器的缓存。那造成访问内存的频率下跌。一种方式就是将图像在档次方向上分为几有些分离处理。每个块的升幅大小需细心选拔,
使得其相应的列直方图不超过缓存的轻重而有能丰硕利用缓存。那种修改的不利点就是扩展的边缘效应。在实践中,这一般能导致处理时间的大幅下落。请留意,在不一样的微机上同时处理那几个块是该算法的一种很粗略的并行算法。

  还说一点,现在一大半的CPU都辅助AVX256了,还是可以够动用AVX进一步加速,如同代码该起来也不是很难,有趣味的仇人可以协调试试。

   
这种优化说实在我不晓得怎么着用代码去完毕。

  可以说,如今那一个速度已经差不离达到了CPU的顶峰了,不过测试过IPP的速度,就如比那几个还要快点,不拔除他利用了AVX,也不排除他选择多核的资源。

  3、多层直方图

  这一个的优化对于BoxBlur来说是第一的一步,可是更主要的是她可以运用在无数场所,比如图像的一部分均方差总结,也可以动用类似的技巧举办加快,两幅图像大的一些平方差也是足以如此优化的,后续我会不难的谈下他在那上边加速的应用。

  在<A
Coarse-to-Fine Algo-rithm for Fast Median Filtering of Image Data With a
Huge Number
of Levels
>一文中显得了多规格直方图是充裕实用的优化手段。其想尽是维持一个平行的较小的直方图,直方图记录了图像的上位数据。例如,对于8位图像,使用两层的直方图很常用,其中高层使用4位,而低层使用全8位数据。习惯上大家分别给她们命名为粗分和剪切直方图。粗分直方图包蕴了16(2^4)个元素,每个元素是呼应的分开直方图高位的积累和。

  源代码下载:https://files.cnblogs.com/files/Imageshop/FastBlur.rar

   
使用多层直方图有四个好处,第四个就是一个钱打二十四个结中值进度的增速。我们得以率先在粗分数据中需找到中值在划分数据之中的职位而不用检查整个256个地方。平均上说那只必要16次而不是128次相比和相加。第一个便宜是有关直方图的相加和相减。当对应的粗分数据为0时,则能够不用统计对应的细分段。当半径
r很时辰,列直方图是稀疏分布的,那几个时候的分支判断是很有必不可少的。

  彩色图工程:https://files.cnblogs.com/files/Imageshop/SSE_Optimization_Demo.rar

   
以上说的很暧昧。举个简单的例子吗,对于3*3的一个核,假设像素的数据值分别为:
100、120、98、77、215、50、243、199、180。对应的要职种类为:6、7、6、4、13、3、15、12、11。低位种类为:4、8、2、13、7、2、3、7、4。
那么粗分直方图的16个要素的值分别为:

图片 9

     
Coarse[3]=1  Coarse[4]=1  Coarse[6]=2  Coarse[7]=1  Coarse[11]=1  Coarse[12]=1  Coarse[13]=1  Coarse[15]=1,其余都为0;

 

  中位数的累加值为3*3/2=5,对粗分直方图进行添加:Coarse[3]+Coarse[4]+Coarse[6]+Coarse[7]=5,获得中值对应的高位段是H=7,由此,中间值就相应在112-128中间,然后再从细分直方图的呼应段中遵从类似的不二法门找到低位的值L,最终取得中值H*16+L。

  不佳意思,图太小,速度为0ms了。

   
4、条件更新核

图片 10

   
粗分和剪切直方图的离别还有一个不明明可是很得力的优化效率。注意到算法的大部分时日都费用更新在核直方图的时抬高或减去列直方图数据,那几个日子随着实时更新粗分直方图而有条件的换代细分直方图而得到下降。

  

   
记得前边说过测算中值的历程是先在粗分数据中搜索中值所在段,然后再从细分数据中找到精确值。对于核的中值,每个列直方图最八只会有2r+1次进献,意味着唯有2r+1个照应的细分段对计量结果有用。那多少个一贯未被利用的段,其相应的分割数据将无需立异。

  

 
 为了促成该意义,大家要求为各样开辟一个记录其最终被更新的地方的列表。当从一个像素移向下个一个像素时,大家立异列直方图以及核直方图的粗分数据。然后按照粗分数据测算出中值再分开数据中所在的段。下一步,按照那几个段上次被更新的职位更新的分割直方图。要是上次更新的地方和方今列的岗位距离2r+1的偏离,那表明旧的职责和眼前地点没有其他交叉。因而大家从0开端更新段。正是经过那种艺术,我们加速了方方面面处理速度。

 

   
交错陈设直方图数据,从而使得相邻列的直方图数据在内存也是相邻的是有益处的。由此,细分直方图数据须求按下述格局布署:段索引、列索引、最终是因素索引。那样,核的划分直方图的翻新就是对一块接二连三的内存的增长了,具体的讲,细分直方图有接近如下的概念形式:int
[,,] Fine= new int
[16,Width,16],其中首个16对应段索引,即像素值的高4位,最后一个16是元素值,即像素的低4位值。因为三维数组的拜访会导致冗余的盘算下标的经过,由此,为了加强速度,应该运用一维数组要么直接用指针访问内存,以一维数组为例,此时的概念应该修改为int
[] Fine=new
int[16*Width*16],对于第X行某个像素值Value,其相应的细分地点则为:
 Fine[(Value>>4) * (Width*16)+X*16+ (Value & 0xF)];

 

   
具体的也许照旧看参考代码更易于通晓。

 

   
参考代码中如下函数:

 

static inline void histogram_sub( const uint16_t x[16], uint16_t y[16] )
{
    int i;
    for ( i = 0; i < 16; ++i ) {
        y[i] -= x[i];
    }
}

 

    第四个指出是把这一个小循环手工展开。第二,我是是用C#编程完毕结果的,C#从未inline,于是我把如此的代码直接举行内嵌到自家的代码中,然则令人惊呆的结果却是调用函数的本子速度却比内嵌的快慢还要快,反汇编函数版本代码如下:

  y[0] -= x[0];
00000000  push        ebp 
00000001  mov         ebp,esp 
00000003  push        esi 
00000004  mov         esi,ecx 
00000006  mov         ecx,edx 
00000008  mov         eax,dword ptr [esi] 
0000000a  sub         dword ptr [ecx],eax 
            y[1] -= x[1];
0000000c  lea         edx,[ecx+4] 
0000000f  mov         eax,dword ptr [esi+4] 
00000012  sub         dword ptr [edx],eax 
            y[2] -= x[2];
00000014  lea         edx,[ecx+8] 
00000017  mov         eax,dword ptr [esi+8] 
0000001a  sub         dword ptr [edx],eax 

...........................其他的减法.....................................

            y[14] -= x[14];
00000074  lea         edx,[ecx+38h] 
00000077  mov         eax,dword ptr [esi+38h] 
0000007a  sub         dword ptr [edx],eax 
            y[15] -= x[15];
0000007c  add         ecx,3Ch 
0000007f  mov         edx,ecx 
00000081  mov         eax,dword ptr [esi+3Ch] 
00000084  sub         dword ptr [edx],eax 
00000086  pop         esi 
        }
00000087  pop         ebp 
00000088  ret 

  关于那几个题材,我的解析是,写成函数的本子,纵然多了几句压栈和出栈的语句外,CPU会丰富利用寄存器来进展操作,而内嵌后,由于变量太多,CPU只可以采纳内存来处理那么些来处理那个赋值的。而寄存器比内存的操作快多了。由此不明白VC的内联函数会不会也有这问题。

   
关于算法的耗时景况,原文给出了一个图纸:

              图片 11

         我本机上安装的是PS
CS4,经过测试,其中间值算法的耗时也是随着用户输入的半径的加码而成线性增添的,可是在小半径的时候仍然比较快的。

 
   对于小半径,我的建议是应用黄算法的方式+多层直方图的艺术来兑现,速度会比本文要更快些。从理论上分析,而唯有在半径大于7时,可选用本文效果达到O(1)复杂度。

 
 图片 12  图片 13  图片 14

          原始图像                  半径=5,百分比=50               半径=40,百分比=50

 
 图片 15  图片 16  图片 17

           
  半径=5,百分比=25                     半径=5,百分比=75  
                         半径=40,百分比=75

     
以一副1024*768的24位真彩色图像为例,进行进程结果报告:

时间比较(ms)
      半径         本文算法        黄算法+多层直方图
      2            506         259
      4            512         323
      6             478         377
      8            469         452
     10            479          515
     20            467         1004
     50            483         2333
     100            525         4947

 

 

 

      

 

 

 

   

 

   
同样,提供个编译好的文本给有趣味商讨该算法的爱人看看效果:

 
  http://files.cnblogs.com/Imageshop/MedianFilter.zip

   
由于上述MedianFilter是用VB6编制的,而VB6不接济指针,只好用数组来操作图像数据,在有点地方会招致功效的惨重丢失,我用C#写的进度(未用三十六线程)是该MedianFilter的速度的三倍到4倍左右,假诺是处理灰度,基本上又可一达到平等大小彩色图像速度的2到3倍。

   
在骨子里运用中,大半径的中值对于去除椒盐噪音是尚未意思的,因为那时图像已经损失了太多立见作用的音讯了。根据我的问询,大半径可以发挥用处的地点有:1、若是您的顺序有和PS一样的选区技术,那么选区的平整这一个职能实在就是对选区数据开展中值处理的长河,那些当然希望之星速度和半径毫无干系。2、一些二值图像的去噪上得以动用,比如一定半径的中值,可以去除一些孤立的块,从而扫除噪声。3、在对一些图像进行艺术处理时,必要大半径的中值。

     2013.10.16 补充:  

     近日上马上学了下C,于是用C语言落成上述进度。而在随想小编的原始代码中,用到SSE2的相干函数举行处理,即有如下的代码:

__inline void  HistgramAdd( unsigned short *x, unsigned short *y  )
{
    *(__m128i*) y = _mm_add_epi16( *(__m128i*) y, *(__m128i*) x );
    *(__m128i*) (y+8) = _mm_add_epi16( *(__m128i*) &y[8], *(__m128i*) &x[8] );
}

__inline void  HistgramSub(unsigned short *x, unsigned short *y )
{
    *(__m128i*) &y[0] = _mm_sub_epi16( *(__m128i*) &y[0], *(__m128i*) &x[0] );
    *(__m128i*) &y[8] = _mm_sub_epi16( *(__m128i*) &y[8], *(__m128i*) &x[8] );
}

  这里_mm_add_epi16是一组Instructions一声令下,编译器在编译的时候会将其编译为SSE对应的汇编代码,而鉴于那几个指令能促成指令级别并行,比如上述_mm_add_epi16方可在同一个指令周期对8个16位数据同时举办拍卖,并且HistgramAdd那几个函数在先后里会大方利用到,由此先后的进程能大幅升高。

   
对于上述测试的1024*768的图片,用C编写DLL,C#调用,同样的机械耗时平均在160ms左右,速度较纯粹的C#的代码快了3倍多,可知SSE的强劲。 

[DllImport("MedianBlur.dll", CallingConvention = CallingConvention.StdCall, CharSet = CharSet.Unicode, ExactSpelling = true)]
private static extern void MedianBlurGray( byte * Src, byte * Dest, int Width, int Height ,int Stride ,int Radius,int Percent);
[DllImport("MedianBlur.dll", CallingConvention = CallingConvention.StdCall, CharSet = CharSet.Unicode, ExactSpelling = true)]
private static extern void MedianBlurRGB( byte * Src, byte * Dest,int Width, int Height ,int Stride ,int Radius,int Percent);

 

    
附件为C#调用C编写的DLL的工程代码文件:http://files.cnblogs.com/Imageshop/MedianBlurTest.rar

      图片 18
     

     
由于_mm_add_epi16是指向16位的数据类型进行的处理,所以中值得半径一般须求不当先128,否则数出现数量溢出等张冠李戴,工程中如此大的半径已经足足应付半数以上场合的。

 

 

 ***************************小编:
laviewpbt   时间: 2013.4.26    联系QQ:  33184777
 转发请保留本行新闻*************************

 

 

 

相关文章