超高速指数模糊算法的实现和优化（10000*10000在100ms左右实现）。

   static inline void _blurinner(guchar* pixel, gint* zR,gint* zG,gint* zB, gint* zA, gint alpha,gint aprec,gint zprec)
  {
    gint R, G, B, A;
    R = *pixel;
    G = *(pixel + 1);
    B = *(pixel + 2);
    A = *(pixel + 3);

    *zR += (alpha * ((R << zprec) - *zR)) >> aprec;
    *zG += (alpha * ((G << zprec) - *zG)) >> aprec;
    *zB += (alpha * ((B << zprec) - *zB)) >> aprec;
    *zA += (alpha * ((A << zprec) - *zA)) >> aprec;

    *pixel       = *zR >> zprec;
    *(pixel + 1) = *zG >> zprec;
    *(pixel + 2) = *zB >> zprec;
    *(pixel + 3) = *zA >> zprec;
  }

static inline void _blurcol(guchar* pixels,
                               gint    width,
                               gint    height,
                               gint    channels,
                               gint    x,
                               gint    alpha,
                               gint    aprec,
                               gint    zprec)
  {
    gint zR;
    gint zG;
    gint zB;
    gint zA;
    gint index;
    guchar* ptr;

    ptr = pixels;

    ptr += x * channels;

    zR = *((guchar*) ptr    ) << zprec;
    zG = *((guchar*) ptr + 1) << zprec;
    zB = *((guchar*) ptr + 2) << zprec;
    zA = *((guchar*) ptr + 3) << zprec;

    for (index = width; index < (height - 1) * width; index += width)
      _blurinner((guchar*) &ptr[index * channels], &zR, &zG, &zB, &zA, alpha,
      aprec, zprec);

    for (index = (height - 2) * width; index >= 0; index -= width)
      _blurinner((guchar*) &ptr[index * channels], &zR, &zG, &zB, &zA, alpha,
      aprec, zprec);
  }

  // pixels   image-data
  // width    image-width
  // height   image-height
  // channels image-channels
  // in-place blur of image ‘img‘ with kernel of approximate radius ‘radius‘
  // blurs with two sided exponential impulse response
  // aprec = precision of alpha parameter in fixed-point format 0.aprec
  // zprec = precision of state parameters zR,zG,zB and zA in fp format 8.zprecb

  void _expblur(guchar* pixels,
                 gint    width,
                 gint    height,
                 gint    channels,
                 gint    radius,
                 gint    aprec,
                 gint    zprec)
  {
    gint alpha;
    gint row = 0;
    gint col = 0;

    if (radius < 1)
      return;

    // calculate the alpha such that 90% of 
    // the kernel is within the radius.
    // (Kernel extends to infinity)
    alpha = (gint) ((1 << aprec) * (1.0f - expf(-2.3f / (radius + 1.f))));

    for (; row < height; row++)
      _blurrow(pixels, width, height, channels, row, alpha, aprec, zprec);

    for (; col < width; col++)
      _blurcol(pixels, width, height, channels, col, alpha, aprec, zprec);

    return;
  }

for (int Y = 0; Y < Height; Y++)
{
    byte *LinePS = Src + Y * Stride;
    byte *LinePD = Dest + Y * Stride;
    int Sum = LinePS[0] << Zprec;
    for (int X = 0; X < Width; X++)　　　　　　//　　从左往右
    {
        Sum += (Alpha * ((LinePS[X] << Zprec) - Sum)) >> Aprec;
        LinePD[X] = Sum >> Zprec;
    }
    for (int X = Width - 1; X >= 0; X--)　　　//　　从右到左
    {
        Sum += (Alpha * ((LinePD[X] << Zprec) - Sum)) >> Aprec;
        LinePD[X] = Sum >> Zprec;
    }
}

int *Buffer = (int *)malloc(Width * sizeof(int));
for (int X = 0; X < Width; X++)        Buffer[X] = Src[X] << Zprec;
for (int Y = 0; Y < Height; Y++)
{
    byte *LinePS = Src + Y * Stride;
    byte *LinePD = Dest + Y * Stride;
    for (int X = 0; X < Width; X++)　　　　　　　　//　　从上到下
    {
        Buffer[X] += (Alpha * ((LinePS[X] << Zprec) - Buffer[X])) >> Aprec;
        LinePD[X] = Buffer[X] >> Zprec;
    }
}
for (int Y = Height - 1; Y >= 0; Y--)　　　　　　//　　从下到上
{
    byte *LinePD = Dest + Y * Stride;
    for (int X = 0; X < Width; X++)
    {
        Buffer[X] += (Alpha * ((LinePD[X] << Zprec) - Buffer[X])) >> Aprec;
        LinePD[X] = Buffer[X] >> Zprec;
    }
}
free(Buffer);

int X = 0;
for (X = 0; X < Width - 8; X += 8)            
{
    //    将8个字节数存入到2个XMM寄存器中
    //    方案1：使用SSE4新增的_mm_cvtepu8_epi32的函数，优点是两行是独立的
    __m128i Dst1 = _mm_cvtepu8_epi32(_mm_cvtsi32_si128((*(int *)(LinePD + X + 0))));    //    _mm_cvtsi32_si128把参数中的32位整形数据移到XMM寄存器的最低32位，其他为清0。
    __m128i Dst2 = _mm_cvtepu8_epi32(_mm_cvtsi32_si128((*(int *)(LinePD + X + 4))));    //    _mm_cvtepu8_epi32将低32位的整形数的4个字节直接扩展到XMM的4个32位中去。
    __m128i Buf1 = _mm_loadu_si128((__m128i *)(Buffer + X + 0));
    __m128i Buf2 = _mm_loadu_si128((__m128i *)(Buffer + X + 4));
    Buf1 = _mm_add_epi32(_mm_srai_epi32(_mm_mullo_epi32(_mm_sub_epi32(_mm_slli_epi32(Dst1, Zprec), Buf1), Alpha128), Aprec), Buf1);
    Buf2 = _mm_add_epi32(_mm_srai_epi32(_mm_mullo_epi32(_mm_sub_epi32(_mm_slli_epi32(Dst2, Zprec), Buf2), Alpha128), Aprec), Buf2);
    _mm_storeu_si128((__m128i *)(Buffer + X + 0), Buf1);
    _mm_storeu_si128((__m128i *)(Buffer + X + 4), Buf2);
    _mm_storel_epi64((__m128i *)(LinePD + X), _mm_packus_epi16(_mm_packs_epi32(_mm_srai_epi32(Buf1, Zprec), _mm_srai_epi32(Buf2, Zprec)), Zero));
}
for (; X < Width; X++)        
{
    Buffer[X] += (Alpha * ((LinePD[X] << Zprec) - Buffer[X])) >> Aprec;
    LinePD[X] = Buffer[X] >> Zprec;
}

void ExpFromLeftToRight_OneLine_SSE(int *Data, int Length, int Radius, int Aprec, int Zprec, int Alpha)
{
    int *LinePD = Data;
    __m128i A = _mm_set1_epi32(Alpha);
    __m128i S1 = _mm_slli_epi32(_mm_load_si128((__m128i *)(LinePD)), Zprec);
    for (int X = 0; X < Length; X++, LinePD += 4)
    {
        S1 = _mm_add_epi32(S1, _mm_srai_epi32(_mm_mullo_epi32(_mm_sub_epi32(_mm_slli_epi32(_mm_load_si128((__m128i *)(LinePD)), Zprec), S1), A), Aprec));
        _mm_store_si128((__m128i *)(LinePD), _mm_srai_epi32(S1, Zprec));
    }
}